Bootstrapped confidence intervals for non-linear mixed effects model

18 Views Asked by At

I'm trying to get bootstrapped confidence intervals for parameter estimates of a non-linear (Gompertz) mixed effects model built in the lme4 package with function nlmer(). I need to be able to graph confidence intervals along the curve, I think with geom_ribbon(). However, I am running into problems and haven't gotten anything to work yet.

I first tried to use the confint() function to get bootstrapped confidence intervals, but that gives me this error: Error: step factor reduced below 0.001 without reducing pwrss

The model runs just fine and seems to converge well. This is my first time using non-linear mixed effects models, which I need to use to deal with repeated sampling of individuals within groups.

Here is some example data:

# create example data with similar structure 

x_data <- seq(1, 14, length.out = 200)
y_data <- 53 * exp(-exp(-0.14 *(x_data - 13))) + rnorm(length(x_data), mean = 0, sd = 0.5)
individ_ID <- c(1:200) # individual ID
group_ID <- c(1:10) # group ID 
group_ID2 <- rep(group_ID, each = 4)
mydata1 <- data.frame(x_data, y_data, individ_ID, group_ID2)

x_data <- seq(15, 28, length.out = 200)
y_data <- 53 * exp(-exp(-0.14 *(x_data - 13))) + rnorm(length(x_data), mean = 0, sd = 0.5)
individ_ID <- c(1:200)
group_ID <- c(1:10)
group_ID2 <- rep(group_ID, each = 4)
mydata2 <- data.frame(x_data, y_data, individ_ID, group_ID2)

mydata3 <- rbind(mydata1, mydata2)

The model is a Gompertz model: Aexp(-exp(-k(d - Ti)))

I set it up like this:


# create function for Gompertz model
gompertz_pp <- ~53*exp(-exp(-k*(d - Ti))) # fixed asymptote

fn.gompertz.pp <- deriv(gompertz_pp,
                        namevec=c("k","Ti"), 
                        function.arg=c("d","k","Ti"))


# create model

startvec <- c(k = 0.15, Ti = 12)

mod <- nlmer(y_data ~ fn.gompertz.pp(x_data, k, Ti) ~ (Ti|group_ID2/individ_ID) + (k|group_ID2/individ_ID), data = mydata3, start = startvec)

boot_result <- bootMer(mod, FUN = fixef, nsim = 100)

The model runs, and now I want need bootstrapped confidence intervals, or something that I can use to graph confidence intervals on a curve.

However, the confint() function failed:

confint(mod1)

and gave me this error message: Error: step factor reduced below 0.001 without reducing pwrss

I then tried to use the bootMer() function:

boot_result <- bootMer(mod1, FUN = fixef, nsim = 100)

Which also failed, with this error message: Error in model.matrix.default(eval(substitute(~foo, list(foo = x[[2]]))), : model frame and formula mismatch in model.matrix()

I am most interested in getting bootMer() to work.

Is there anywhere that I can improve my set up to get bootMer() or confint() to work? Or are there any other methods? I need confidence intervals that account for the random effect variation and which I can plot along a curve.

0

There are 0 best solutions below