【发布时间】:2020-02-26 06:24:01
【问题描述】:
我正在使用mle2 来估计非线性模型的参数,并且我想要估计参数估计周围的误差(标准误差)。同样,我想使用该模型然后使用 newdata 进行预测,并且在此过程中的几个步骤中我遇到了问题(错误)。
这是数据:
table<- "ac.temp performance
1 2.17 47.357923
4 2.17 234.255317
7 2.17 138.002633
10 2.17 227.545902
13 2.17 28.118072
16 9.95 175.638448
19 9.95 167.392218
22 9.95 118.162747
25 9.95 102.770622
28 9.95 191.874867
31 16.07 206.897159
34 16.07 74.741619
37 16.07 127.219884
40 16.07 208.231559
42 16.07 89.544612
44 20.14 314.946107
47 20.14 290.994063
50 20.14 243.322497
53 20.14 192.335133
56 20.14 133.841776
58 23.83 139.746673
61 23.83 224.135993
64 23.83 126.726493
67 23.83 246.443386
70 23.83 163.019896
83 28.04 4.614154
84 28.04 2.851866
85 28.04 2.935584
86 28.04 153.868415
87 28.04 103.884295
88 30.60 0.000000
89 29.60 0.000000
90 30.30 0.000000
91 29.90 0.000000
92 30.80 0.000000
93 28.90 0.000000
94 30.00 0.000000
95 30.20 0.000000
96 30.40 0.000000
97 30.70 0.000000
98 27.90 0.000000
99 28.60 0.000000
100 28.60 0.000000
101 30.40 0.000000
102 30.60 0.000000
103 29.70 0.000000
104 30.50 0.000000
105 29.30 0.000000
106 28.90 0.000000
107 30.40 0.000000
108 30.20 0.000000
109 30.10 0.000000
110 29.50 0.000000
111 31.00 0.000000
112 30.60 0.000000
113 30.90 0.000000
114 31.10 0.000000"
perfdat<- read.table(text=table, header=T)
首先,我必须为我的非线性模型设置几个固定参数,以了解动物在温度方面的表现
pi = mean(subset(perfdat, ac.temp<5)$performance)
ti = min(perfdat$ac.temp)
定义x 变量(温度)
t = perfdat$ac.temp
为非线性模型创建函数
tpc = function(tm, topt, pmax) {
perfmu = pi+(pmax-pi)*(1+((topt-t)/(topt-tm)))*(((t-ti)/(topt-ti))^((tm-ti)/(topt-tm)))
perfmu[perfmu<0] = 0.00001
return(perfmu)
}
创建负对数似然函数
LL1 = function(tm, topt, pmax, performance=perfdat$performance) {
perfmu = tpc(tm=tm, topt=topt, pmax=pmax)
loglike = -sum(dnorm(x=performance, mean=perfmu, log=TRUE))
return(loglike)
}
使用mle2 - maximum likelihood estimator 的模型性能
m1<- mle2(minuslogl = LL1, start=list(tm=15, topt=20, pmax=150), control=list(maxit=5000))
summary(m1)
这给了我参数估计值,但不是warning message: In sqrt(diag(object@vcov)) : NaNs produced 的误差估计值(标准错误)。但是,参数估计很好,可以让我做出有意义的预测。
Plot of non-linear curve using parameter estimates
我尝试了许多不同的优化器和方法,并得到了关于无法计算标准的相同错误。错误,通常带有关于无法反转粗麻布的警告。或者我对我的参数的估计非常不稳定,而且没有任何意义。
如果我使用:
confint(m1)
我的每个参数都有 95% 的间隔,但我无法将它们合并到可以用来制作如下图的预测方法中,我使用 nls 模型和 predict() 制作了该图:
non-linear model with error graphed
如果我通过将模型公式嵌入mle2 函数来重新创建我的mle2() 模型:
tpcfun<- function(t, tm.f, topt.f, pmax.f) {
perfmu = pi+(pmax.f-pi)*(1+((topt.f-t)/(topt.f-tm)))*(((t-ti)/(topt.f-ti))^((tm.f-ti)/(topt.f-tm.f)))
perfmu[perfmu<0] = 0.00001
return(perfmu)
}
m2<- mle2(performance ~ dnorm(mean=-sum(log(tpcfun(t=ac.temp, tm.f, topt.f, pmax.f))), sd=1), method='L-BFGS-B', lower=c(tm.f=1, topt.f=1, pmax.f=1), control=list(maxit=5000, parscale=c(tm.f=1e1, topt.f=1e1, pmax.f=1e2)), start=list(tm.f=15, topt.f=20, pmax.f=150), data=perfdat)
summary(m2)
我的参数得到了不合理的估计值,但我仍然没有得到错误估计值。
我的问题是,任何人都可以看出我的模型(模型函数和似然函数)有什么问题,或者我做错了什么吗?我有一种感觉,我可能写错了似然函数,但我尝试了各种分布和不同的方式,但我可能完全搞砸了。
或者有没有一种方法可以让我对我的参数的误差进行估计,以便我可以使用它们在图表中可视化我的模型预测的误差?
谢谢,
赖利
PS。我决定制作一个只有点估计值的图表,然后是趋势线没有错误,但我想在每个点估计值上放置 95% CI 的条形图,但是 confint() 给了我非常小的 CI '甚至没有出现在图表上,因为它们比我使用的点字符小,哈。
【问题讨论】:
标签: r predict mle hessian non-linear