10

x下の図に示すように、時間（軸単位）とともに変化する周波数値があります。いくつかの正規化の後、これらの値は、いくつかの分布の密度関数のデータ ポイントとして表示される場合があります。

Q:これらの周波数ポイントがワイブル分布からのものであると仮定すると、そこから分布パラメーター を推測するためTに、どうすれば最適なワイブル密度関数をポイントに適合させることができますか?T

sample <- c(7787,3056,2359,1759,1819,1189,1077,1080,985,622,648,518,
            611,1037,727,489,432,371,1125,69,595,624)

plot(1:length(sample), sample, type = "l")
points(1:length(sample), sample)

ここに画像の説明を入力

更新します。誤解を招かないように、もう少し説明を加えたいと思います。時間（軸単位）とともに変化する周波数値がxあると言うことは、次のようなデータがあることを意味します。

  • 7787 の価値の実現 1
  • 価値 2 の 3056 の実現
  • 2359 の価値 3 の実現 ... など

私の目標（私が思うに、間違ったもの）に向けた何らかの方法は、これらの実現のセットを作成することです：

# Loop to simulate values 
set.values <- c()
for(i in 1:length(sample)){
  set.values <<- c(set.values, rep(i, times = sample[i]))
}

hist(set.values)
lines(1:length(sample), sample)
points(1:length(sample), sample)

ここに画像の説明を入力

で使用fitdistrしますset.values:

f2 <- fitdistr(set.values, 'weibull')
f2

それが間違った方法だと思う理由と、より良い解決策を探しているのはなぜRですか?

  • 上記の分布フィッティングアプローチでは、分布からの私の認識の完全なset.valuesセットであると想定されていますT

  • 私の最初の質問では、密度曲線の最初の部分のポイントを知っています-そのテールがわからず、テール（および密度関数全体）を推定したい

4

3 に答える 3

3

最初にすべてのポイントで試してください

最初のポイントを落として 2 回目のトライoptimこれは、ボックス内の値のセットに制約された最適な値を見つけるため に使用する前のように、より良い試みです(呼び出しでlowerとupperベクトルによって定義されoptimます)。ワイブル分布形状パラメーターに加えて、最適化の一部として x と y をスケーリングすることに注意してください。したがって、最適化する 3 つのパラメーターがあります。

残念なことに、すべてのポイントを使用すると、ほとんどの場合、制約ボックスの端に何かが見つかります。これは、ワイブルがすべてのデータに適していない可能性があることを示しています. 問題は 2 つの点です。それらは大きすぎます。最初のプロットでは、すべてのデータに当てはめようとしたことがわかります。

最初の 2 点を削除して、残りの点だけを合わせると、はるかによく適合します。これは2 番目のプロットで確認できます。これは適切だと思います。いずれにせよ、制約ボックスの内部での極小値です。

library(optimx)
sample <- c(60953,7787,3056,2359,1759,1819,1189,1077,1080,985,622,648,518,
            611,1037,727,489,432,371,1125,69,595,624)
t.sample <- 0:22

s.fit <- sample[3:23]
t.fit <- t.sample[3:23]

wx <- function(param) { 
  res <- param[2]*dweibull(t.fit*param[3],shape=param[1])
  return(res)
} 
minwx <- function(param){
  v <- s.fit-wx(param)
  sqrt(sum(v*v))
}

p0 <- c(1,200,1/20)
paramopt <- optim(p0,minwx,gr=NULL,lower=c(0.1,100,0.01),upper=c(1.1,5000,1))

popt <- paramopt$par
popt
rms <- paramopt$value
tit <- sprintf("Weibull - Shape:%.3f xscale:%.1f  yscale:%.5f rms:%.1f",popt[1],popt[2],popt[3],rms)

plot(t.sample[2:23], sample[2:23], type = "p",col="darkred")
lines(t.fit, wx(popt),col="blue")
title(main=tit)
于 2015-05-03T10:33:36.337 に答える
3

こちら で説明されているように、最尤パラメーターを直接計算できます。

# Defining the error of the implicit function
k.diff <- function(k, vec){
  x2 <- seq(length(vec))
  abs(k^-1+weighted.mean(log(x2), w = sample)-weighted.mean(log(x2), 
                                                            w = x2^k*sample))
}

# Setting the error to "quite zero", fulfilling the equation
k <- optimize(k.diff, vec=sample, interval=c(0.1,5), tol=10^-7)$min

# Calculate lambda, given k
l <- weighted.mean(seq(length(sample))^k, w = sample)

# Plot
plot(density(rep(seq(length(sample)),sample)))
x <- 1:25
lines(x, dweibull(x, shape=k, scale= l))
于 2015-05-06T07:43:02.180 に答える
1

データがワイブル分布からのものであると仮定すると、次のように形状とスケールのパラメーターを推定できます。

sample <- c(7787,3056,2359,1759,1819,1189,1077,1080,985,622,648,518,
        611,1037,727,489,432,371,1125,69,595,624)
 f<-fitdistr(sample, 'weibull')
 f

分散ワイブルかどうかわからない場合は、ks.test を使用することをお勧めします。これは、データが仮説分布からのものかどうかをテストします。データの性質に関する知識があれば、選択したいくつかの分布をテストして、どれが最も効果的かを確認できます。

あなたの例では、これは次のようになります。

 ks = ks.test(sample, "pweibull", shape=f$estimate[1], scale=f$estimate[2])
 ks

p 値は有意ではないため、データがワイブル分布からのものであるという仮説を棄却しません。

更新: ワイブルまたは指数関数のいずれかのヒストグラムは、データとよく一致しているように見えます。指数分布の方が適合すると思います。パレート分布は別のオプションです。 

f<-fitdistr(sample, 'weibull')
z<-rweibull(10000, shape= f$estimate[1],scale= f$estimate[2])
hist(z)

f<-fitdistr(sample, 'exponential')
z = rexp(10000, f$estimate[1]) 
hist(z)
于 2015-05-03T09:37:35.257 に答える