R语言 极大似然法求带积分函数最优参数 报错non-finite function value

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

各位大牛们好，我在R中尝试用极大似然法（maxLik）求一个带积分（integrate）函数的最优参数值，运行后发现报了“non-finite function value ”的错误。


逐步排查后，发现是下面的简化的代码里出了问题，用Rstudio自带的debug功能尝试找了下，也没有找到具体的问题。恳请各位看看是哪里有问题呢？多谢多谢！

具体代码如下：

1. integrand <− function(x,para) {exp(5+para*x)}
   
2. 
   
3. integrand2 <− function(para1) {integrate( integrand, lower = -2 , upper = 2, para=para1)$value }
   
4. 
   
5. MLE <- maxLik(logLik=integrand2,start=0.3)

复制代码

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

关键词：function finite value 极大似然法 极大似然

代码没骗你，exp(5+para*x)这个函数的定积分是
(e^(5+para*2)-e^(5-para*2))/para，把para看作自变量，这个函数的曲线类似一个开口向上的抛物线，没有最大值。

你把第一行的函数换一个试试，maxLik()是可以出结果的，说明maxLik(), integrate()这些函数都没有问题。

换句话说，你的函数找不到最优参数

llb_321 发表于 2021-3-4 20:09
换句话说，你的函数找不到最优参数

谢谢你的回复，确实是我简化的函数设定的问题，我把简化函数改了，但还是会有问题，可以麻烦你再看一下吗？我的目的是估计“para”的值，使U_Cur/（U_Cur+U_Tgt）最大，其中U_Cur和U_Tgt是exp(A)/(exp(A)+exp(B))的形式。这是我重新设定的函数：

1. library("maxLik")
   
2. 
   
3. integrand <− function(x,para)
   
4. {
   
5.   # print(para)
   
6.   
   
7.   # dv, dv_T, dn, dn_T are the actual input data, here only gives one set of value for easy debugging
   
8.   dv <- -1.4527280  
   
9.   dv_T <- 0.2038909
   
10.   dn <- 17.5174383
    
11.   dn_T <- 26.6993828
    
12.   
    
13.   beta10 = para[1];
    
14.   beta11 = para[2];
    
15.   beta12 = para[3];
    
16.   alpha1 = para[4];
    
17.   beta20 = para[5];
    
18.   beta21 = para[6];
    
19.   beta22 = para[7];
    
20.   alpha2 = para[8];
    
21.   
    
22.   U_Cur = exp(beta10 + beta11*dv + beta12*dn + alpha1*x);
    
23.   U_Tgt = exp(beta20 + beta21*dv_T + beta22*dn_T + alpha2*x);
    
24.   
    
25.   # print(U_Cur)
    
26.   # print(U_Tgt)
    
27.   
    
28.   P <- (U_Cur/(U_Cur+U_Tgt))
    
29.   
    
30.   # print(P)
    
31.   
    
32.   return(P)
    
33. }
    
34. 
    
35. 
    
36. integrand2 <− function(para1) {integrate( integrand, lower = -2 , upper = 2, para=para1)$value }
    
37. 
    
38. MLE<-maxLik(logLik=integrand2,start=c(0.5, 0.2, 0.1, 0.3, -0.5, 0.2, 0.1, 0.6))

复制代码

我通过调试，发现报错的原因是U_Cur或U_Tgt有时可以是inf或0。我把U_Cur/（U_Cur+U_Tgt）加个log再变换，或者变成把U_Cur=exp(A)/(exp(A)+exp(B))变成1-(exp(B)/(exp(A)+exp(B)))等方法都试了。还是会报错。请问还有别的什么办法吗？

我试了新的函数，也用你的初始参数做了函数图像和积分图像，然后把你的函数又作了简化，用一个参数时exp(x)/(exp(x)+exp(a*x))，似乎有结果。改成用两个参数时，报错。我也不知道为什么。抱歉，能力有限。
类似exp(x)/(exp(x)+exp(a*x))的函数，似乎只在某个区间内可以做定积分的，结果基本上是区间上限的值。
猜，可能还是参数对函数积分结果的影响无法确定。

llb_321 发表于 2021-3-5 20:32
我试了新的函数，也用你的初始参数做了函数图像和积分图像，然后把你的函数又作了简化，用一个参数时exp(x ...

嗯嗯，好的，还是谢谢你，麻烦了。我后来直接给exp(A)/(exp(A)+exp(B))限定了范围，如果出现(-)inf，我就赋值为(-)2^1000，不过运行结果还是会报错"The integral is probably divergent"。我目前打算直接用GA算法求参数值，但就是没办法得到参数的t检验结果了。。再次谢谢你，耽误你时间和精力了，祝好~