Hi Raeanne,
gamm fits using PQL, which doesn't always converge, however looking at
your data, a couple of things stand out...
1. Trying a Tweedie seemed to give nicer residual plots than the
negative binomial, and also converges with gamm. e.g.
SB.gam3<-gam(count~offset(vol_offset)+
s(Depth_m, by=StnF, bs="cs")+StageF*RegionF+s(StnF,bs="re"),
family=Tweedie(p=1.4),
data=Zoop_2011[Zoop_2011$SpeciesF=='SB',],method="ML")
SB.gam3
plot(fitted(SB.gam3)^.5,residuals(SB.gam3))
qq.gam(SB.gam3)
... I just searched over p=1.2-1.6 in steps of .1 to find the optimum
value for the Tweedie p parameter according to the ML score. see ?Tweedie
2. If you only need a random effect for StnF, and don't have to have the
spatial autocorrelation structure then you can add the random effect in
gam, as above - this will be much more efficient than using gamm
(doesn't use PQL). See ?random.effects.
best,
Simon
On 07/12/12 12:00, Raeanne Miller wrote:
Hi there,
My question is about the 'theta' parameter in specification of a NB GAMM.
I have fit a GAM with an optimum structure of:
SB.gam4<-gam(count~offset(vol_offset)+
s(Depth_m, by=StnF, bs="cs")+StageF*RegionF,
family=negbin(1, link=log),
data=Zoop_2011[Zoop_2011$SpeciesF=='SB',])
However, this GAM shows heterogeneity in the residuals, and the model fit is
not great. I would like to add a random effect to this GAM, for station (StnF),
and a spatial correlation structure, for Depth_m. When trying to fit an
equivalent GAMM to start with, I used the formula:
SB.gamm1<-gamm(count~offset(vol_offset)+
s(Depth_m,by=StnF,bs="cs")+StageF*RegionF,
family=negative.binomial(theta=1.41), niterPQL=50,
data=Zoop_2011[Zoop_2011$SpeciesF=='SB',])
The theta=1.41 came from SB.gam4, using the theta.ml() command. However, this
model does not converge, even with 5000 iterations. When I increase theta to 6,
the model converges, but I'm not sure why. I have also tried theta=c(1,10),
which converges, but gives me lots of warnings:
Warning messages:
1: In y + .Theta :
longer object length is not a multiple of shorter object length
Can anyone tell me why this might be?
The model also does not ever converge if I replace family=negative.binomial()
with family=negbin(). I've looked at the help files for GAM and GAMM, and I'm
not sure what the difference between these is?
I have attached code to create my dataset, for those interested.
With many thanks,
Raeanne
The Scottish Association for Marine Science (SAMS) is registered in Scotland as
a Company Limited by Guarantee (SC009292) and is a registered charity (9206).
SAMS has an actively trading wholly owned subsidiary company: SAMS Research
Services Ltd a Limited Company (SC224404). All Companies in the group are
registered in Scotland and share a registered office at Scottish Marine
Institute, Oban Argyll PA37 1QA. The content of this message may contain
personal views which are not the views of SAMS unless specifically stated.
Please note that all email traffic is monitored for purposes of security and
spam filtering. As such individual emails may be examined in more detail.
______________________________________________
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
______________________________________________
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.