非线性混合模型

Question

作为R 初学者，我正在尝试修复一个包含随机因素的模型。公式为：

Temp ~ a - (b * exp(-c *rate))

其中Temp 是温度，rate 是变化的量度（时间/温度）。为了构建模型并获取初始参数，我使用 nlme 包：

data1<-groupedData((Temp~rate|Year), data=data)
fm1<-nlme(Temp ~ a - (b * exp(-c *rate)), data=data, fixed=Temp~rate, start=c(a=8.10,b=7.24,c=0.5))
    Error in eval(expr, envir, enclos) : object 'a' not found`

我也试过这个：

`fm100<-selfStart(~a-(b*exp(-c*rate)),
function(mCall, data, LHS)
{
xy<-sortedXyData(mCall[["x"]], LHS, data)
tmp<-coef(lm(Temp~rate, data=data),
value<-c(exp(tmp[1],temp[2])
getInitial=c("a","b","c"))
}`

Error: unexpected symbol in:
"value<-c(exp(tmp[1],temp[2]) getInitial"

也许是一个简单的问题，但我还没有发现任何有用的东西。

以下是数据：

Temp<-c(9,11,11,11,11,12,12,12,12,12,12,12,12,12,12,12,12,13,13,13,13,13,13,13,14,14,17.1,17.6,17.6,19.6,20.6,21.3,21.3,22.8,22.8,22.8,22.8,22.8,22.8,22.8,22.8,22.8,22.8,22.8,23.5,23.8,11.4,11.4,11.4,11.4,11.7,11.7,12.6,13.6,13.6,14.6,14.6,14.6,14.6,15.7,15.7,15.7,16.1,16.1,16.7,16.7,11.6,12.6,12.6,12.6,14.5,14.5,14.7,15.8,15.8,15.8,15.8,16,16,16,16,16,16,16,16,16,16)
Rate<-c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1)
Year<-c(2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011) 
data<-data.frame(Temp, Rate, Year)

提前谢谢你！

Answer 1

首先，Rate 在rate 之间有一些拼写错误。

其次，我认为模型在当前形式下被过度参数化了，或许可以尝试以下类似的方法？

library(nlme)
Temp <- c(9,11,11,11,11,12,12,12,12,12,12,12,12,12,12,12,12,13,13,13,13,13,13,13,14,14,17.1,17.6,17.6,19.6,20.6,21.3,21.3,22.8,22.8,22.8,22.8,22.8,22.8,22.8,22.8,22.8,22.8,22.8,23.5,23.8,11.4,11.4,11.4,11.4,11.7,11.7,12.6,13.6,13.6,14.6,14.6,14.6,14.6,15.7,15.7,15.7,16.1,16.1,16.7,16.7,11.6,12.6,12.6,12.6,14.5,14.5,14.7,15.8,15.8,15.8,15.8,16,16,16,16,16,16,16,16,16,16)
Rate <- c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,0.0417,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1)
Year <- c(2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2006,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2010,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011,2011) 
data <- data.frame(Temp, Rate, Year)

fm1 <- nlme(Temp ~ a * exp(-b * Rate), 
    data = data, 
    fixed = a + b ~ 1, 
    random = a + b ~ 1 | Year, 
    start = c(a = 20.0, b = 0.5))

第三，我更改了起始值以更好地拟合数据，看看Temp 与Rate 的图：

f <- function(Rate) {
    a <- 20.0
    b <- 0.5
    a * exp(-b * Rate)
}

plot(Temp ~ Rate, data = data)
curve(f, add = TRUE, lwd = 2, col = "red")

最后，我猜测了你的随机效应的结构，我认为?nlme 中的fixed 和random 参数是一个很好的参考，所以你可以自己看看这些是如何工作和选择的正确的结构。

希望这有用！