【问题标题】:R Nonlinear Least Squares (nls) Model FittingR 非线性最小二乘 (nls) 模型拟合
【发布时间】:2012-06-29 08:56:35
【问题描述】:

我正在尝试将数据的 G 函数中的信息拟合到以下数学模式:y = A / ((1 + (B^2)*(x^2))^((C+1 )/2)) 。该图的形状可以在这里看到:

http://www.wolframalpha.com/input/?i=y+%3D+1%2F+%28%281+%2B+%282%5E2%29*%28x%5E2%29%29%5E%28%282%2B1%29%2F2%29%29

这是我一直在做的一个基本示例:

data(simdat)

library(spatstat)

simdat.Gest <- Gest(simdat) #Gest is a function within spatstat (explained below)

Gvalues <- simdat.Gest$rs

Rvalues <- simdat.Gest$r

GvsR_dataframe <- data.frame(R = Rvalues, G = rev(Gvalues))

themodel <- nls(rev(Gvalues) ~ (1 / (1 + (B^2)*(R^2))^((C+1)/2)), data = GvsR_dataframe, start = list(B=0.1, C=0.1), trace = FALSE)

“Gest”是“spatstat”库中的一个函数。它是 G 函数或最近邻函数,它显示了独立轴上粒子之间的距离,以及在从属轴上找到最近邻粒子的概率。因此,它从 y=0 开始,在 y=1 处达到饱和点。

如果您绘制 simdat.Gest,您会注意到曲线是“s”形的,这意味着它从 y = 0 开始并在 y = 1 结束。出于这个原因,我反转了向量 Gvalues,它是因变量。因此,信息处于适合上述模型的正确方向。

您可能还注意到我已经自动设置了 A = 1。这是因为 G(r) 总是在 1 处饱和,所以我没有费心将它保留在公式中。

我的问题是我不断收到错误。对于上面的示例,我收到此错误:

Error in nls(rev(Gvalues) ~ (1/(1 + (B^2) * (R^2))^((C + 1)/2)), data = GvsR_dataframe,  : 
  singular gradient

我也遇到了这个错误:

Error in nls(Gvalues1 ~ (1/(1 + (B^2) * (x^2))^((C + 1)/2)), data = G_r_dataframe,  : 
  step factor 0.000488281 reduced below 'minFactor' of 0.000976562

我不知道第一个错误来自哪里。然而,我认为第二个发生是因为我没有为 B 和 C 选择合适的起始值。

我希望有人可以帮助我找出第一个错误的来源。另外,选择起始值以避免第二个错误的最有效方法是什么?

谢谢!

【问题讨论】:

  • 当您输入 simdat.Gest &lt;- Gest(simdat) 时,您是在告诉我们您有一个名为 Gest 的函数。但是你没有给我们。我认为使用rnorm 创建测试数据集并不难,尽管理想情况下我们会得到前二十行,但我们确实需要知道Gest 在做什么。
  • 您在该公式中使用的名称 ('Gvalues') 与 data.frame ('G') 中使用的名称不同
  • 那是我的错误,我现在已经编辑了帖子。 “Gest”是“spatstat”库中的一个函数。 Gest 是最近邻函数,它显示独立轴上粒子之间的距离,以及在从属轴上找到最近邻粒子的概率。因此,它从 y=0 开始,在 y=1 处达到饱和点。
  • 另外,我尝试过使用 nls.lm,它也让我很伤心。

标签: r model plot nls model-fitting


【解决方案1】:

如前所述,您的问题很可能是起始值。您可以使用两种策略:

  1. 使用蛮力查找起始值。请参阅包 nls2 了解执行此操作的函数。
  2. 尝试合理猜测起始值。 根据您的值,可以对模型进行线性化。

G = (1 / (1 + (B^2)*(R^2))^((C+1)/2))

ln(G)=-(C+1)/2*ln(B^2*R^2+1)

如果 B^2*R^2 很大,则变为大约。 ln(G) = -(C+1)*(ln(B)+ln(R)),这是线性的。

如果 B^2*R^2 接近 1,它大约是。 ln(G) = -(C+1)/2*ln(2),是常数。

(请检查错误,由于足球比赛,昨晚很晚。)

在提供其他信息后

编辑: 数据看起来像遵循累积分布函数。如果它像鸭子一样嘎嘎叫,那它很可能是一只鸭子。事实上,?Gest 声明了一个 CDF 估计值。

library(spatstat)
data(simdat)
simdat.Gest <- Gest(simdat)
Gvalues <- simdat.Gest$rs
Rvalues <- simdat.Gest$r
plot(Gvalues~Rvalues)

#let's try the normal CDF
fit <- nls(Gvalues~pnorm(Rvalues,mean,sd),start=list(mean=0.4,sd=0.2))
summary(fit)
lines(Rvalues,predict(fit))
#Looks not bad. There might be a better model, but not the one provided in the question.

【讨论】:

  • 这很有帮助,谢谢!我可能会探索其他模型,并尝试优化我之前使用的模型,但这是一个很好的起点。
  • Roland,事实证明,使用 CDF 对我的信息进行建模要简单得多,但并没有提供我希望的那么多有用的信息。我最初的问题中的函数,如果我能够用数据对其进行建模,将提供更多有用的信息。你知道我可以使用什么策略来找到我上面提到的错误的根源吗?
  • 试图拟合一个不能很好地描述数据“形状”的模型通常会带来优化问题。这基本上就是错误消息告诉您的内容。非线性回归模型的选择应该(我几乎说必须)基于物理或数学考虑。您选择模型的理由是什么?我怀疑您正在尝试解决一些未知的其他问题,使用其他工具可以更好地解决。
猜你喜欢
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2012-08-29
  • 2013-07-28
  • 2017-02-24
  • 1970-01-01
相关资源
最近更新 更多