R nls奇异梯度

Question

我已经尝试在这个主题上搜索其他线程,但没有一个修复程序对我有用.我有一个自然实验的结果,我想显示一个符合指数分布的事件的连续出现次数.我的R shell贴在下面

f <- function(x,a,b) {a * exp(b * x)}
> x
 [1]  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
[26] 26 27
> y
 [1] 1880  813  376  161  100   61   31    9    8    2    7    4    3    2    0
[16]    1    0    0    0    0    0    1    0    0    0    0    1
> dat2
    x    y
1   1 1880
2   2  813
3   3  376
4   4  161
5   5  100
6   6   61
7   7   31
8   8    9
9   9    8
10 10    2
11 11    7
12 12    4
13 13    3
14 14    2
> fm <- nls(y ~ f(x,a,b), data = dat2, start = c(a=1, b=1)) 
Error in numericDeriv(form[[3L]], names(ind), env) : 
  Missing value or an infinity produced when evaluating the model
> fm <- nls(y ~ f(x,a,b), data = dat2, start = c(a=7, b=-.5)) 
Error in nls(y ~ f(x, a, b), data = dat2, start = c(a = 7, b = -0.5)) : 
  singular gradient
> fm <- nls(y ~ f(x,a,b), data = dat2, start = c(a=7,b=-.5),control=nls.control(maxiter=1000,warnOnly=TRUE,minFactor=1e-5,tol=1e-10),trace=TRUE) 
4355798 :   7.0 -0.5
Warning message:
In nls(y ~ f(x, a, b), data = dat2, start = c(a = 7, b = -0.5),  :
  singular gradient

请原谅糟糕的格式,首先发布在这里.x包含直方图的区间,y包含该直方图中每个区间的出现次数.dat2在14处截止,因为0计数箱会甩掉指数回归,我真的只需要适应那些前14位.那些数量超过14的箱子我有生理上的理由相信它们是特殊的.我最初得到的问题是无穷大,我没有得到,因为没有一个值为0.在给出了不同的帖子所提出的合适的起始值后,我得到了奇异的梯度误差.我看到的唯一其他帖子有更多变量,我尝试增加迭代次数但是没有成功.任何帮助表示赞赏.一个

Answer 1

1)线性化以获得起始值您需要更好的起始值:

# starting values
fm0 <- nls(log(y) ~ log(f(x, a, b)), dat2, start = c(a = 1, b = 1))

nls(y ~ f(x, a, b), dat2, start = coef(fm0))

赠送:

Nonlinear regression model
  model: y ~ f(x, a, b)
   data: x
        a         b 
4214.4228   -0.8106 
 residual sum-of-squares: 2388

Number of iterations to convergence: 6 
Achieved convergence tolerance: 3.363e-06

1a)同样,我们可以lm通过写作来获得初始值

y ~ a * exp(b * x)

如

y ~ exp(log(a) + b * x)

并记录两者的日志以获得log(a)和b中的模型线性:

log(y) ~ log(a) + b * x

可以使用lm以下方法解决:

fm_lm <- lm(log(y) ~ x, dat2)
st <- list(a = exp(coef(fm_lm)[1]), b = coef(fm_lm)[2])
nls(y ~ f(x, a, b), dat2, start = st)

赠送:

Nonlinear regression model
  model: y ~ f(x, a, b)
   data: dat2
       a        b 
4214.423   -0.811 
 residual sum-of-squares: 2388

Number of iterations to convergence: 6 
Achieved convergence tolerance: 3.36e-06

1b)我们也可以通过重新参数化来实现它.在这种情况下,只要我们根据参数变换转换初始值,a = 1和b = 1就可以工作.

nls(y ~ exp(loga + b * x), dat2, start = list(loga = log(1), b = 1))

赠送:

Nonlinear regression model
  model: y ~ exp(loga + b * x)
   data: dat2
  loga      b 
 8.346 -0.811 
 residual sum-of-squares: 2388

Number of iterations to convergence: 20 
Achieved convergence tolerance: 3.82e-07

所以b如图所示,a = exp(loga)= exp(8.346)= 4213.3

2)线性 另一种更容易的可能性是使用alg="plinear"在这种情况下线性进入的参数不需要起始值.在这种情况下,b=1问题的起始值似乎就足够了.

nls(y ~ exp(b * x), dat2, start = c(b = 1), alg = "plinear")

赠送:

Nonlinear regression model
  model: y ~ exp(b * x)
   data: dat2
        b      .lin 
  -0.8106 4214.4234 
 residual sum-of-squares: 2388

Number of iterations to convergence: 11 
Achieved convergence tolerance: 2.153e-06

Answer 2

请检查 minpack.lm 包中的 nlsLM 函数。这是 nls 的更强大版本，可以处理残差平方和为零的数据。

https://www.r-bloggers.com/a-better-nls/