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 で終わることを意味します。変数。したがって、情報は上記のモデルに適合する正しい方向にあります。

また、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

最初のエラーがどこから来ているのか、私には手がかりがありません。しかし、2 つ目は、B と C に適切な開始値を選択しなかったために発生したと思います。

最初のエラーがどこから来たのかを誰かが見つけてくれることを願っていました。また、2 番目のエラーを回避するために開始値を選択する最も効果的な方法は何ですか?

ありがとう!

4

1 に答える 1

3

前述のように、問題はおそらく開始値です。使用できる戦略は 2 つあります。

  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 に答える