[R-sig-ME] nlmer model definition and starting values

[Sebastian P. Luque]

Hi,

I can't understand how to specify starting values for the following
nonlinear growth model.  The fake data produced below is designed to
closely match a real case scenario.  The model is for the following
Gompertz growth parameterization:

y=Linf * exp(-b * exp(-k * x))

where y is some size variable, and Linf, b, and k are parameters
indicating asymptotic size, size at birth, and a growth constant,
respectively.  I've hacked a self-start function for it.  Each row in
the fake data represents a single individual, with each individual
measured only once.  The data correspond to 4 species and 2 sexes, with
each combination having its own set of parameters, and some normal
heteroscedastic error added.

I would like to start with the assumption that b is a random effect,
while species and sex have fixed effects on Linf and k.


---<--------------------cut here---------------start------------------->---
"gompertz.fun" <- function(x, Linf, b, k) {Linf * exp(-b * exp(-k * x))}

## Set parameters for, say, 4 species and 2 sexes

grps <- expand.grid(species=LETTERS[1:4], sex=c("M", "F"))
parms <- c("Linf", "b", "k")
parms.mat <- matrix(c(432.3, 411.5, 389.7, 448.2,
                      381.5, 360.5, 338.1, 385.7,
                      1.16, 0.59, 0.73, 0.61,
                      7.56, 0.49, 0.51, 0.51,
                      0.15, 0.08, 0.07, 0.09,
                      0.29, 0.14, 0.11, 0.14),
                    ncol=3, dimnames=list(paste(grps[, 1], grps[, 2]), parms))

## Generate random x for each combination of species and sex, and y
## based on the parameters, with different y normal random noise
## (heteroscedastic) per combination
set.seed(123)
x <- rep(seq(1, 60, 3), 8) + runif(160, min=0, max=2)
fake <- data.frame(grps[rep(seq(nrow(grps)), each=20), ], x)
fake <- with(fake, fake[order(species, sex, x), ])
y <- apply(fake, 1, function(k) {
    grp <- paste(k[1], k[2])
    parms <- parms.mat[rownames(parms.mat) %in% grp, ]
    gompertz.fun(as.numeric(k[3]), parms[1], parms[2], parms[3])
})
## Add the increasing error
fake$y <- rnorm(160, mean=y, sd=seq(5, 30, length.out=20))
## Quick visualization
xyplot(y ~ x | species, data=fake, groups=sex, aspect=0.7,
       scales=list(alternating=1, tck=c(0.75, 0), rot=c(0, 90)))
---<--------------------cut here---------------end--------------------->---


How could I write the nlmer call for such a model, and how could I
anticipate the order of starting values to supply.  Thanks in advance
for any feedback on this.



Cheers,

-- 
Seb