0

I have a data set like this example data, which includes three sites. I want to model the trend of "cumSR" along "timeInt". I did it for each site using the following model ("cumSR~b/(cc+(timeInt)^(-z))"). Now I want to use the same model to get only one curve for these three sites. I think that I need to use nonlinear mixed model, right? I have tried "nlme" function in "nlme" package and "nlmer" function in "lme4" function, but haven't figured out. Could anyone help on this? Thanks so much.

    dat = data.frame(cumSR=c(17, 25, 34, 43, 46, 49, 51, 52, 55, 57, 59, 60, 18, 26, 34, 38, 41, 43, 45, 49, 51, 54, 
                     57, 58, 59, 23, 29, 33, 37, 38, 40, 44, 46, 47),
             timeInt=c(2.5, 7.5, 12.5, 17.5, 22.5, 27.5, 32.5, 42.5, 47.5, 52.5, 62.5, 67.5, 2.5, 7.5, 12.5, 17.5, 
                       22.5, 27.5, 32.5, 37.5, 42.5, 47.5, 52.5, 57.5, 62.5, 2.5, 7.5, 12.5, 17.5, 22.5, 37.5, 42.5, 47.5, 57.5),
             site=c("A", "A", "A", "A", "A", "A", "A", "A", "A", "A", "A", "A", "B", "B", "B", "B", "B", "B","B", 
                    "B", "B", "B", "B", "B", "B", "C","C","C","C","C","C","C","C","C"))


    par(mfcol=c(1,3))
    dat1 <- dat[which(dat$site=="A"),]
    SAR.nls.logist <- nls(cumSR~b/(cc+(timeInt)^(-z)), data=dat1,start = list(b=10,cc=.1,z=.5))
    plot(dat1$timeInt,dat1$cumSR, type="b",col=1)
    curve(coef(SAR.nls.logist)[1]/(coef(SAR.nls.logist)[2]+(x)^(-coef(SAR.nls.logist)[3])), col = "pink", type="l",add = TRUE,lwd=2)

    dat1 <- dat[which(dat$site=="B"),]
    SAR.nls.logist <- nls(cumSR~b/(cc+(timeInt)^(-z)), data=dat1,start = list(b=10,cc=.1,z=.5))
    plot(dat1$timeInt,dat1$cumSR, type="b",col=1)
    curve(coef(SAR.nls.logist)[1]/(coef(SAR.nls.logist)[2]+(x)^(-coef(SAR.nls.logist)[3])), col = "pink", type="l",add = TRUE,lwd=2)

    dat1 <- dat[which(dat$site=="C"),]
    SAR.nls.logist <- nls(cumSR~b/(cc+(timeInt)^(-z)), data=dat1,start = list(b=10,cc=.1,z=.5))
    plot(dat1$timeInt,dat1$cumSR, type="b",col=1)
    curve(coef(SAR.nls.logist)[1]/(coef(SAR.nls.logist)[2]+(x)^(-coef(SAR.nls.logist)[3])), col = "pink", type="l",add = TRUE,lwd=2)
Jian Zhang
  • 1,173
  • 5
  • 15
  • 20
  • 2
    Can you add the `nlme` and `nlmer` code that you tried, along with a brief explanation of how exactly it didn't work? – joran Nov 26 '13 at 02:55
  • 2
    with only three sites it probably won't be feasible to use mixed models (you will be trying to estimate a variance from three points). I would suggest just fitting the pooled model and comparing with the separate models (consider the `?nlsList` function too) – Ben Bolker Nov 26 '13 at 03:31
  • Thanks for the reply. My real data have more than 3 sites. The data I listed is just an example. – Jian Zhang Nov 26 '13 at 18:53
  • Some codes I tried to use "nlme" and "nlmer". I don't think that I use them correctly. library(nlme) nm1 <- nlme(cumSR~b/(cc+(timeInt)^(-z)),data = dat, fixed=cumSR + timeInt ~1, random = site ~ 1, start = c(b=10,cc=.1,z=.5)) library(lme4) nm2 <- nlmer(cumSR~b/(cc+(timeInt)^(-z)) ~ site|1, data=dat, start = c(b=10,cc=.1,z=.5)) – Jian Zhang Dec 16 '13 at 05:45
  • Have a look at this example: http://rpackages.ianhowson.com/cran/agridat/man/vold.longterm.html Maybe you want something like this (but it didn't converge): `m1 <- nlme(cumSR ~b/(cc+(timeInt)^(-z)), groups=~site, data = dat, fixed=list(b~1, cc~1, z~1), random = b+cc+z~1, start = c(b=18, cc= -.02, z=.2), control=list(maxIter=100))` – Kevin Wright Jun 11 '15 at 14:46

0 Answers0