[R] Help to solve modeling problem with gamm

[Abby Spurdle]

Sorry, I need to apologize.

My statement that "In general, linear models work well for biological
growth curves with small to medium sized data", may not be correct.
Also, I read through you code too quickly, my bad...

If you want to set weights (not sure if that's a good idea or not
here), I'd recommend you set them to a ***numeric vector***.

Your code is trying to set the weights to more complex objects (from
the nlme package), and it's possible that the gamm function is trying
to coerce them to numeric vectors, and that may be causing problems.

For nontrivial weights, you may want to compute weights (ensuring that
they are a numeric vector) in a separate step, and then fit the model
in a second step.


On 2/28/20, José Antonio García Pérez <garci95 using hotmail.com> wrote:
> 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 using 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.
>