- html - 出于某种原因,IE8 对我的 Sass 文件中继承的 html5 CSS 不友好?
- JMeter 在响应断言中使用 span 标签的问题
- html - 在 :hover and :active? 上具有不同效果的 CSS 动画
- html - 相对于居中的 html 内容固定的 CSS 重复背景?
对于这个数据集:
dat = structure(list(x = c(5L, 5L, 5L, 5L, 10L, 10L, 10L, 10L, 15L,
15L, 15L, 15L, 17L, 17L, 17L, 17L, 20L, 20L, 20L, 20L, 20L, 20L,
20L, 20L, 22L, 22L, 22L, 22L, 24L, 24L, 24L, 24L, 25L, 25L, 25L,
25L, 27L, 27L, 27L, 27L, 30L, 30L, 30L, 30L, 35L, 35L, 35L, 35L),
y = c(2.2, 2.2, 1.95, 1.9, 4.1, 3.95, 3.75, 3.4, 5.15, 4.6,
4.75, 5.15, 3.7, 4.1, 3.9, 3.5, 7, 6.7, 6.7, 6.95, 4.95, 6, 6.45,
6.4, 7, 4.45, 6.15, 6.4, 7, 6.6, 6.7, 7, 4.5, 4.7, 5.75, 4.35,
5.4, 5.15, 5.7, 5.7, 0, 0, 0.5, 0, 0, 0, 0, 0)), .Names = c("x", "y"),
row.names = c(6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L,
15L, 16L, 17L, 34L, 35L, 36L, 37L, 18L, 19L, 20L, 21L, 38L, 39L,
40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 22L, 23L, 24L,
25L, 50L, 51L, 52L, 53L, 26L, 27L, 28L, 29L, 30L, 31L, 32L, 33L),
class = "data.frame")
beta.reg<-function(x, Yopt,Tmin,Topt,Tmax, b1) {
Yopt*((x-Tmin)/(Topt-Tmin))^(b1*(Topt-Tmin)/(Tmax-Topt))*((Tmax-x) / (Tmax-Topt)) ^ b1
}
mod <- nls(y ~ beta.reg(x, Yopt,Tmin,Topt,Tmax, b1), data=dat,
start=c(Yopt=6, Tmin=0.1, Topt=24, Tmax=30, b1=1),
control=nls.control(maxiter=800))
Error en numericDeriv(form[[3L]], names(ind), env) : Missing value or an infinity produced when evaluating the model
rnorm<-(10)
y <- c(20,60,70,49,10)
rnorm<-(10)
y <- c(20,60,70,49,10)
dat<-data.frame(x = rep(c(15,20,25,30,35), times=5),
rep = as.factor(rep(1:5, each=5)),
y = c(y+rnorm(5), y+rnorm(5),y+rnorm(5),y+rnorm(5),y+rnorm(5)))
R version 3.1.1 (2014-07-10)
Platform: x86_64-pc-linux-gnu (64-bit)
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] nlme_3.1-118 latticeExtra_0.6-26 RColorBrewer_1.0-5 lattice_0.20-29
loaded via a namespace (and not attached):
[1] grid_3.1.1 tools_3.1.1
最佳答案
这里有很多问题,我怀疑它是否可以在 SO 帖子中充分涵盖,但这应该可以帮助您入门。
首先,看起来您想要 Tmax < max(dat$x)
,例如 < 35。这会导致问题,因为然后 Tmax - x < 0
用于 x
的某些值,并且当您尝试将负数提高到幂(在公式的第二项中)时,您会得到NA
的。这是错误消息的原因。
其次,非线性模型的收敛取决于模型公式和数据,因此该过程收敛于一组数据而不是另一组数据的事实完全无关紧要。
第三,非线性建模迭代地最小化作为参数函数的残差平方和。如果 RSS 表面具有 局部最小值 ,并且您的 start
接近于 1,则算法会找到它。但只有 全局最小值 才是真正的解决方案。你的问题有很多很多局部最小值。
第四,nls(...)
默认使用高斯牛顿法。众所周知,高斯牛顿在移动参数(添加到预测器或从预测器中减去的参数,因此在您的情况下 Tmin
和 Tmax
)时不稳定。幸运的是,minpak.lm
包实现了 Levenberg Marquardt 方法,它在这些条件下更加稳定。该包中的 nlsLM(...)
函数使用与 nls(...)
相同的调用序列,并返回 nls
类型的对象,因此该类对象的所有方法也能正常工作。用那个。
第五,非线性回归(实际上都是最小二乘回归)的一个基本假设是残差是正态分布的。因此,您必须使用 Q-Q 图验证任何解决方案。
第六,你的模型有一组反常的特征。当 Tmin -> -Inf
模型中的第一项接近 1
时。事实证明,与 Tmin
小于 min(dat$x)
的任何其他值相比,这产生的 RSS 更低,因此所有算法都倾向于将 Tmin
驱动为较大的负值。你可以很容易地看到这一点,如下所示:
library(minpack.lm)
mod <- nlsLM(y ~ beta.reg(x, Yopt,Tmin,Topt,Tmax, b1), data=dat,
start=c(Yopt=6,Tmin=0,Topt=24,Tmax=50, b1=1),
control=nls.lm.control(maxiter=1024,maxfev=1024))
coef(summary(mod))
# Estimate Std. Error t value Pr(>|t|)
# Yopt 6.347019 0.2919686 21.73870235 8.055342e-25
# Tmin -155.530098 2204.0011003 -0.07056716 9.440694e-01
# Topt 21.157545 0.6702713 31.56564484 2.240134e-31
# Tmax 35.000000 11.4838614 3.04775537 3.933164e-03
# b1 3.321326 9.1844548 0.36162468 7.194035e-01
sum(residuals(mod)^2)
# [1] 50.24696
par(mfrow=c(1,2))
plot(y~x,dat)
with(as.list(coef(mod)),curve(beta.reg(x, Yopt,Tmin,Topt,Tmax, b1),add=TRUE))
qqnorm(residuals(mod))
Tmin
和
b1
的估计都很差,而且
Tmin
的值在物理上没有意义,这是数据问题,而不是拟合问题。
Tmin
、
Tmax
和
b1
进行网格搜索来看到这一点(省略
Yopt
和
Topt
以节省时间,因为无论起点如何,这些参数都得到了很好的估计)。
init <- c(Yopt=6, Topt=24)
grid <- expand.grid(Tmin= seq(0,4,len=100),
Tmax= seq(35,100,len=10),
b1 = seq(1,10,len=10))
mod.lst <- apply(grid,1,function(gr){
nlsLM(y ~ beta.reg(x, Yopt,Tmin,Topt,Tmax, b1), data=dat,
start=c(init,gr),control=nls.control(maxiter=800)) })
rss <- sapply(mod.lst,function(m)sum(residuals(m)^2))
mod <- mod.lst[[which.min(rss)]] # fit with lowest RSS
coef(summary(mod))
# Estimate Std. Error t value Pr(>|t|)
# Yopt 6.389238 0.2534551 25.208557840 2.177168e-27
# Topt 22.636505 0.5605621 40.381798589 7.918438e-36
# Tmin 35.000002 104.6221159 0.334537316 7.396005e-01
# Tmax 36.234602 133.4987344 0.271422809 7.873647e-01
# b1 -41.512912 7552.0298633 -0.005496921 9.956395e-01
sum(residuals(mod)^2)
# [1] 34.24019
plot(y~x,dat)
with(as.list(coef(mod)),curve(beta.reg(x, Yopt,Tmin,Topt,Tmax, b1),add=TRUE))
qqnorm(residuals(mod))
x
之外没有为
(Tmin,Tmax)
定义。由于您有数据输出到
x=35
,拟合算法永远不会产生
Tmax < 35
(如果它收敛)。一种解决此问题的方法会稍微改变您的模型函数,使其在该范围之外裁剪为 0。 (不过,我不知道根据您的问题的物理性质,这是否合法......)。
beta.reg<-function(x, Yopt,Tmin,Topt,Tmax, b1) {
ifelse(x>Tmax,0,
ifelse(x<Tmin,0,
Yopt*((x-Tmin)/(Topt-Tmin))^(b1*(Topt-Tmin)/(Tmax-Topt))*((Tmax-x) / (Tmax-Topt)) ^ b1
))
}
coef(summary(mod))
# Estimate Std. Error t value Pr(>|t|)
# Yopt 6.1470413 0.21976766 27.970636 3.202940e-29
# Tmin -52.8172658 184.16899439 -0.286787 7.756528e-01
# Topt 23.0777898 0.63750721 36.200045 7.638121e-34
# Tmax 30.0039413 0.02529877 1185.984187 1.038918e-98
# b1 0.5966129 0.32439982 1.839128 7.280793e-02
sum(residuals(mod)^2)
# [1] 28.10144
par(mfrow=c(1,2))
plot(y~x,dat)
with(as.list(coef(mod)),curve(beta.reg(x, Yopt,Tmin,Topt,Tmax, b1),add=TRUE))
qqnorm(residuals(mod))
qqline(residuals(mod))
b1
的估计要好得多(并且与早期模型函数的估计非常不同)。残差仍然不正常,但在这种情况下,我想检查数据是否有异常值。
关于r - nls - 收敛错误,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/26681385/
在 Eclipse 中运行我的程序时出现错误ORA-12705:无法访问 NLS 数据文件或指定的环境无效。 到目前为止我尝试过的事情: 我已经尝试在注册表的\HKEY_LOCAL_MACHINE\S
在 R 中工作,我正在尝试使用 nls() 将一些数据拟合到以下模型: y ~ c - a * exp(-b * x) 我的数据: x <- c(8, 8, 10, 10, 10, 10, 12, 1
我尝试搜索有关此主题的其他线程,但没有一个修复对我有用。我有一个自然实验的结果,我想显示符合指数分布的事件连续发生的次数。我的 R shell 粘贴在下面 f x [1] 1 2 3 4
对于这个数据集: dat = structure(list(x = c(5L, 5L, 5L, 5L, 10L, 10L, 10L, 10L, 15L, 15L, 15L, 15L, 17L, 17L
我正在使用 NLS.initializeMessages(BUNDLE_NAME, XYZ.class) 加载常量。 当“BUNDLE_NAME”和“XYZ.class”存在于同一个项目中时,它可以工
使用线性模型函数 lm() 多项式公式可以包含这样的快捷符号: m <- lm(y ~ poly(x,3)) 这是一个快捷方式,可以让用户不必创建 x^2 和 x^3 变量或将它们输入到公式中,例如
想象这两个列表如下: seq<-seq(1,16) L1<-list(A=seq, B=seq, C=seq, D=seq) v1=c(11,15,17,19,21,22,24,25,26,27,28
我想拟合一个线性平台 (nls) 模型,该模型将高度描述为年龄的函数,并且我想测试区域间模型的任何参数是否存在显着差异。 这是我目前所拥有的: # Create data df1 head (df1
我有一个输入日期(比如 sysdate),我想获取上一个星期一的日期。我试过 select trunc(sysdate, 'D') from dual; 但它依赖于 NLS。此外,我不想通过名称检查结
我正在尝试使用 R 中的 nls 拟合矩形双曲线。 curve.nlslrc = nls(photolrc ~ (1/(2*theta))*(AQY*PARlrc+Am-sqrt((AQY*PARlr
我正在使用 nls 拟合一些指数数据. 我正在使用的代码是: fit <- nls(y ~ expFit(times, A, tau, C), start = c(A=100, tau=-3, C=0
尝试拟合我的数据时遇到问题。 这是数据: x = c(1, 1.071519305, 1.148153621, 1.230268771, 1.318256739, 1.412537545, 1.513
我正在尝试制作 nls适合一个有点复杂的表达式,其中包括两个积分,其中两个拟合参数在其上限中。 我得到了错误 "Error in nlsModel(formula, mf, start, wts) :
有没有办法限制 NLS 系数在 R 中的取值范围?我知道我的数据应该存在的曲线形状;然而,NLS 无法通过产生 < 1 的功率系数来产生这样的曲线。 从本质上讲,我正在尝试为一组幼树茎(树苗)数据生成
我想将以下函数拟合到我的数据中: f(x) = Offset+Amplitudesin(FrequencyT+Phase), 或根据 Wikipedia : f(x) = C+alphasin(ome
我不明白为什么我不能为这些数据使用 nls 函数。 我尝试了很多不同的起始值,但总是出现相同的错误。 这是我一直在做的事情: expFct2 = function (x, a, b,c) { a*
我一直在研究一个曲线拟合脚本,该脚本将 3 个指数修改的高斯 (EMG) 拟合到卷积曲线。我的基函数类似于高斯分布,但包括第三个参数(前两个是 mu 和 sigma ),它确定函数的指数分量的权重。
我们的 java 类调用 PLSQL proc,它以由 NLS_DATE_FORMAT 定义的默认格式返回日期。我们的应用程序为国际化设置了自己的语言环境,但我希望日期格式保持为“DD-MON-RR”
谁能给我一个很好的解释,说明参数“算法”在 R 中的 nls 函数中的作用? 另外,公式是如何工作的?我知道它使用了 tilda,但我真的找不到对它的实际解释。 另外,起始值有多重要?我是否需要尝试多
我正在尝试使用 nls 将高斯峰拟合到 R 中的密度图.当我使用以下等式时: fit :2:0: unexpected end of input 1: ~ ^ 谁能指出我哪里出错了? 最佳答案 试
我是一名优秀的程序员,十分优秀!