[R-sig-ME] gamm: problems with corCAR1

Karel Viaene Karel.Viaene at UGent.be
Thu Oct 6 15:52:31 CEST 2011


Dear all,

I?m analyzing this dataset containing time measures (Week),  
contaminant concentrations (Treatment) and diversity indices (richness  
for example). The diversity indices were calculated from abundance  
data of phytoplankton in a mesocosmos experiment. The abundances were  
determined weekly for two replicates (Replicate) per treatment.
I?m looking for the effects of time (Week) and contaminant levels  
(Treatment) on diversity indices (e.g. richness)
Initial analysis with GAM models showed temporal autocorrelation. So  
now I?m trying to fit this gamm (gamm1):

gamm1    <- gamm(richness~
s(Week,by=as.numeric(Treatment=="0"),k=6,fx=FALSE) +
s(Week,by=as.numeric(Treatment=="0.5"),k=6,fx=FALSE) +
s(Week,by=as.numeric(Treatment=="5"),k=6,fx=FALSE) +
s(Week,by=as.numeric(Treatment=="15"),k=6,fx=FALSE) +
s(Week,by=as.numeric(Treatment=="50"),k=6,fx=FALSE) +
s(Week,by=as.numeric(Treatment=="150"),k=6,fx=FALSE) +
s(Treatment,k=6,fx=FALSE) + factor(Treatment),
correlation=corCAR1(form=~Week|factor(Treatment),data=indices,family=gaussian)

I seem to be having difficulties with the correlation structure. An  
initial error occured because replicates were taken at the same time:

Error in Initialize.corCAR1(X[[2L]], ...) :
Covariate must have unique values within groups for corCAR1 objects

I solved this by selecting one replicate but is there another solution  
for this?
Moreover, when analyzing the data of one replicate, I received  
following error:

Error in MEestimate(lmeSt, grps) :
Singularity in backsolve at level 0, block 1

I have no idea how to solve this. It seems to be related with the  
complexity of the model because no error occurred when running a  
simpler gamm (gamm2):

gamm2 <- gamm(richness~
s(Week,k=6,fx=FALSE) + factor(Treatment),
correlation=corCAR1(form=~Week|conc.f), data=test,family=gaussian)

I?ve added my dataset. Any help would be well appreciated!

With kind regards,
Karel




More information about the R-sig-mixed-models mailing list