party for prediction [REPOST]

This was an exceptionally helpful answer, I can only thank you again.
I have plenty of avenues ahead where I was worried before I was
getting trapped in a dead end. If all else fails, the idea of using
anova is brilliant. Thank you!

Ed

On 14 October 2012 18:36, Achim Zeileis <Achim.Zeileis at uibk.ac.at> wrote:

On Sun, 14 Oct 2012, Ed wrote:

First up, thanks hugely for your response. I've been beating my head
against this!

On 14 October 2012 16:51, Achim Zeileis <Achim.Zeileis at uibk.ac.at> wrote:

I'm not sure what you mean by "integral vector". If you want to apply the
approach to hundreds of thousands of observations, I gues that these are
categorical (maybe even binary?) but maybe not...

I'm sorry I can't go into the details of the data, I would if I could.
z are categorical variables represented as integers, mostly ordered,
but not all. I've tried fitting them as integers, as well as ordered,
but O don't think it made a huge difference.

The tests performed for categorical partitioning variables are rather
different from the tests for numerical partitioning variables. If all of the
variables are categorical, this may not be immediately obvious, but the
factor coding should be more appropriate (especially if the number of levels
is small or moderate).

If I recall correctly, we kept linearModel as simple as we did to save as
much time as possible. This can be particularly important when one of the
partitioning variables has many possible splits and the linearModel has
to
be fitted thousands of times.

I can appreciate that, but maybe having an alternative linearModel which
will predict when the fit is degenerate would be worth including? I'm happy
to contribute what I have, although it's pretty obvious stuff (and probably
done suboptimally since I'm not much of an R coder at this point). For me at
least, even with huge datasets, the speed of party is quite good; it's
getting a better result that's the problem.

As I explained in my last e-mail. In your situation this does not solve the
problem completely because subsequent the tests are also not adapted to
this. Setting the empirical estimating functions to zero for non-identified
coefficients might alleviate the problem but is not really a clean solution.

Also, mob() assesses the stability of all coefficients of the model in
all
nodes during partitioning. If any of the coefficients is not identified,
this would have to be excluded from all subsequent parameter stability
tests
in that node (and its child nodes). This is currently not provided for in
mob().

Would pretending the coefficients were fit at 0 fool mob into doing
something moderately meaningful here?

The coefficients are not looked at during fitting, only the estfun(). This
would have to be set to 0.

If not, I would try to hack the code, but I'm honestly at something of
a loss as to how to modify it and feed the results back into my
interpreter. I have bytecode installed; I downloaded the source, but I
haven't squared the circle of modifying the source and installing the
result. I will check out the docs on writing extensions you suggest.

Writing (or modifying) R packages and installing them under Windows is
pretty standard and well documented. The pointers I gave you should
hopefully get you started.

The second problem I have is that I get "Cholesky not positive definite"
errors at some nodes. I guess this is because of numerical error and
degeneracy in the covariance matrix? Any thoughts on how to avoid having
this happen would be welcome; it is ignorable though for now.

This comes from the parameter stability tests and might be a result of an
unidentified (or close to unidentified) model fit.

This is a great help to know. I improved my results quite considerably
with aggressive scaling of everything (scaling the response and all
the predictors to lie between 0 an 1). That deepened my tree by a
factor of two or so (say depth 3 to 7) and improved the quality of fit
substantially. Is there any way I can engage a more numerically robust
Cholesky in mob?

No, I don't think that this is conceivable with the way this is implemented
at the moment. Instead of the currently implemented tests, one could in
principle use the likelihood ratio test which is invariant against parameter
transformations and doesn't need the covariance matrix.

With hundreds of thousands of observations, you would need some
additional
pruning strategy anyway. Significance test-based splitting will probably
overfit because tiny differences in the coefficients will be picked up at
such large sample sizes.

I'm okay with overfitting, honestly. At the moment it is underfitting
by quite a large amount I think (the quality of the predictions on the
training set is not very high). The problem really is there is so much
going on the data, but the "noise" level is probably very low. I
wouldn't be surprised if my data was accurate to 5 or 6 s.f.

OK. Maybe with so little noise some other splitting strategy (not based on
significance tests) would be better?

Furthermore, computationally the extensive search over all possible
splits
might be too burdensome with this many observations.

I have plenty of compute time/power and RAM, though R seems to be
running single threaded. But even on a few million observations, it's
still pretty fast and doesn't use more than 30 or so gig of memory. If
it takes a day and requires 150gig of RAM, that is absolutely fine,
even over that would be viable though less optimal.

With this setup, you may consider writing your own partitioning algorithm
using the same type of ideas as MOB. Instead of using the parameter
stability tests, you could use plain likelihood ratio tests or ANOVAs to
decide about splitting further.

If d is the data in the current node, y and x are response/regressors, and
f1 to fn are your categorical partitioning variables, you could do

m0 <- lm(y ~ x, data = d)
m1 <- lm(y ~ f1 * x, data = d)
...
mn <- lm(y ~ fn * x, data = d)

and then

anova(m0, m1)
...
anova(m0, mn)

and then choose the lowest p-value for splitting. You could then also
parallelize the computations in the daughter nodes.

None of this is readily coded in party though...

Hence, using some subsampling strategy might not be the worst thing.

I've tried this at various degrees, but the data is really very
complicated with not a lot of error. I'm trying to encourage party to
fit more closely, which I thought more data might encourage. At the
moment I'm a long way from a clean fit. I have subsampled at various
levels down to 1%, and although that increases the depth of the tree
and quality of fit, it still doesn't give a very good quality fit and
can encourage it to overlook obvious aspects of the training set.

I guess what I really want to know is:
(a) has anyone else had this problem, and if so how did they overcome
it?

We have had non-identified model fits in binary GLMs (with quasi-complete
separation) where we then set estfun() to all zero so that partitioning
stops. But I don't think that such a strategy helps here.

I've considered using rpart() to partition into cells of constant
gradient, then fitting linear models myself to the cells. This is my
next thought. I'm pretty sure partitioning over linear regression is
the way forward for the data we have. I tried mars and glm but there
are good reasons to think they're less reasonable, even though the fit
wasn't particularly poor. I'm not particularly wedded to party's
approach except that it looked like it immediately returned what we
needed, and with some degree of "optimality" into the bargain.

Possibly boosting linear models may be another route to go.

(b) is there any way to get a line or stack trace out of a try()
without source modification?

Not sure, I don't know any off the top off my head.

I guess I really will have to bite the bullet and try to figure out
how to install modified libraries. Thanks.

I think that's easier compared to solving the conceptual problems (thinking
about how to best partition degenerate linear models etc.).
;-)

Good luck!
Best wishes,

Z

(c) failing all of that, does anyone know of an alternative to mob
that does the same thing; for better or worse I'm now committed to
recursive partitioning over linear models, as per mob?

If your partitioning variables are particularly simple (e.g., all binary)
you could exploit that and it may be easier to write a custom function
for
your particular data. Then likelihood-ratio tests (rather than LM-type
tests) would also be easier to apply in case of unidentified parameters.

But if there are partitioning variables with different measurement
scales,
then this will not be that simple...

Unfortunately each partitioning variable is essentially a state
indicator, taking values say 0,...,R where R is different for each
component. I'm not a stats expert either; I've spent some time with
the party manuals and papers, but I wouldn't be confident of
implementing something like it in the time available to me (though if
I have to I will, but that wouldn't be a good situation to be in).

(d) failing all of this, does anyone have a link to a way to rebuild, or
locally modify, an R package (preferably windows, but anything would
do)?

Have a look at the "Writing R Extensions" manual and the R for Windows
FAQ.

Will do.

Thank you very much for your responses, I really appreciate it.

Best wishes,

Ed