廿TT

譬如水怙牛過窓櫺 頭角四蹄都過了 因甚麼尾巴過不得

交代再生過程のある時刻における状態がポアソン分布に従うか否かシミュレーションしてみた.

状況

なんらかの機械を設置したとする.

この機械はある程度時間が経過すると故障する.

修理にはまたある程度の時間がかかる.

つまりこの機械は,

  • 状態1:稼働中
  • 状態2:修理中

と, 2 つの状態を交互に繰り返す.

f:id:abrahamcow:20150107212609p:plain

ある時刻 o にこの機械をみたとき, この機械が稼働中ならば成功, 修理中なら失敗とみなす.

この成功の回数はポアソン分布に従うか知りたい.

シミュレーション

指数分布に従って生起するイベントの発生回数はポアソン分布に従うことが知られているので, ここではワイブル分布を仮定することにした.

ワイブル分布 - Wikipedia

  • 故障が発生するまでにかかる時間は shape パラメータ = 3 のワイブル分布
  • 修理が完了するまでにかかる時間は shape パラメータ = 1/3 のワイブル分布

scale パラメータは両者ともに 1 とした.

時刻 o に機械の状態を調べることを 100 回繰り返した場合の成功の回数求めるシミュレーションを, 1,000 回繰り返した.

時刻 o は一様乱数でランダムに振った.

it <- 1000
suc <- numeric(100)
suc_C <- numeric(it)
RVs <- numeric(20000)
for(j in 1:it){
  for(i in 1:100){
    RVs[2*(1:10000)-1] <- rweibull(10000,shape=3) #稼働中; 奇数
    RVs[2*(1:10000)] <- rweibull(10000,shape=1/3) #修理中; 偶数
    tim <- cumsum(RVs)
    o <- runif(1,min(tim),max(tim))
    suc[i] <- which.max(tim>o) %% 2
  }
  suc_C[j] <- sum(suc)
}

tab1 <-table(suc_C)
x1 <-as.numeric(names(tab1))
bp1<- barplot(tab1)
l1 <-mean(suc_C)
points(bp1,dpois(x1,l1)*it)

f:id:abrahamcow:20150107215126p:plain

棒グラフが観測度数, ○印が理論度数.

みた感じ, ポアソン分布で近似して問題なさそうである.

ただし, 奥村先生のページ(カイ2乗検定)を参考に,

> tab1
suc_C
  3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25  26 
  1   1   3  12  27  42  63  80 108 129 108 106  87  79  45  42  31  18   8   5   1   1   2   1 

ポアソン分布に従うという帰無仮説についてカイ二乗検定をやってみると,

> chisq.test(tab1)

	Chi-squared test for given probabilities

data:  tab1
X-squared = 996.848, df = 23, p-value < 2.2e-16

p 値めっちゃ小さいのでこれはポアソン分布とは言えない.

二項分布だろうか.

bp1<- barplot(tab1)
points(bp1,dbinom(x1,prob=l1/100,size=100)*it)

f:id:abrahamcow:20150107222620p:plain

しかし, 二項検定の結果 p 値はわりと一様に分布している.

pvs <- sapply(1:it,function(i)binom.test(suc_C[i], n=100,p=l1/100)$p.value)
hist(pvs)

f:id:abrahamcow:20150107223525p:plain

二項分布でもない?

検定の仕方まちがえてるのだろうか?

検定と区間推定

今後の課題

正確な分布を調べるにはどうしたらいいのだろう?

また,

例えばPISA「盗難事件」問題ほかで挙げた6年間の事故発生件数

> x = c(762,792,795,794,849,834)
が独立なポアソン分布に従うという帰無仮説についてカイ2乗検定をしてみましょう。Rでは次のようにして簡単に計算できます。

> chisq.test(x)

Chi-squared test for given probabilities

data: x
X-squared = 6.2329, df = 5, p-value = 0.2842

カイ2乗検定

とあるが, chisq.test(x) でなんで「ポアソン分布に従うという帰無仮説についてのカイ二乗検定」になるのかよくわかってない.

追記

やはり検定の仕方がまちがっているとのこと。お恥ずかしい限りです。

ありがとうございます。大変助かりました。

この裏RjpWikiというブログにはいつも大変勉強させて頂いております。

そして恐れ多いことに奥村先生にも言及リンク頂きました。

この分布はポアソン?

ありがとうございます。がんばって勉強します。

7.2節 適合度検定