我正在使用R中的nleqslv程序包来求解非线性方程组。
require(nleqslv)
x <- c(6,12,18,24,30)
NMfun1 <- function(k,n) {
y <- rep(NA, length(k))
y[1] <- -(5/k[1])+sum(x^k[2]*exp(k[3]*x))+2*sum(k[4]*x^k[2]*exp(-k[1]*x^k[2]*exp(k[3]*x)+k[3]*x)/(1-k[4]*exp(-k[1]*x^k[2]*exp(k[3]*x))))
y[2] <- -sum(log(x))-sum(1/(k[2]+k[3]*x))+sum(k[1]*x^k[2]*exp(k[3]*x)*log(x))+2*sum(k[1]*k[4]*exp(-k[1]*x^k[2]*exp(k[3]*x)+k[3]*x)*log(x)/(1-k[4]*exp(-k[1]*x^k[2]*exp(k[3]*x))))
y[3] <- -sum(x/(k[2]+k[3]*x))+sum(k[1]*x^(k[2]+1)*exp(k[3]*x))-sum(x)+2*sum(k[4]*x^k[2]*exp(-k[1]*x^k[2]*exp(k[3]*x)+k[3]*x)/(1-k[4]*exp(-k[1]*x^k[2]*exp(k[3]*x))))
y[4] <- -(5/(1-k[4]))+2*sum(exp(-k[1]*x^k[2]*exp(k[3]*x))/(1-k[4]*exp(-k[1]*x^k[2]*exp(k[3]*x))))
return(y)
}
kstart <- c(0.05, 0, 0.35, 0.9)
NMfun1(kstart)
nleqslv(kstart, NMfun1, control=list(btol=.0001),method="Newton")所获得的k的估计值为：0.04223362 -0.08360564 0.14216026 0.37854908但是k的估计值要大于零。
1条答案
按热度按时间ovfsdjhp1#
好的。所以你想要大于0的真实的解，如果它们存在的话。创建一个新函数，在将输入参数传递给
NMfun1之前对其进行平方。然后使用nleqslv包中的searchZeros函数来搜索解。像这样这段代码的最后4个非注解行的输出是
你会看到包含在对象
z中的解有一个负的元素，并且是平方的。从这个实验中可以看出，你的系统有一个正的解。