Hi Will,
Your edf interpretation is not quite right. The smooths are subject to a
centring constraint for identifiability reasons, and this removes a
degree of freedom, so EDF=1 corresponds to a straight line fit.
On your second point and third points. anova(gamb1.1,gamb1.2,
test=Chisq) conducts the crude approximate GLRT described in the
?anova.gam. It will not produce sensible results for comparing two
models that are essentially identical (nor, although it's not the issue
here, does it do well for testing random effects for equality to zero).
mgcv generally tries not to second guess the user, so it won't stop you
trying to compare two essentially identical models, but that doesn't
mean that the results will be meaningful.
Basically, the results of anova(gamb1.1,gamb1.2, test=Chisq) will give
reasonable ball-park p-values when both models can be reasonably well
approximated by some un-penalized model (approximately matching the
effective degrees of freedom of the penalized version), and the null
model at least has a lower EDF than the alternative (better still, one
would condition on the larger model smoothing parameter estimates, to
force proper nesting). The test is in any case conditional on the
smoothing parameter estimates. This can lead to the p-values being
somewhat too low when smoothing psrameters are highly uncertain (as they
may be at small sample sizes). Note that ML or REML smoothing parameter
estimation (see gam argument 'method') tends to lead to best p-value
performance in simulations.
Generally speaking, the p-values reported by single model calls to
`anova.gam' or `summary.gam' are better justified than the approximate
GLRT (but prior to 1.7-19, only for terms which cannot be penalized away
altogether). Even then, what is reported is conditional on the smoothing
parameter estimates, so can be too low when these are not well identified.
best,
Simon
On 14/07/12 17:41, wshadish wrote:
We are using GAM in mgcv (Wood), relatively new users, and wonder if anyone
can advise us on a problem we are encountering as we analyze many short time
series datasets. For each dataset, we have four models, each with intercept,
predictor x (trend), z (treatment), and int (interaction between x and z).
Our models are
Model 1: gama1.1 - gam(y~x+z+int, family=quasipoisson) ##no smooths
Model 2: gama1.2 - gam(y~x+z+s(int, bs=cr), family=quasipoisson) ##smooth
the interaction
Model 3: gama1.3 - gam(y~s(x, bs=cr)+z+int, family=quasipoisson) ##smooth
the trend
Model 4: gama1.4 - gam(y~s(x, bs=cr)+z+s(int, bs=cr),
family=quasipoisson) ##smooth trend and interaction
We have three questions. One question is simple. We occasionally obtain edf
=1 and Ref.df=1 for some smoothed predictors (x, int). Because Wood says
that edf can be interpreted roughly as functional form (quadratic, cubic
etc) + 1, this would imply x^0 functional form for the predictor, and that
doesn't make a lot of sense. Does such a result for edf and rdf indicate a
problem (e.g., collinearity) or any particular interpretation?
The other two questions concern which model fits the data best. We do look
at the usual various fit statistics (R^2, Dev, etc), but our question
concerns using the anova function to do model comparisons, e.g.,
anova(gama2.1,gama2.2, test=Chisq).
1. Is there research on the power of the model comparison test? Anecdotally,
the test seems to reject the null even in cases that would appear to have
only small differences. These are not hugely long time series, ranging from
about 17 to about 49, so we would not have thought them to yield large
power.
2. More important, in a few cases, we are getting a result that looks like
this:
anova(gamb1.1,gamb1.2, test=Chisq)
Analysis of Deviance Table
Model 1: y ~ x + z + int
Model 2: y ~ x + z + s(int, bs = cr)
Resid. Df Resid. Dev Df Deviance P(|Chi|)
130 36.713
230 36.713 1.1469e-05 1.0301e-05 6.767e-05 ***
We are inclined to think that the significance p value here is simply a
result of rounding error in the computation of the df difference and
deviance difference, and that we should treat this as indicating the models
are not different from each other. Has anyone experienced this before? Is
our interpretation reasonable?
Thanks to anyone who is able to offer advice.
Will Shadish
--
View this message in context:
http://r.789695.n4.nabble.com/GAM-Chi-Square-Difference-Test-tp4636523.html
Sent from the R help mailing list archive at Nabble.com.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Simon Wood, Mathematical Science, University of Bath BA2 7AY UK
+44 (0)1225 386603 http://people.bath.ac.uk/sw283
__
