読者です 読者をやめる 読者になる 読者になる

廿TT

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

指数分布を仮定した右打ち切りデータの平均の区間推定(ってこれでいいのかな?)

指数分布の平均の信頼区間

指数分布 - Wikipedia

信頼水準を(1-\alpha)とした指数分布の平均の信頼区間は
\displaystyle \left[  \frac{2n}{\chi^2_{2n, 1- \frac{\alpha}{2}}} \hat{\mu}, \frac{2n}{\chi^2_{2n, \frac{\alpha}{2}}} \hat{\mu}  \right]
で求められる。

このことは、\hat{\mu} がふつうの標本平均の場合は、指数分布からのランダム標本の和がガンマ分布に従うことを使って証明できる。

等を参照。

入門・演習 数理統計

入門・演習 数理統計

一方、観測に右打ち切りがある場合、指数分布の母平均の最尤推定量は
\displaystyle  \hat{\mu} = \left(\sum^{n}_{i=1} x_i\right) / \left(\sum^{n}_{i=1} d_i \right)
である。

ここで、 d_i は打ち切りのとき 0、打ち切りでないとき 1 の値を取るインジケータ。

ああああ: R による生存時間分析(指数分布を仮定したパラメトリックモデル) を参照。

観測に右打ち切りがある場合も、やはり指数分布の平均の信頼区間は
\displaystyle \left[  \frac{2n}{\chi^2_{2n, 1- \frac{\alpha}{2}}} \hat{\mu}, \frac{2n}{\chi^2_{2n, \frac{\alpha}{2}}} \hat{\mu}  \right]
でオッケーなんだろうか?

直感的には大丈夫そうな気がしてちょっと考えたんだけど、ちゃんとした証明はできなかった。

ケーススタディ

Gehan データ(Rと生存時間分析(1) を参照)で指数分布を仮定した右打ち切りデータの平均と信頼区間を求めてみる。

信頼水準は95%にした。

library(MASS)

gehan_c <-subset(gehan,treat=="control")
gehan_6 <-subset(gehan,treat=="6-MP")

exp_CI <- function(time,cens,alpha=0.05){
  n = length(time)
  mu_hat = sum(time)/sum(cens)
  num = 2 * n * mu_hat
  data.frame(mu_hat, lower= num/ qchisq((1-alpha/2),df=2*n), upper= num/ qchisq(alpha/2,df=2*n))
}

df1 <- rbind(
exp_CI(time =gehan_c$time, cens=gehan_c$cens),
exp_CI(time =gehan_6$time, cens=gehan_6$cens)
)

df1$treat <- c("control","6-MP")

library(ggplot2)
theme_set(theme_bw(base_size=15))
ggplot(df1, aes(x=treat,y=mu_hat))  +
  geom_bar(stat="identity", alpha=0.6)+
  labs(y=expression(hat(mu))) + 
  geom_errorbar(aes(ymin=lower,ymax=upper), colour="red3")
treat \hat{\mu} lower upper
control 8.666667 5.892184 14.00072
6-MP 39.888889 27.119154 64.43921

f:id:abrahamcow:20141124022653p:plain
(6-MP 投与群のほうがはっきり平均生存時間が長い。)

お願い

甚だ不躾なお願いではございますが、ちゃんとした証明ができる方いらっしゃいましたらご教示願えませんでしょうか。

参考になりそうな文献

鎌倉稔成『生存時間・再発事象分析:理論と応用』(pdf)

生存時間解析

生存時間解析

9 命題の成立確率と信頼区間