[R] NLME question...

[Jeffrey Harring]

I would like to fit a nonlinear model with "nlme."

The model is a nonlinear growth model with five time points: y = X*b + 
e, where design matrix X is defined as

X= |  1   0  |
      |  1   1  |
      |  1   k  |
      |  1  2k |
      |  1  3k |

and parameter vector b = (b0, b1). And where "k" is a parameter to be 
estimated. Of course I also want to estimate the intercept (b0), slope 
(b1). Error variances (e) with 5 free parameters and random effects 
covariance matrix  (2x2: for b0 and b1).

I am not certain how to get the design matrix in the code so that the 
"nlme" procedure will be able to
estimate parameter "k" as well as the intercept and slope.

I have tried the following:

meanfunc <- function(x,b0,b1,k){
meangrad <- array(0,c(length(x),2),list(NULL,c("b0","b1")))

for(i in 1:length(x)){
    if(x[i]==0){
     err[i] <- b0 + b1*x[i]
       meangrad[i,"b0"] <- 1
       meangrad[i,"b1"] <- x[i]}
    if(x[i]==1){
     err[i] <- b0 + b1*x[i]
       meangrad[i,"b0"] <- 1
       meangrad[i,"b1"] <- x[i]}
    if(x[i]==2){
     err[i] <- b0 + b1*(x[i]-1)*k
     meangrad[i,"b0"] <- 1
       meangrad[i,"b1"] <- (x[i]-1)*k}
    if(x[i]==3){
       err[i] <- b0 + b1*2*(x[i]-1)*k
     meangrad[i,"b0"] <- 1
       meangrad[i,"b1"] <- (x[i]-1)*k}
    if(x[i]==4){
       err[i] <- b0 + b1*3*(x[i]-1)*k
     meangrad[i,"b0"] <- 1
       meangrad[i,"b1"] <- (x[i]-1)*k}
}
#  compute analytical derivatives
   attr(err,"gradient") <- meangrad
   err
}

data.mlfit <- nlme(y ~ meanfunc(x,b0,b1,k),
   fixed=list(b0 ~ 1, b1 ~ 1, k ~ 1),
   random=list(b0 ~ 1, b1 ~ 1),
   groups = ~subj,
   data=nd,
   start=list(fixed=c(300,50,1)),
   method="ML",verbose=T)

I get the following error message:

Error in meangrad[i, "b1"] <- (x[i] - 1) * k :
  number of items to replace is not a multiple of replacement length
In addition: Warning messages:
1: In err[i] <- b0 + b1 * x[i] :
  number of items to replace is not a multiple of replacement length
2: In err[i] <- b0 + b1 * x[i] :
  number of items to replace is not a multiple of replacement length
3: In err[i] <- b0 + b1 * (x[i] - 1) * k :
  number of items to replace is not a multiple of replacement length

Any suggestions or insights would be greatly appreciated.

Thank you,
Jeff

-- 
**********************************************************
Jeffrey R. Harring, Assistant Professor
Department of Measurement, Statistics & Evaluation (EDMS)
1230 Benjamin Building
University of Maryland
College Park, MD 20742-1115

Phone: 	301.405.3630
Fax: 	301.314.9245
Email: 	harring at umd.edu
Web:  	http://www.education.umd.edu/EDMS/fac/Harring/webpage.html