1

我正在尝试将数据的 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 选择合适的起始值。

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

谢谢！

4

1 回答 1

3

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

  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.
于 2012-06-28T07:55:20.897 回答