I conducted an experiment where earthworms were subjected to two treatments,
with and without herbicide in the soil. Biomass measurements were taken every
12 days for 398 days and the biomass growth curves as a function of time were
plotted.
There was clearly a non-linear growth pattern such that an additive mixed
effects model was proposed to model the behavior of biomass as a function of
time and treatments.
When plotting the residuals a clear cone-shaped pattern was observed, therefore
a series of additive models were proposed sequentially to deal with violations
of the assumption of homogeneity. Below we can see the models with the
following names: M.1; M.2; M.3; M.4
lmc <- lmeControl (niterEM = 5000, msMaxIter = 1000)
f1 <- formula (Biomass ~ Treat + s (Time, by = Treat))
M.1 <-gamm (f1, random = list (fcajita = ~ 1), method = "REML", control = lmc,
data = Acorticis)
#This first model uses the experimental box factor (i.e. fcajita) as the random
element of the model. This random effects model assumes homogeneity between the
experimental boxes and within them over time
M.2 <-gamm (f1, random = list (fcajita = ~ 1), method = "REML", control = lmc,
data = Acorticis, weights = varIdent (form = ~ 1 | fcajita))
#This second model assumes heterogeneity between boxes, but homogeneity within
each box over time
M.3 <- gamm (f1, random = list (fcajita = ~ 1), method = "REML", control = lmc,
data = Acorticis, weights = varExp (form = ~ Time10))
#The third model assumes homogeneity between boxes but heterogeneity within
each box over time
Finally, we decided to model the heterogeneity using the 'varComb' function in
order to combine the variances where the model allows heterogeneity between the
experimental boxes and heterogeneity within the experimental boxes over time:
M.4 <- gamm (f1, random = list (fcajita = ~ 1), data = Acorticis, method =
"REML", control = lmc, weights = varComb (varIdent (form = ~ 1 | fcajita),
varPower (form = ~ Time10)))
The first three models executed perfectly and the following values of the AIC
indicator were obtained:
> AIC (M.1 $ lme, M.2 $ lme, M.3 $ lme)
df AIC
M.1 8 379.6464
M.2 15 309.5736
M.3 9 310.4828
Unfortunately, the execution of the M.4 model failed and the following error
message was obtained:
Error in environment (attr (ret $ lme $ modelStruct $ varStruct, "formula"))
<-. GlobalEnv:
attempt to set an attribute on NULL
A final model I tried was M5:
M.5 <- gamm(f1, random = list(fcajita =~ 1), data = Acorticis, method = "REML",
control = lmc, weights = varComb(varIdent(form = ~1|fcajita), varExp(form =~
Time10|fcajita)))
and this time I got the following error message:
Error in lme.formula(y ~ X - 1, random = rand, data = strip.offset(mf), :
nlminb problem, convergence error code = 1
message = function evaluation limit reached without convergence (9)
Además: Warning message:
In logLik.reStruct(object, conLin) :
Singular precision matrix in level -1, block 1
My question is: Could someone help me fix these problems to run the M.4 and M.5
models?
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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.
