Skip to content

nlme starting values are not the correct length

3 messages · Lara Reichmann, PIKAL Petr

Lara Reichmann
# 
Dear R community,

I am trying to fit an nlme model where I want to estimate the fixed effects of two treatments on the parameters on the following equation Photo~(a*(1-exp(-c*PARi/a)))-b
I was able to fit a simple model without covariates following the method described in Mixed-Effects Methods and Classes for S and S-PLUS, version 3.0, but when I add the covariates, I get the error " starting values for the fixed component are not the correct length"

My data has the following structure "Subject" "Species" "Fert" "Photo" "PARi" , where several "Photo" measurements where taken on the same subject by changing "PARi", 4 Species levels and 2 Fert levels, there are 31 Subjects (one missing value), and 323 observations

DATA extract

Subject	Species	Fert	Photo	PARi
1	bb	1	22.5389	1499.3307
1	bb	1	21.881	1248.913
1	bb	1	21.2862	999.3387
1	bb	1	20.5836	799.9308
1	bb	1	19.3758	601.1412
1	bb	1	15.5915	399.815
1	bb	1	8.7978	200.1087
1	bb	1	4.4347	99.686
1	bb	1	2.0387	49.7842
1	bb	1	-1.4854	0.0576
2	sw	0	6.782	1500.5337
2	sw	0	7.1432	1249.2749
2	sw	0	7.3319	1000.9891
2	sw	0	7.5848	799.1752
2	sw	0	7.1882	599.5544
2	sw	0	6.809	399.988
2	sw	0	5.3877	198.7574
2	sw	0	3.5104	100.7115
2	sw	0	0.8856	50.7015
2	sw	0	-1.121	0.0569
3	jg	1	16.0827	2000.4941
3	jg	1	16.0236	1501.1957
3	jg	1	16.3818	1248.9551
3	jg	1	16.7815	1499.6414
3	jg	1	17.175	2000.6851
3	jg	1	16.6529	1000.2707
3	jg	1	15.7987	799.676
3	jg	1	15.5437	598.9409
3	jg	1	11.7683	400.7715
3	jg	1	4.89	200.7468
3	jg	1	4.1294	100.9664
3	jg	1	1.6008	50.9254
3	jg	1	-0.89	0.5347
4	sw	1	25.2889	2000.1454
4	sw	1	24.7284	1499.6191
4	sw	1	24.3637	1249.7523
4	sw	1	23.3523	1000.0944
4	sw	1	21.6057	800.2209
4	sw	1	18.8926	599.7022
4	sw	1	14.6598	398.9366
4	sw	1	7.7182	201.5697
4	sw	1	3.4775	100.5139
4	sw	1	1.169	49.7045
4	sw	1	-1.3558	1.6914
5	jg	0	6.1626	2000.9351
5	jg	0	7.5573	1499.6581
5	jg	0	7.7129	1249.5073
5	jg	0	7.442	1000.7276
5	jg	0	7.5135	799.1286
5	jg	0	7.1559	599.5568
5	jg	0	6.8161	400.3576
5	jg	0	4.0097	199.7442
5	jg	0	2.7202	101.1253
5	jg	0	1.0746	51.1787
5	jg	0	-0.5913	0.975



This works so far:

lightresponse<-groupedData(Photo~PARi|Subject,data=lightr,outer = ~ Species * Fert,labels = list(x = "PAR", y = "CO2 uptake rate"),units = list(x = "(photon s-1)", y = "(umol/m?2 s)"))  
Photo.resp<-function(PARi,A,B,C)A*(1-exp(-C*PARi/A))-B
Photo.resp<-deriv ((~A *(1-exp(-C*PARi/A))-B),c("A","B","C"),function(PARi,A,B,C){})
OUTPUT
Nonlinear mixed-effects model fit by maximum likelihood
  Model: Photo ~ Photo.resp(PARi, A, B, C) 
  Data: lightresponse 
  Log-likelihood: -494.5926
  Fixed: A + B + C ~ 1 
          A           B           C 
24.89334793  1.77983637  0.06499634 

Random effects:
 Formula: list(A ~ 1, B ~ 1, C ~ 1)
 Level: Subject2
 Structure: General positive-definite, Log-Cholesky parametrization
         StdDev      Corr       
A        10.67382785 A     B    
B         0.52572012 1.000      
C         0.01433605 0.371 0.384
Residual  0.71900020            

Number of Observations: 323
Number of Groups: 31 

##Now, I want to test the effect of Species and Fert, I don't fully understand how to modify the start vector, as I tried several options and no one seems to be correct. Do the number of levels in each factor matter? In that case 4 Species and 2 Fert levels, I would need 6 initial parameters x 3? This didn't work either
Error in nlme.formula(model = Photo ~ Photo.resp(PARi, A, B, C), fixed = A +  : 
  starting values for the fixed component are not the correct length

I hope someone out there has the answer!
Thanks!!!

Lara

Lara G. Reichmann
Postdoctoral Fellow
USDA-ARS
808 E Blackland Rd
Temple, TX 76502
1 day later
PIKAL Petr
# 
Hi

see inline
A aqm not an expert in nlme models but it seems that you have only three fixed parameters A,B,C and want to feed them with 12 starting values. This is probably the reason for error.

Maybe you shall also search in R-sig-mixed-models

Regards
Petr
Lara Reichmann
# 
Hi Petr,
Thanks for the email, the number of parameters is not 3 because I am trying to fit a parameter for each treatment level, therefore the start vector changes depending on the fixed statement (see below). 

I finally got some of these to work....but for some reason I cannot run a complete model with both fixed effects and interactions (Species*Fert), I still couldn't get the syntax for the start vector! 

Thanks again,
Lara

All these work:
<unfert <- subset(lightresponse, Fert==0)
<model.cov4<-nlme(Photo~A*(1-exp(-lC*PARi/A))-B,fixed=A+B+lC~Species,random=A+B+lC~1,data=unfert,start=c(24.89,0,0,0,1.78,0,0,0,0.06,0,0,0))
<plot(augPred(model.cov4,outer=T))
#Here is the summary of the non-linear mixed model:
<summary(model.cov4)

<fert <- subset(lightresponse, Fert==1)
<model.cov5<-nlme(Photo~A*(1-exp(-lC*PARi/A))-B,fixed=A+B+lC~Species,random=A+B+lC~1,data=fert,start=c(24.89,0,0,0,1.78,0,0,0,0.06,0,0,0))
<plot(augPred(model.cov5))
#Here is the summary of the non-linear mixed model:
<summary(model.cov5)

<model.cov6<-nlme(Photo~A*(1-exp(-lC*PARi/A))-B,fixed=A+B+lC~Fert,random=A+B+lC~1,data=lightresponse,start=c(24.89,0,1.78,0,0.06,0))
<plot(augPred(model.cov6))
#Here is the summary of the non-linear mixed model:
<summary(model.cov6)

This one I am interested in but cannot figure out the start vector!

<model.nlme<-nlme(Photo~A*(1-exp(-lC*PARi/A))-B,fixed=A+B+lC~Species*Fert,random=A+B+lC~1,data=lightresponse,start=c(?????)

Also, once I have the summary of effects, how do I do contrasts for specific species differences?
I understand from this output that the p value is comparing each species with the baseline species, but how do I compare jg with sw for example? Is there a simple "Tukey" or family contrast for nlme method?

Fixed effects: A + B + lC ~ Species 
                   Value Std.Error  DF   t-value p-value
A.(Intercept)  24.132181  3.649618 281  6.612249  0.0000
A.Speciesjg    -4.233378  5.160757 281 -0.820302  0.4127
A.Specieslb     2.862573  5.161393 281  0.554612  0.5796
A.Speciessw     4.964391  5.344752 281  0.928835  0.3538
B.(Intercept)   1.523276  0.262774 281  5.796899  0.0000
B.Speciesjg     0.115693  0.372653 281  0.310457  0.7564
B.Specieslb     0.469510  0.372018 281  1.262062  0.2080
B.Speciessw     0.476496  0.380826 281  1.251220  0.2119
lC.(Intercept)  0.065695  0.005050 281 13.009079  0.0000
lC.Speciesjg   -0.009383  0.007155 281 -1.311311  0.1908
lC.Specieslb    0.009004  0.007178 281  1.254371  0.2107
lC.Speciessw   -0.002268  0.007369 281 -0.307801  0.7585
On Dec 6, 2012, at 12:52 AM, PIKAL Petr <petr.pikal at precheza.cz> wrote: