nls уникальная градиентная матрица - параметры посадки в верхних пределах интеграла

Question

Я пытаюсь сделать nls подходящим для немного сложного выражения, которое включает в себя два интеграла с двумя параметрами соответствия в их верхних пределах.

я получил ошибку

"Ошибка в nlsModel (формула, М.Ф., начало, WTS): сингулярное градиент матрицы на начальных оценок параметров".

Я уже искал в предыдущих ответах, но не помог. Инициализация параметров выглядит нормально, я попытался изменить параметры, но никто не работает. Если у моей функции есть только один интеграл, все работает очень хорошо, но при добавлении второго интегрального члена просто получена ошибка. Я не считаю, что функция переопределена, поскольку я выполнял другие приемы с гораздо большим количеством параметров, и они работали. Ниже я написал список с некоторыми данными.

Минимально Примером может служить следующее:

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) 

Я не могу понять, что произойдет. Мне нужно выполнить эту подгонку для трех интегральных компонентов, но даже для двух я имею эту проблему. Я очень ценю вашу помощь. Спасибо.

Answer 1

Вы можете попробовать некоторые другие оптимизаторы:

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 

Левенберга-Марквардт сходится тоже, но 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