[R-sig-ME] nlmer model specification

James Widman jwidman at mi.nmfs.gov
Wed Dec 22 20:04:05 CET 2010


Hi Ben,
   Thanks for your help, I'm closer but not there yet - see below.
On 12/21/2010 04:16 PM, Ben Bolker wrote:
> -----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:

   I'm assuming the equations were supposed to be Asym+a+b. I still 
don't see a fixed affect for foodtype.
>   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
>  nlmer(Weight..g. ~ SSgompertz(week, Asym, a, b) ~

+ Asym+b+a|foodtype/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) ~ Asym + b + a | foodtype/tankf/FishIDf

    Data: scupgrowth

   AIC  BIC logLik deviance

  4251 4383  -2104     4207

Random effects:

  Groups                   Name Variance   Std.Dev.   Corr

  FishIDf:(tankf:foodtype) Asym 3.1676e+01 5.62811999

                           b    3.4431e-04 0.01855561  1.000

                           a    2.7675e-02 0.16635946 -1.000 -1.000

  tankf:foodtype           Asym 3.9691e-07 0.00063001

                           b    5.4414e-07 0.00073766 -0.994

                           a    2.1861e-05 0.00467557 -0.996  1.000

  foodtype                 Asym 6.2047e-01 0.78769773

                           b    1.9456e-05 0.00441087  1.000

                           a    9.6219e-04 0.03101925 -1.000 -1.000

  Residual                      5.9801e-01 0.77331039

Number of obs: 2925, groups: FishIDf:(tankf:foodtype), 2925; tankf:foodtype, 8; foodtype, 2

Fixed effects:

       Estimate Std. Error t value

Asym 85.091422   1.569495   54.22

a     5.302438   0.031940  166.01

b     0.856355   0.003421  250.30

Correlation of Fixed Effects:

   Asym   a

a -0.131

b  0.647 -0.717

Warning message:

In mer_finalize(ans) : iteration limit reached without convergence (9)


>
>>> 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-----




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