統計学 学習ノート

はじめに

「試験の途中でまだ動いている電球のデータは、捨ててしまっていいの？」

電球の寿命を調べるために、10個を同時に点灯させたとします。全部が切れるまで待つと、何か月もかかります。そこで4個が切れた時点で試験を止めました。残りの6個は「まだ切れていない」ということしか分かりません。この6個を捨てて4個だけで平均をとると、丈夫な電球を無視したことになるので、寿命を短く見積もってしまいます。

別の場面も考えてみましょう。ある店の会員の来店回数を調べたら、平均は2回でした。ただし、この記録に載るのは1回以上来た人だけです。入会しただけで1回も来ていない人は、記録のどこにも現れません。

実は、この2つの例は、同じ「データが足りない」でも性質がまったく違います。前者はいることが分かっていて値が足りないデータで、後者はいることすら分からないデータです。だから、同じ公式では扱えません。正しい分析の鍵は、「このデータが手に入る確率はいくらか」という尤度を、欠け方に合わせて書くことです。

そこで登場するのが、今回の主役 欠測、打ち切り、トランケーション の3つです。

前回の記事では、順位だけを使う検定を導きました。今回の記事では、3種類の不完全さの区別、欠測の尤度の分解、打ち切りの尤度と正確な信頼区間、トランケーションの条件付き分布、そして尤度の寄与の書き分けについて解説します。

3種類の不完全さを見分ける

はじめに、不完全データの3つの種類について解説します。

上限が150kgの体重計で、健康診断をしていると考えてください。不完全なデータの生まれ方は、次の3通りです。

【体重計で見る3つの不完全さ】

1. 欠測： 記録用紙が濡れて、数字が読めなくなった人がいる
2. 打ち切り： 150kgを超える人は「150kg以上」とだけ記録される
3. トランケーション： 150kgを超える人は会場に入れず、記録に何も残らない

欠測（欠損、missing）は、値があるはずなのに手に入らないことです。準1級で学んだ空欄のことです。

打ち切り（censoring）は、個体は観測されたが、値はある範囲に入ることしか分からないことです。

トランケーション（切断、truncation）は、値がある範囲に入る個体が、標本にまったく現れないことです。

種類	個体がいると分かるか	値について分かること
欠測	分かる	何も分からない（他の変数は分かる）
打ち切り	分かる	値が入る範囲は分かる
トランケーション	分からない	標本にいる個体は値も分かる

要するに、打ち切りとトランケーションの違いは、「人数が分かるかどうか」です。打ち切りでは「150kg以上の人が3人いた」と言えますが、トランケーションでは何人いたか分かりません。

せっかくなので、別の例も見てみましょう。

【実際の場面での区別】

• 血液検査で検出限界の0.1未満は「0.1未満」と記録する： 打ち切り（小さい側なので左打ち切り）
• 研究の終了時に生きていた患者の生存時間： 打ち切り（大きい側なので右打ち切り）
• 調査の開始前に閉店した店が、名簿に載っていない： トランケーション（左トランケーション）
• 合格者の入試の点数だけが手元にある： トランケーション

試験のポイント：問題文で「人数が分かるか」を確かめます。人数が分かれば打ち切り、分からなければトランケーションです。

以上が、3種類の不完全さについての解説です。

欠測：観測データの尤度を分解する

次に、欠測があるときの尤度について説明します。

準1級で学んだように、欠測メカニズムがMAR（観測値だけで欠測の確率が決まる）なら、メカニズムを尤度から外せます。これを無視可能と呼びました。ここでは、無視可能なときに最尤推定量を実際に導きます。

2つの変数を考えます。$x$ は $n$ 人全員が観測され、$y$ は最初の $r$ 人だけが観測されたとします。欠測の確率は $x$ だけで決まるとします。同時密度を $f(x, y) = f(x; \phi)\, f(y \mid x; \psi)$ と分解すると、観測データの尤度は次のとおりです。

$$L = \prod_{i=1}^{r} f(x_i; \phi)\, f(y_i \mid x_i; \psi) \prod_{i=r+1}^{n} \int f(x_i; \phi)\, f(y \mid x_i; \psi)\, dy$$

後ろの積分は、$y$ の密度を積分して1になるので、$f(x_i; \phi)$ だけが残ります。つまり、次のように2つに分かれます。

$$L = \underbrace{\prod_{i=1}^{n} f(x_i; \phi)}_{\text{全員の } x} \times \underbrace{\prod_{i=1}^{r} f(y_i \mid x_i; \psi)}_{r \text{ 人の回帰}}$$

$\phi$ と $\psi$ が別々の母数なら、それぞれを最大化すればよいことになります。2人の係が別々の帳簿をつけるように、全員分の $x$ の帳簿と、$r$ 人分の回帰の帳簿に分かれたのです。

$(x, y)$ が2変量正規分布なら、$E(y \mid x) = \alpha + \beta x$ です。$\mu_x$ は全員の平均 $\bar{x}$ で推定し、$\alpha, \beta$ は $r$ 人の回帰で推定します。$\mu_y = \alpha + \beta\mu_x$ なので、最尤推定量の不変性から、次の式が得られます。

$$\hat{\mu}_y = \bar{y}_r + \hat{\beta}(\bar{x} - \bar{x}_r)$$

$\bar{x}_r, \bar{y}_r$ は $r$ 人の平均です。言葉で言うと、「観測できた人の平均を、$x$ の偏りの分だけ回帰で直した値」です。

例題：MARのもとでの平均の最尤推定

100人の中間テストの平均は60点だった。期末テストは中間テストの点によって受験の確率が変わり、受験した40人の中間テストの平均は66点、期末テストの平均は70点だった。受験者での回帰係数は $\hat{\beta} = 0.8$ である。$(x, y)$ を2変量正規分布とする。〔1〕$\mu_y$ の最尤推定値を求めよ。〔2〕受験者の平均70点をそのまま使うことの問題点を述べよ。

【解答】

〔1〕$\hat{\mu}_y = 70 + 0.8 \times (60 - 66) = 65.2$ である。

〔2〕受験者は中間テストの点が高い人に偏っている（66点 > 60点）。欠測は $x$ だけに依存するのでMARであり、$y$ の平均をそのまま使うと、$x$ の偏りの分だけ過大になる。

一方で、受験の確率が期末テストの点そのものに依存するMNARなら、この分解は成り立ちません。欠測の確率が $y$ を含むので、積分の外に出せないからです。メカニズムのモデルを同時に書く必要があり、その妥当性はデータから確かめられません。

以上が、欠測の尤度についての説明です。

打ち切り：「ここまでは大丈夫だった」を尤度に入れる

続いて、打ち切りの理論について解説します。

打ち切りには、値の位置による分類と、打ち切りの決まり方による分類があります。

【打ち切りの種類】

• 右打ち切り： 値が $c$ より大きいことだけが分かる
• 左打ち切り： 値が $c$ より小さいことだけが分かる
• 区間打ち切り： 値が $(a, b]$ に入ることだけが分かる
• I型打ち切り： 前もって決めた時刻 $c$ で観察を止める
• II型打ち切り： 前もって決めた個数 $r$ の故障が起きたら止める
• ランダム打ち切り： 個体ごとの打ち切り時刻が確率変数である

右打ち切りの尤度を導きましょう。寿命 $T$ の密度を $f$、生存関数を $S$ とします。打ち切り時刻 $C$ は $T$ と独立で、密度 $g$、生存関数 $\bar{G}$ をもつとします。観測できるのは、$Y = \min(T, C)$ と、$\delta = 1\{T \le C\}$ です。

$\delta = 1$ の個体は、時刻 $y_i$ で故障し、それまで打ち切られていません。$\delta = 0$ の個体は、時刻 $y_i$ で打ち切られ、それまで故障していません。したがって、次のようになります。

$$L = \prod_i \{f(y_i)\bar{G}(y_i)\}^{\delta_i}\{g(y_i)S(y_i)\}^{1 - \delta_i} = \prod_i f(y_i)^{\delta_i} S(y_i)^{1 - \delta_i} \times \prod_i \bar{G}(y_i)^{\delta_i} g(y_i)^{1 - \delta_i}$$

ここがポイントです。$g$ と $\bar{G}$ が寿命の母数を含まないなら、後ろの積は定数です。この仮定を 無情報な打ち切り と言います。これは、欠測の無視可能性と同じ発想です。逆に、「体調が悪い人ほど途中で来なくなる」ような打ち切りでは、この分解は成り立ちません。

指数分布なら、準1級で見たように、$d = \sum \delta_i$ として $\hat{\lambda} = d/\sum y_i$ です。1級では、情報量がどれだけ減るかまで計算します。対数尤度 $\ell = d\log\lambda - \lambda\sum y_i$ を2回微分すると、$-\ell'' = d/\lambda^2$ です。時刻 $c$ のI型打ち切りでは $E(d) = n(1 - e^{-\lambda c})$ なので、次のとおりです。

$$I(\lambda) = \frac{n(1 - e^{-\lambda c})}{\lambda^2}$$

完全なデータの $n/\lambda^2$ に比べて、情報は $1 - e^{-\lambda c}$ 倍に減ります。例えば、$c$ が平均寿命に等しいとき、$\lambda c = 1$ なので63%です。打ち切りは、故障の起きた数の割合だけ情報を残す、と言えます。

1. II型打ち切りの正確な分布

II型打ち切りでは、指数分布のとき正確な信頼区間が作れます。$n$ 個のうち先に故障した $r$ 個の時刻を $t_{(1)} < \cdots < t_{(r)}$ とします。順序統計量の同時密度は、残りの $n - r$ 個が $t_{(r)}$ より長生きすることも考えると、次の形です。

$$\frac{n!}{(n-r)!} \prod_{i=1}^{r} \lambda e^{-\lambda t_{(i)}} \cdot e^{-\lambda(n-r)t_{(r)}} = \frac{n!}{(n-r)!} \lambda^{r} e^{-\lambda T}$$

$T = \sum_{i=1}^{r} t_{(i)} + (n - r)t_{(r)}$ を 総試験時間 と呼びます。みんなの稼働時間の合計です。最尤推定量は、$\hat{\lambda} = r/T$ です。

ここで、$t_{(0)} = 0$ として、$W_i = (n - i + 1)(t_{(i)} - t_{(i-1)})$ と変換します。$\sum_i W_i = T$ で、ヤコビアンは定数なので、密度は $\lambda^r e^{-\lambda\sum w_i}$ に比例します。つまり、$W_1, \dots, W_r$ は独立に指数分布に従います。

$i$ 番目の故障を待っている間は $n - i + 1$ 個が動いていて、指数分布には無記憶性があるので、全員が毎回まっさらな状態から競争し直す、というイメージです。指数分布の和はガンマ分布なので、次が成り立ちます。

$$2\lambda T \sim \chi^2_{2r}$$

例題：II型打ち切りの信頼区間

電球10個を同時に点灯し、4個目が切れた時点で試験を止めた。切れた時刻は 12, 20, 31, 45（時間）だった。寿命は指数分布とする。〔1〕総試験時間と平均寿命 $\theta = 1/\lambda$ の最尤推定値を求めよ。〔2〕$\theta$ の95%信頼区間を求めよ。ただし $\chi^2_{0.975}(8) = 2.18$、$\chi^2_{0.025}(8) = 17.53$ とする。

【解答】

〔1〕$T = 12 + 20 + 31 + 45 + 6 \times 45 = 378$ である。$\hat{\theta} = T/r = 378/4 = 94.5$ 時間である。

〔2〕$2T/\theta \sim \chi^2_8$ だから、$2.18 < 756/\theta < 17.53$ を解いて、次を得る。

$$\frac{756}{17.53} = 43.1 < \theta < \frac{756}{2.18} = 346.8$$

切れた4個の平均は27時間で、推定値の94.5時間とは大きく違います。わずか4回の故障では区間も広くなります。

2. 指数分布以外の場合

分布を決めないときも、同じ尤度の発想が使えます。故障の時刻 $t_1 < t_2 < \cdots$ で、時刻 $t_j$ の直前に $n_j$ 個が動いていて、そのうち $d_j$ 個が故障したとします。ハザードを $h_j$ とすると、尤度は $\prod_j h_j^{d_j}(1 - h_j)^{n_j - d_j}$ です。各因子は二項分布の形なので、$\hat{h}_j = d_j/n_j$ です。これを $\hat{S}(t) = \prod_{t_j \le t}(1 - \hat{h}_j)$ に入れたものが、準1級のカプラン・マイヤー推定量でした。つまり、ノンパラメトリックな最尤推定量です。

補足（ここは読み飛ばしてOK）： 打ち切りを欠測の一種と考えて、EMアルゴリズムで解くこともできます。指数分布なら、打ち切りの値を $E(T \mid T > c) = c + 1/\lambda$ で埋めて平均をとる操作をくり返します。収束先は $1/\lambda = \sum y_i/d$ で、直接の最尤推定値と一致します。

試験のポイント：総試験時間には、打ち切られた個体の時間を入れます。$2\lambda T \sim \chi^2_{2r}$ の自由度は個数 $n$ ではなく、故障数 $r$ の2倍です。

以上が、打ち切りについての解説です。

トランケーション：見えない人の分だけ割り直す

さらに、トランケーションの理論について説明します。

値が集合 $A$ に入るときだけ標本に現れるなら、手元のデータの分布は、条件付き分布です。

$$f_A(t) = \frac{f(t)}{P(T \in A)}, \qquad t \in A$$

分子は元の密度と同じで、見えない部分の確率で割って全体を1にそろえています。ピザの一部が切り取られたとき、残りを全体と考え直すようなものです。

1. 左トランケーションと遅延エントリー

老人ホームの入居者の寿命を調べると、入居前に亡くなった人は、記録に現れません。これは、入居した年齢 $v_i$ での左トランケーションで、遅延エントリー とも呼ばれます。個体ごとの尤度の寄与は、$f(t_i)/S(v_i)$ です。

指数分布なら、この寄与は $\lambda e^{-\lambda(t_i - v_i)}$ で、$\hat{\lambda} = n/\sum(t_i - v_i)$ となります。無記憶性のおかげで、入ってからの時間だけを見ればよいのです。例えば、$(v, t)$ が (2, 5), (0, 4), (3, 9), (1, 2), (4, 6) なら、$\sum(t_i - v_i) = 16$ で、平均の推定値は3.2です。$t$ をそのまま平均すると5.2となり、長生きした人が選ばれた分だけ、長く見えてしまいます。

2. 切断正規分布

入試の合格者のように、点数が $c$ 以上の人だけが標本に残る場合を考えます。$X \sim N(\mu, \sigma^2)$、$\alpha = (c - \mu)/\sigma$ として、$\phi'(z) = -z\phi(z)$ を使って積分すると、次の公式が得られます。

$$E(X \mid X > c) = \mu + \sigma\lambda(\alpha), \qquad \lambda(\alpha) = \frac{\phi(\alpha)}{1 - \Phi(\alpha)}$$

$\lambda(\alpha)$ は 逆ミルズ比 と呼ばれます。分散は $\sigma^2\{1 + \alpha\lambda(\alpha) - \lambda(\alpha)^2\}$ です。

具体的に見てみましょう。受験者の点数が $N(60, 10^2)$ で、合格点が70点なら、$\alpha = 1$、$\lambda(1) = 0.242/0.159 = 1.525$ です。合格者の平均は75.3点、標準偏差は4.5点になります。合格者だけを見ると、平均は高く、ばらつきは半分以下に見えるのです。

同じ場面を、打ち切りと比べましょう。70点未満の人が全員「70点」と記録された場合、人数は全員分あるので、平均は $0.841 \times 70 + 0.159 \times 75.25 = 70.8$ 点です。同じ体重計でも、人数が分かるかどうかで、答えがこれだけ変わります。

例題：零切断ポアソン分布

ある店の会員のうち、1回以上来店した人だけが記録に現れる。来店回数はポアソン分布に従うと仮定する。記録された120人の平均来店回数は2.0回だった。〔1〕記録される来店回数の確率関数を書け。〔2〕$\lambda$ の尤度方程式を導け。〔3〕$\hat{\lambda}$ を、$\lambda = 2$ から始めるニュートン法で求め、記録に現れない会員の人数を推定せよ。

【解答】

〔1〕$P(X = k \mid X \ge 1) = \dfrac{e^{-\lambda}\lambda^{k}}{k!\,(1 - e^{-\lambda})}$、$k = 1, 2, \dots$ である。

〔2〕対数尤度は $\ell = \sum x_i \log\lambda - n\lambda - n\log(1 - e^{-\lambda}) + \text{定数}$ である。微分すると、次を得る。

$$\ell'(\lambda) = \frac{\sum x_i}{\lambda} - n - \frac{ne^{-\lambda}}{1 - e^{-\lambda}} = \frac{n\bar{x}}{\lambda} - \frac{n}{1 - e^{-\lambda}} = 0$$

よって、尤度方程式は $\bar{x} = \lambda/(1 - e^{-\lambda})$ である。

〔3〕$h(\lambda) = \lambda/(1 - e^{-\lambda}) - 2$ の根を求める。$\lambda = 2$ から始めると、1.606、1.594 と収束し、$\hat{\lambda} = 1.594$ である。$\hat{P}(X = 0) = e^{-1.594} = 0.203$ だから、全会員の推定値は $120/(1 - 0.203) = 150.6$ 人で、記録に現れない会員は約31人である。

尤度方程式の右辺は $E(X \mid X \ge 1)$ です。「標本平均 = 切断された分布の平均」という形になっています。そのまま $\hat{\lambda} = 2.0$ とするのは、見えない0回の人を無視した誤りです。

試験のポイント：トランケーションの尤度は、分母の $P(T \in A)$ が母数を含むことが本質です。分母を忘れると、尤度方程式が変わります。

以上が、トランケーションについての説明です。

尤度の寄与を書き分ける

最後に、どの場面でも使える尤度の組み立て方を見ていきましょう。

ここまでの結論は単純です。個体ごとに「手に入った情報の確率」を書いて掛けるだけです。料理の材料を、1人ずつレシートに書き出すようなものです。

観測の種類	分かっていること	尤度の寄与
完全な観測	$T = t$	$f(t)$
右打ち切り	$T > c$	$S(c)$
左打ち切り	$T \le c$	$1 - S(c)$
区間打ち切り	$a < T \le b$	$S(a) - S(b)$
左トランケーション	$T = t$（ただし $T > v$）	$f(t)/S(v)$
左トランケーションと右打ち切り	$T > c$（ただし $T > v$）	$S(c)/S(v)$
右トランケーション	$T = t$（ただし $T \le v$）	$f(t)/\{1 - S(v)\}$
欠測（MAR）	観測された部分のみ	欠測部分を積分した密度

例題：いろいろな不完全データの尤度

寿命が指数分布に従う5個の個体について、次の情報を得た。① 時刻3で故障、② 時刻5でまだ動いている、③ 時刻1にはすでに故障していた、④ 時刻2から4の間に故障した、⑤ 時刻2より長く動いたものだけが調査の対象になる条件で、時刻6で故障。〔1〕尤度を書け。〔2〕対数尤度の微分を0とおく式を導け。〔3〕$\hat{\lambda}$ の値は約0.26である。平均寿命の推定値を述べよ。

【解答】

〔1〕それぞれの寄与を掛けて、次を得る。

$$L(\lambda) = \lambda e^{-3\lambda} \cdot e^{-5\lambda} \cdot (1 - e^{-\lambda}) \cdot (e^{-2\lambda} - e^{-4\lambda}) \cdot \lambda e^{-4\lambda}$$

〔2〕$e^{-2\lambda} - e^{-4\lambda} = e^{-2\lambda}(1 - e^{-2\lambda})$ だから、$\ell = 2\log\lambda - 14\lambda + \log(1 - e^{-\lambda}) + \log(1 - e^{-2\lambda})$ である。微分して、次を得る。

$$\frac{2}{\lambda} - 14 + \frac{e^{-\lambda}}{1 - e^{-\lambda}} + \frac{2e^{-2\lambda}}{1 - e^{-2\lambda}} = 0$$

〔3〕$1/0.26 \approx 3.85$ である。

逆に、尤度を正しく書かないと、結果は決まった向きに偏ります。平均10の指数分布を時刻8で打ち切るシミュレーションで確かめましょう。

【扱いを間違えたときの偏り】

• 打ち切りの個体を捨てる： 約3.5。長生きの個体を捨てるので短くなる
• 打ち切りの時刻を故障の時刻として平均する： 約5.5。やはり短くなる
• 正しい最尤推定（故障数で割る）： 約10.0
• 左トランケーションを無視する： 入り口の時刻を $(0, 10)$ の一様分布としたとき約14.2。長生きが選ばれるので長くなる

細かい公式を忘れてしまっても問題ありません。ここでは、「個体ごとに、手に入った情報の確率を書いて掛ける」という点だけ押さえておきましょう。

以上が、尤度の寄与の書き分けについての解説です。

今回のまとめ

それでは、今回の内容を振り返ります。

今回は不完全なデータの分析を、尤度から紹介しました。データが欠けたときに必要なのは、特別な公式ではありません。何が分かっていて、何が分からないのかを言葉にして、その確率を式にするだけです。欠け方を無視した計算は、決まった向きに偏ります。

今回のポイントを振り返りましょう。

振り返り

• 人数が分かれば打ち切り、人数すら分からなければトランケーションである。
• MARの欠測では、尤度が全員の分と観測できた人の分に分かれ、回帰で平均を補正できる。
• 無情報な打ち切りでは、尤度は $\prod f^{\delta}S^{1-\delta}$ になる。
• II型打ち切りの指数分布では、総試験時間について $2\lambda T \sim \chi^2_{2r}$ が成り立つ。
• トランケーションの尤度は、見えない分の確率で割った条件付き分布である。

キーワードのおさらい

• 欠測（欠損）: 値があるはずなのに手に入らないこと
• 打ち切り: 個体は観測されたが、値が入る範囲しか分からないこと
• トランケーション: ある範囲の値の個体が標本に現れないこと。切断とも言う

さて、今回のニュートン法やEMアルゴリズムは、尤度を数値で解く方法でした。しかし、尤度の式が書けても解けない問題や、統計量の分布が導けない問題もあります。そんなときは、コンピュータでくじを何度も引いて、分布そのものを作ってしまう手があります。

次回は、シミュレーションについて紹介します。

それでは、次の記事へ進んでいきましょう。