[R-sig-ME] nlmer model specification

[Ben Bolker]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

On 10-12-21 02:43 PM, James Widman wrote:
> Happy Holidays to All,
> 
>  I'm trying to use nlmer to model some data from an experiment I ran and
> am having some trouble specifying the model correctly. Would appreciate
> any insight you can provide.
> 
> I am examining the growth of fish fed two diets over 9 weeks - with
> weekly measurements. Thirty fish were placed in each of 8 tanks, 4 tanks
> fed foodtype NRD and 4 tanks fed foodtype OTO.
> Unbalanced data due to mortality.
> I'm trying to fit a nlmer to the data and want to make sure I have
> specified the model correctly. I can provide additional information if
> needed.
> This is what I came up with and have tried a number of alternatives but
> don't seem to  be progressing.  I realize the model will probably be
> simplified.
> 
> foodtype should be a fixed effect.
> 

  Hmmm.  It's not obvious to me from looking at the ?lmer page how one
would go about fitting this sort of 'split-plot' design in lmer: I see
how one would allow all three variables to vary among tanks and among
fish within tanks:

 nlmer(Weight..g. ~ SSgompertz(week, Asym, a, b) ~
Asym+b+c|tankf/FishIDf, ...)

  If there weren't a problem with specifying a fixed effect with only
two types, we could do

nlmer(Weight..g. ~ SSgompertz(week, Asym, a, b) ~
Asym+b+c|foodtype/tankf/FishIDf, ...)

  but that probably won't work well.

  My next try, unless someone on the list has a better solution, would
probably be to try this in AD Model Builder.

  (I just got a copy of Madsen and Thyregod "Intro to general and
generalized linear models", which has an interesting section on
extending the orange-tree analysis (?Orange) using Laplace
approximations coded around nlme and using AD Model Builder ...)

  Ben Bolker


>> str(scupgrowth)
> 'data.frame':   2925 obs. of  11 variables:
>  $ StartDate  : Factor w/ 1 level "2008/07/29": 1 1 1 1 1 1 1 1 1 1 ...
>  $ SampleDate : Factor w/ 10 levels "2008/07/29","2008/08/05",..: 1 1 1
> 1 1 1 1 1 1 1 ...
>  $ FishID     : int  1 2 3 4 5 6 7 8 9 10 ...
>  $ Length..mm.: int  38 39 34 33 36 40 39 34 32 35 ...
>  $ Weight..g. : num  0.64 0.81 0.47 0.48 0.62 0.78 0.7 0.47 0.47 0.62 ...
>  $ Tank       : int  1 1 1 1 1 1 1 1 1 1 ...
>  $ Days       : int  0 0 0 0 0 0 0 0 0 0 ...
>  $ FishIDf    : Factor w/ 2925 levels "1","2","3","4",..: 1 2 3 4 5 6 7
> 8 9 10 ...
>  $ foodtype   : Factor w/ 2 levels "NRD","OTO": 1 1 1 1 1 1 1 1 1 1 ...
>  $ tankf      : Factor w/ 8 levels "1","2","3","4",..: 1 1 1 1 1 1 1 1 1
> 1 ...
>  $ week       : num  0 0 0 0 0 0 0 0 0 0 ...
>> (nlm1 <- nlmer(Weight..g. ~ SSgompertz(week, Asym, a, b) ~ foodtype +
> tankf + tankf|FishIDf,
> +  scupgrowth, start = c(Asym = 72.70484, a = 5.34048, b = 0.83792)))
> Nonlinear mixed model fit by the Laplace approximation
> Formula: Weight..g. ~ SSgompertz(week, Asym, a, b) ~ foodtype + tankf
> +      tankf | FishIDf
>    Data: scupgrowth
>   AIC  BIC logLik deviance
>  6010 6249  -2965     5930
> Random effects:
>  Groups   Name        Variance   Std.Dev.  Corr
>  FishIDf  foodtypeOTO 4.7489e-06 0.0021792
>           tankf2      2.8193e-04 0.0167909 -0.829
>           tankf3      1.4213e-04 0.0119219  0.001  0.110
>           tankf4      2.0968e-04 0.0144804 -0.913  0.759  0.056
>           tankf5      1.8160e-04 0.0134760 -0.023 -0.007 -0.035  0.019
>           tankf6      3.7620e-04 0.0193958 -0.800  0.704  0.014  0.731
> -0.260
>           tankf7      2.0907e-04 0.0144592  0.023 -0.195 -0.062 -0.022
> -0.038
>           tankf8      3.4983e-04 0.0187038 -0.958  0.701 -0.033  0.870 
> 0.222
>  Residual             1.1874e+00 1.0896851
> 
>  -0.021
>   0.702 -0.107
> 
> Number of obs: 2925, groups: FishIDf, 2925
> 
> Fixed effects:
>       Estimate Std. Error t value
> Asym 90.150991   5.205204   17.32
> a     5.160145   0.042640  121.02
> b     0.868886   0.003198  271.73
> 
> Correlation of Fixed Effects:
>   Asym  a
> a 0.748
> b 0.970 0.576
> 
> When I tried the following I get an error.
>   (nlm1 <- nlmer(Weight..g. ~ SSgompertz(week, Asym, a, b) ~ (1 |
> foodtype) + tankf + tankf|FishIDf,
> +  scupgrowth, start = c(Asym = 72.70484, a = 5.34048, b = 0.83792)))
> Error in model.matrix.default(eval(substitute(~expr, list(expr =
> x[[2]]))),  :
>   model frame and formula mismatch in model.matrix()
> 


-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.10 (GNU/Linux)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org/

iEYEARECAAYFAk0RGUIACgkQc5UpGjwzenPswwCfUo5zcC9OBCKRTwwoUuGBrr0k
8xEAmwa6mpJaboeDXTZrod8jqqDbBXj4
=VpL8
-----END PGP SIGNATURE-----