nls奇异梯度矩阵 – 拟合积分上限的参数

2023年5月2日 952次阅读

我试图使nls适合一个有点复杂的表达式,包括两个积分,其中两个拟合参数在其上限.

我收到了错误

“Error in nlsModel(formula, mf, start, wts) : singular gradient
matrix at initial parameter estimates”.

我已经在之前的答案中搜索过,但没有帮助.参数初始化似乎没问题,我试图改变参数但没有工作.如果我的函数只有一个积分,一切都很好,但是当添加第二个积分项时,就会得到错误.我不相信该功能是过度参数化的,因为我已经执行了更多参数的其他拟合并且它们起作用.下面我写了一些包含一些数据的列表.

最小的例子如下：

integrand <- function(X) {
  return(X^4/(2*sinh(X/2))^2)
}

fitting = function(T1, T2, N, D, x){
  int1 = integrate(integrand, lower=0, upper = T1)$value
  int2 = integrate(integrand, lower=0, upper = T2)$value
  return(N*(D/x)^2*(exp(D/x)/(1+exp(D/x))^2
)+(448.956*(x/T1)^3*int1)+(299.304*(x/T2)^3*int2))
}
fit = nls(y ~ fitting(T1, T2, N, D, x),
start=list(T1=400,T2=200,N=0.01,D=2))

——>作为参考,适用的适合如下：

integrand <- function(X) {
  return(X^4/(2*sinh(X/2))^2)
}
fitting = function(T1, N, D, x){
  int = integrate(integrand, lower=0, upper = T1)$value
  return(N*(D/x)^2*(exp(D/x)/(1+exp(D/x))^2 )+(748.26)*(x/T1)^3*int)
}
fit = nls(y ~ fitting(T1 , N, D, x), start=list(T1=400,N=0.01,D=2))

——->数据说明问题：

dat<- read.table(text="x       y
0.38813 0.0198
0.79465 0.02206
1.40744 0.01676
1.81532 0.01538
2.23105 0.01513
2.64864 0.01547
3.05933 0.01706
3.47302 0.01852
3.88791 0.02074
4.26301 0.0256
4.67607 0.03028
5.08172 0.03507
5.48327 0.04283
5.88947 0.05017
6.2988  0.05953
6.7022  0.07185
7.10933 0.08598
7.51924 0.0998
7.92674 0.12022
8.3354  0.1423
8.7384  0.16382
9.14656 0.19114
9.55062 0.22218
9.95591 0.25542", header=TRUE)

我无法弄清楚发生了什么.我需要对三个完整的组件执行此操作,但即使是两个我也有这个问题.我非常感谢你的帮助.谢谢.

最佳答案 您可以尝试其他一些优化器：

fitting1 <- function(par, x, y) {
  sum((fitting(par[1], par[2], par[3], par[4], x) - y)^2)
}

library(optimx)
res <-  optimx(c(400, 200, 0.01, 2),
       fitting1,
       x = DF$x, y = DF$y,
       control = list(all.methods = TRUE))

print(res)

#                  p1       p2          p3         p4         value fevals gevals niter convcode kkt1 kkt2 xtimes
#BFGS        409.7992 288.6416  -0.7594461   39.00871  1.947484e-03    101    100    NA        1   NA   NA   0.22
#CG          401.1281 210.9087  -0.9026459   20.80900  3.892929e-01    215    101    NA        1   NA   NA   0.25
#Nelder-Mead 414.6402 446.5080  -1.1298606 -227.81280  2.064842e-03     89     NA    NA        0   NA   NA   0.02
#L-BFGS-B    412.4477 333.1338  -0.3650530   37.74779  1.581643e-03     34     34    NA        0   NA   NA   0.06
#nlm         411.8639 333.4776  -0.3652356   37.74855  1.581644e-03     NA     NA    45        0   NA   NA   0.04
#nlminb      411.9678 333.4449  -0.3650271   37.74753  1.581643e-03     50    268    48        0   NA   NA   0.07
#spg         422.0394 300.5336  -0.5776862   38.48655  1.693119e-03   1197     NA   619        0   NA   NA   1.06
#ucminf      412.7390 332.9228  -0.3652029   37.74829  1.581644e-03     45     45    NA        0   NA   NA   0.05
#Rcgmin            NA       NA          NA         NA 8.988466e+307     NA     NA    NA     9999   NA   NA   0.00
#Rvmmin            NA       NA          NA         NA 8.988466e+307     NA     NA    NA     9999   NA   NA   0.00
#newuoa      396.3071 345.1165  -0.3650286   37.74754  1.581643e-03   3877     NA    NA        0   NA   NA   1.02
#bobyqa      410.0392 334.7074  -0.3650289   37.74753  1.581643e-03   7866     NA    NA        0   NA   NA   2.07
#nmkb        569.0139 346.0856 282.6526588 -335.32320  2.064859e-03     75     NA    NA        0   NA   NA   0.01
#hjkb        400.0000 200.0000   0.0100000    2.00000  3.200269e+00      1     NA     0     9999   NA   NA   0.01

Levenberg-Marquardt也会收敛,但是当nlsLM尝试从结果创建一个nls模型对象时会失败,因为渐变矩阵是单数的：

library(minpack.lm)
fit <- nlsLM(y ~ fitting(T1, T2, N, D, x),
          start=list(T1=412,T2=333,N=-0.36,D=38), data = DF, trace = TRUE)
#It.    0, RSS = 0.00165827, Par. =        412        333      -0.36         38
#It.    1, RSS = 0.00158186, Par. =    417.352    329.978    -0.3652     37.746
#It.    2, RSS = 0.00158164, Par. =    416.397    330.694  -0.365025    37.7475
#It.    3, RSS = 0.00158164, Par. =    416.618    330.568  -0.365027    37.7475
#It.    4, RSS = 0.00158164, Par. =    416.618    330.568  -0.365027    37.7475
#Error in nlsModel(formula, mf, start, wts) : 
#  singular gradient matrix at initial parameter estimates
点赞