I have a data set which collected biomass of a species from 7 trees repeatedly. I assumed the biomass would change as a Gaussian peak function through time. Since it's non-linear with nested repeated measure and I would like to add tree group as a random effect, I assume a non-linear mixed effect model should do the job?
Here is the sample data: (can be downloaded from here https://www.dropbox.com/s/i5vm2fasx75shp0/tree.RData?dl=0 )
Here is the code:
>tree
time id Biomass
1 41 A1 7.6
2 41 A2 0.0
3 41 A3 71.1
4 57 A4 26.9
5 57 A5 52.1
6 57 A6 101.1
7 72 A1 0.0
8 72 A2 0.0
9 72 A3 34.0
10 83 A4 40.9
11 83 A5 195.4
12 83 A6 167.0
13 97 A1 17.6
14 97 A2 12.7
15 97 A3 12.4
16 111 A4 266.5
17 111 A5 139.6
18 111 A6 256.2
19 127 A1 111.4
20 127 A2 35.8
21 127 A3 72.9
22 149 A4 159.5
23 149 A5 305.5
24 149 A6 366.4
25 159 A1 19.5
26 159 A2 57.5
27 161 A3 205.6
28 174 A4 257.3
29 174 A5 166.2
30 175 A6 374.3
31 187 A1 159.3
32 187 A2 54.6
33 187 A3 136.4
34 204 A4 256.1
35 205 A5 423.3
36 204 A6 237.8
37 216 A1 67.4
38 216 A2 205.6
39 216 A3 316.1
40 232 A4 128.5
41 233 A5 20.8
42 233 A6 67.8
fm1<-Biomass~time+rd
#run model
model<-nlme(fm1, data=tree,
fixed=time~1,
random=rd~1,
groups=~id,
start=c(time=43))
#model summary
summary(model)
Is my formula correct? It seemed overly simple that I've been scratching my head to find where it indicated the Gaussian function? Since package nlme was designed to "Fit and compare Gaussian linear and nonlinear mixed-effects models" (as indicated in the package manual), I think variables were assumed to be based on the Gaussian distribution? Plus it did provide sigma estimate:
> model$sigma
[1] 104.1907
Sorry if it looks like a dumb question. This is my first time trying to fit such model. Thanks
