ジョン・ラチン(2020)『医薬データのための統計解析』 問題6.4 解答例

公開日: 更新日:

【2022年11月4週】 【A000】生物統計学 【A092】ロジスティック回帰分析

この記事をシェアする
  • B!
サムネイル画像

本稿は、ジョン・ラチン(2020)『医薬データのための統計解析』の「問題6.4」の自作解答例です。積二項尤度にもとづく対数リスク比のロジスティック回帰モデルに関する問題です。

なお、閲覧にあたっては、以下の点にご注意ください。

  • スマートフォンやタブレット端末でご覧の際、数式が見切れている場合は、横にスクロールすることができます。
  • 曝露(発症)状況を表す右下の添え字は、「0」である場合($n_0,\pi_0$ など)や「2」である場合($n_2,\pi_2$ など)がありますが、どちらも「非曝露群(コントロール群)」を表しています。
  • 漸近的な性質を用いる際は、①中心極限定理が成り立つ、②漸近分散を推定する際に、母数をその一致推定量で置き換えることができるということが成り立つと仮定しています。
  • デルタ法を用いる際、剰余項(2次の項)が漸近的に無視できる($0$に確率収束する)と仮定しています。
  • 上述の参考書では、標準正規分布の上側 $100\alpha\%$ 点を $Z_{1-\alpha}$ と表記していますが、本サイトでは、$Z_\alpha$ としています。そのため、参考書に載っている式の形式と異なる部分があります。
  • 著作権の関係上、問題文は、掲載しておりません。上述の参考書をお持ちの方は、お手元にご用意してご覧ください。
  • この解答例は、筆者が自作したものであり、公式なものではありません。あくまでも参考としてご覧いただければ幸いです。

問題6.4.1:スコア関数①

積二項尤度の基本式より、 \begin{align} L \left(\pi_1,\pi_2\right)={}_{n_1}C_a \cdot \pi_1^a \left(1-\pi_1\right)^{n_1-a} \cdot {}_{n_2}C_b \cdot \pi_2^b \left(1-\pi_2\right)^{n_2-b} \end{align} 対数尤度関数の定義 $l \left(\theta,\boldsymbol{x}\right)=\log{L \left(\theta,\boldsymbol{x}\right)}$ より、 \begin{align} l \left(\pi_1,\pi_2\right)&=\log{ \left[{}_{n_1}C_a \cdot \pi_1^a \left(1-\pi_1\right)^{n_1-a} \cdot {}_{n_2}C_b \cdot \pi_2^b \left(1-\pi_2\right)^{n_2-b}\right]}\\ &=\log{{}_{n_1}C_a}+a\log{\pi_1}+ \left(n_1-a\right)\log{ \left(1-\pi_1\right)}+\log{{}_{n_2}C_b}+b\log{\pi_2}+ \left(n_2-b\right)\log{ \left(1-\pi_2\right)} \end{align} $\log{{}_{n_1}C_a}+\log{{}_{n_2}C_b}=C$ とすると、各セルどうしの関係、対数リスクモデルの仮定より、 \begin{gather} n_1-a=c n_2-b=d a+b=m_1 \log{\pi_1}=\alpha+\beta \log{\pi_2}=\alpha \end{gather} よって、 \begin{align} l \left(\pi_1,\pi_2\right)&=a\log{\pi_1}+c\log{ \left(1-\pi_1\right)}+b\log{\pi_2}+d\log{ \left(1-\pi_2\right)}+C\\ &=a\log{ \left(e^{\alpha+\beta}\right)}+c\log{ \left(1-e^{\alpha+\beta}\right)}+b\log{ \left(e^\alpha\right)}+d\log{ \left(1-e^\alpha\right)}+C\\ &=a \left(\alpha+\beta\right)+c\log{ \left(1-e^{\alpha+\beta}\right)}+b\alpha+d\log{ \left(1-e^\alpha\right)}+C\\ &= \left(a+b\right)\alpha+a\beta+c\log{ \left(1-e^{\alpha+\beta}\right)}+d\log{ \left(1-e^\alpha\right)}+C\\ &=m_1\alpha+a\beta+c\log{ \left(1-e^{\alpha+\beta}\right)}+d\log{ \left(1-e^\alpha\right)}+C \end{align} パラメータ $\alpha$ のスコア関数は、この対数尤度関数を $\alpha$ で偏微分して、 \begin{align} U_\alpha \left(\boldsymbol{\theta}\right)&=m_1-c \cdot \frac{e^{\alpha+\beta}}{1-e^{\alpha+\beta}}-d \cdot \frac{e^\alpha}{1-e^\alpha}\\ &=m_1-\frac{c\pi_1}{1-\pi_1}-\frac{d\pi_2}{1-\pi_2} \end{align} $\blacksquare$

問題6.4.2:スコア関数②

同様に、パラメータ $\beta$ のスコア関数は、この対数尤度関数を $\beta$ で偏微分して、 \begin{align} U_\beta \left(\boldsymbol{\theta}\right)&=a-c \cdot \frac{e^{\alpha+\beta}}{1-e^{\alpha+\beta}}\\ &=a-\frac{c\pi_1}{1-\pi_1} \end{align} $\blacksquare$

問題6.4.3:最尤推定量

最尤推定量の定義 $U \left(\hat{\theta}\right)=0$ より、 \begin{gather} U_\alpha \left(\hat{\alpha},\hat{\beta}\right)=m_1-\frac{c\pi_1}{1-\pi_1}-\frac{d\pi_2}{1-\pi_2}=0\tag{1}\\ U_\beta \left(\hat{\alpha},\hat{\beta}\right)=a-\frac{c\pi_1}{1-\pi_1}=0\tag{2} \end{gather} この連立方程式を解くと、 \begin{gather} m_1-a-\frac{de^{\hat{\alpha}}}{1-e^{\hat{\alpha}}}=0\\ \frac{e^{\hat{\alpha}}}{1-e^{\hat{\alpha}}}=\frac{b}{d}\\ e^{\hat{\alpha}}=\frac{b}{d}-\frac{b}{d}e^{\hat{\alpha}}\\ \left(1+\frac{b}{d}\right)e^{\hat{\alpha}}=\frac{b}{d}\\ \left(\frac{d+b}{d}\right)e^{\hat{\alpha}}=\frac{b}{d}\\ e^{\hat{\alpha}}=\frac{b}{d} \cdot \frac{d}{d+b}\\ e^{\hat{\alpha}}=\frac{d}{d+b}=\frac{b}{n_2}\\ \hat{\alpha}=\log{\frac{b}{n_2}} \end{gather} \begin{gather} \frac{e^{\hat{\alpha}+\hat{\beta}}}{1-e^{\hat{\alpha}+\hat{\beta}}}=\frac{a}{c}\\ e^{\hat{\alpha}+\hat{\beta}}=\frac{a}{c}-\frac{a}{c}e^{\hat{\alpha}+\hat{\beta}}\\ \left(1+\frac{a}{c}\right)e^{\hat{\alpha}+\hat{\beta}}=\frac{a}{c}\\ e^{\hat{\alpha}+\hat{\beta}}=\frac{a}{c} \cdot \frac{c}{a+c}=\frac{a}{n_1}\\ e^{\hat{\beta}}\frac{b}{n_2}=\frac{a}{n_1}\\ e^{\hat{\beta}}=\frac{an_2}{bn_1}\\ \hat{\beta}=\log{\frac{an_2}{bn_1}}=\log{\widehat{RR}} \end{gather} $\blacksquare$

問題6.4.4:期待情報行列

ヘッセ行列の成分の定義と $c=n_1 \left(1-\pi_1\right),d=n_2 \left(1-\pi_2\right)$ より、 \begin{align} H_\alpha \left(\boldsymbol{\theta}\right)&=\frac{\partial^2}{\partial\alpha^2}l \left(\boldsymbol{\theta}\right)\\ &=-c \cdot \frac{e^{\alpha+\beta}}{ \left(1-e^{\alpha+\beta}\right)^2}-d \cdot \frac{e^\alpha}{ \left(1-e^\alpha\right)^2}\\ &=-c\frac{\pi_1}{ \left(1-\pi_1\right)^2}-d\frac{\pi_2}{ \left(1-\pi_2\right)^2}\\ &=-n_1 \cdot \frac{\pi_1}{1-\pi_1}-n_2 \cdot \frac{\pi_2}{1-\pi_2} \end{align} \begin{align} H_\beta \left(\boldsymbol{\theta}\right)&=\frac{\partial^2}{\partial\beta^2}l \left(\boldsymbol{\theta}\right)\\ &=-c \cdot \frac{e^{\alpha+\beta}}{ \left(1-e^{\alpha+\beta}\right)^2}\\ &=-c\frac{\pi_1}{ \left(1-\pi_1\right)^2}\\ &=-n_1 \cdot \frac{\pi_1}{1-\pi_1} \end{align} \begin{align} H_{\alpha\beta} \left(\boldsymbol{\theta}\right)&=\frac{\partial^2}{\partial\alpha\partial\beta}l \left(\boldsymbol{\theta}\right)\\ &=c \cdot \frac{e^{\alpha+\beta}}{ \left(1-e^{\alpha+\beta}\right)^2}\\ &=-c\frac{\pi_1}{ \left(1-\pi_1\right)^2}\\ &=-n_1 \cdot \frac{\pi_1}{1-\pi_1} \end{align} 観測情報行列の定義 $\boldsymbol{i} \left(\boldsymbol{\theta}\right)=-\boldsymbol{H} \left(\boldsymbol{\theta}\right)$ と二項分布における期待情報行列との関係 $\boldsymbol{i} \left(\boldsymbol{\theta}\right)=\boldsymbol{I} \left(\boldsymbol{\theta}\right)$ から、$\psi_i=\frac{\pi_i}{1-\pi_i}$ とすると、 \begin{align} \boldsymbol{I} \left(\boldsymbol{\theta}\right)&=\boldsymbol{i} \left(\boldsymbol{\theta}\right)\\ &=- \left[\begin{matrix}H_\alpha \left(\boldsymbol{\theta}\right)&H_{\alpha\beta} \left(\boldsymbol{\theta}\right)\\H_{\alpha\beta} \left(\boldsymbol{\theta}\right)&H_\beta \left(\boldsymbol{\theta}\right)\\\end{matrix}\right]\\ &= \left[\begin{matrix}n_1\psi_1+n_2\psi_2&n_1\psi_1\\n_1\psi_1&n_1\psi_1\\\end{matrix}\right] \end{align} $\blacksquare$

問題6.4.5:漸近分散・共分散行列

期待情報行列の成分は、 \begin{align} \boldsymbol{I} \left(\boldsymbol{\theta}\right)= \left[\begin{matrix}\boldsymbol{I}_\alpha \left(\boldsymbol{\theta}\right)&\boldsymbol{I}_{\alpha\beta} \left(\boldsymbol{\theta}\right)\\\boldsymbol{I}_{\alpha\beta} \left(\boldsymbol{\theta}\right)^T&\boldsymbol{I}_\beta \left(\boldsymbol{\theta}\right)\\\end{matrix}\right] \end{align} このとき、推定量の共分散行列は、$\boldsymbol{V} \left(\hat{\boldsymbol{\theta}}\right)={\boldsymbol{I} \left(\boldsymbol{\theta}\right)}^{-1}$ となる。期待情報行列の逆行列は、 \begin{align} {\boldsymbol{I} \left(\boldsymbol{\theta}\right)}^{-1}=\frac{\boldsymbol{1}}{ \left|\boldsymbol{I} \left(\boldsymbol{\theta}\right)\right|} \left[\begin{matrix}\boldsymbol{I}_\beta \left(\boldsymbol{\theta}\right)&-\boldsymbol{I}_{\alpha\beta} \left(\boldsymbol{\theta}\right)^T\\-\boldsymbol{I}_{\alpha\beta} \left(\boldsymbol{\theta}\right)&\boldsymbol{I}_\alpha \left(\boldsymbol{\theta}\right)\\\end{matrix}\right] \end{align} 2×2行列の行列式の公式 $ \left|\boldsymbol{A}\right|=a_{11}a_{22}-a_{12}a_{21}$ より、 \begin{align} \left|\boldsymbol{I} \left(\boldsymbol{\theta}\right)\right|&= \left(n_1\psi_1+n_2\psi_2\right)n_1\psi_1- \left\{n_1\psi_1\right\}^2\\ &= \left(n_1\psi_1+n_2\psi_2-n_1\psi_1\right)n_1\psi_1\\ &=n_1\psi_1n_2\psi_2 \end{align} したがって、 \begin{align} \boldsymbol{V} \left(\hat{\alpha},\hat{\beta}\right)&={\boldsymbol{I} \left(\boldsymbol{\theta}\right)}^{-1}\\ &=\frac{1}{n_1\psi_1n_2\psi_2} \left[\begin{matrix}n_1\psi_1&-n_1\psi_1\\-n_1\psi_1&n_1\psi_1+n_2\psi_2\\\end{matrix}\right]\\ &= \left[\begin{matrix}\frac{1}{n_2\psi_2}&-\frac{1}{n_2\psi_2}\\-\frac{1}{n_2\psi_2}&\frac{n_1\psi_1+n_2\psi_2}{n_1\psi_1n_2\psi_2}\\\end{matrix}\right] \end{align} したがって、推定量の共分散行列の成分は、 \begin{align} V \left(\hat{\alpha}\right)&= \left[{\boldsymbol{I} \left(\boldsymbol{\theta}\right)}^{-1}\right]_\boldsymbol{\alpha}\\ &=\frac{1-\pi_2}{n_2\pi_2} \end{align} \begin{align} V \left(\hat{\beta}\right)&= \left[{\boldsymbol{I} \left(\boldsymbol{\theta}\right)}^{-1}\right]_\boldsymbol{\beta}\\ &=\frac{1}{n_1\psi_1}+\frac{1}{n_2\psi_2}\\ &=\frac{1-\pi_1}{n_1\pi_1}+\frac{1-\pi_2}{n_2\pi_2} \end{align} \begin{align} \mathrm{Cov} \left(\hat{\alpha},\hat{\beta}\right)&= \left[{\boldsymbol{I} \left(\boldsymbol{\theta}\right)}^{-1}\right]_{\boldsymbol{\alpha},\boldsymbol{\beta}}\\ &=-\frac{1-\pi_2}{n_2\pi_2} \end{align} $\blacksquare$

問題6.4.6:対数リスク比の推定分散

この一致推定量は、${\hat{\pi}}_1=p_1=e^{\hat{\alpha}+\hat{\beta}}=\frac{an_2}{bn_1},{\hat{\pi}}_2=p_2=e^{\hat{\alpha}}=\frac{b}{n_2}$ を代入して、 \begin{align} \hat{V} \left(\hat{\beta}\right)&=\hat{V} \left(\log{\mathrm{\widehat{RR}}}\right)\\ &=\frac{1-p_1}{n_1p_1}+\frac{1-p_2}{n_2p_2}\\ &=\frac{1-p_1}{a}+\frac{1-p_2}{b}\\ &=\frac{1}{a}-\frac{p_1}{a}+\frac{1}{b}-\frac{p_2}{b} \end{align} $a=n_1p_1,b=n_2p_2$ なので、 \begin{align} \hat{V} \left(\log{\mathrm{\widehat{RR}}}\right)=\frac{1}{a}-\frac{1}{n_1}+\frac{1}{b}-\frac{1}{n_2} \end{align} $\blacksquare$

問題6.4.7:帰無仮説のもとでの共通の母比率の最尤推定量

対数尤度の式より、帰無仮説 $\beta=0$ における対数尤度は、 \begin{align} l \left(\pi_1,\pi_2\right)&=m_1\alpha+a\beta+c\log{ \left(1-e^{\alpha+\beta}\right)}+d\log{ \left(1-e^\alpha\right)}\\ l \left(\hat{\alpha},0\right)&=\log{L \left(\hat{\alpha},0\right)}\\ &=m_1\alpha+a \cdot 0+c\log{ \left(1-e^{\alpha+0}\right)}+d\log{ \left(1-e^\alpha\right)}\\ &=m_1\alpha+c\log{ \left(1-e^\alpha\right)}+d\log{ \left(1-e^\alpha\right)}=m_1\alpha+m_2\log{ \left(1-e^\alpha\right)} \end{align} パラメータ $\alpha$ のスコア関数は、この対数尤度関数を $\alpha$ で偏微分すると、 \begin{align} U_\alpha \left(\boldsymbol{\theta}\right)=m_1-m_2\frac{e^\alpha}{1-e^\alpha} \end{align} 尤度方程式 $U_\alpha \left(\boldsymbol{\theta}\right)=0$ を解くことによって得られる $\alpha$ の最尤推定量は、 \begin{gather} m_1-m_2\frac{e^{{\hat{\alpha}}_0}}{1-e^{{\hat{\alpha}}_0}}=0\\ \frac{e^{{\hat{\alpha}}_0}}{1-e^{{\hat{\alpha}}_0}}=\frac{m_1}{m_2}\\ \left(1+\frac{m_1}{m_2}\right)e^{{\hat{\alpha}}_0}=\frac{m_1}{m_2}\\ e^{{\hat{\alpha}}_0}=\frac{m_1}{m_2} \cdot \frac{m_2}{m_1+m_2}\\ e^{{\hat{\alpha}}_0}=\frac{m_1}{N}\\ {\hat{\alpha}}_0=\log{\frac{m_1}{N}} \end{gather} したがって、共通の発症確率 $\pi_1=\pi_2=\pi$ の最尤推定量は、 \begin{align} \hat{\pi}=\frac{m_1}{N} \end{align} $\blacksquare$

問題6.4.8:帰無仮説のもとでのスコア方程式

$\beta$ に関するスコア関数 $U_\beta \left(\boldsymbol{\theta}\right)=a-c \cdot \frac{e^{\alpha+\beta}}{1-e^{\alpha+\beta}}$ は、帰無仮説 $\beta=0$ においては、 \begin{align} U_\beta \left(\boldsymbol{\theta}\right)&=a-n_1 \cdot \frac{e^{{\hat{\alpha}}_0}}{1+e^{{\hat{\alpha}}_0}}\\ &=a-c \cdot \frac{\frac{m_1}{N}}{1-\frac{m_1}{N}}\\ &=a-c \cdot \frac{m_1}{m_2} \end{align} $\blacksquare$

問題6.4.9:帰無仮説のもとでの推定情報行列

また、帰無仮説のもとで $\beta$ に対応する推定情報量の逆行列の成分は、 \begin{align} V \left(\hat{\beta}\right)&= \left[{\boldsymbol{I} \left(\boldsymbol{\theta}\right)}^{-1}\right]_\boldsymbol{\beta}\\ &=\frac{1-\hat{\pi}}{n_1\hat{\pi}}+\frac{1-\hat{\pi}}{n_2\hat{\pi}}\\ &=\frac{1-\hat{\pi}}{\hat{\pi}} \left(\frac{1}{n_1}+\frac{1}{n_2}\right)\\ &=\frac{1-\frac{m_1}{N}}{\frac{m_1}{N}} \left(\frac{1}{n_1}+\frac{1}{n_2}\right)\\ &=\frac{m_2}{m_1} \left(\frac{1}{n_1}+\frac{1}{n_2}\right)\\ &=\frac{m_2}{m_1} \cdot \frac{n_1+n_2}{n_1n_2}\\ &=\frac{m_2}{m_1} \cdot \frac{N}{n_1n_2} \end{align} $\blacksquare$

問題6.4.10:有効スコア検定

結果として有効スコア検定 $\chi^2=\frac{ \left\{U_\beta \left(\boldsymbol{\theta}\right)\right\}^2}{\hat{V} \left(\hat{\beta}\right)}$ は、 \begin{align} \chi^2= \left(a-c \cdot \frac{m_1}{m_2}\right)^2 \left(\frac{m_1}{m_2} \cdot \frac{n_1n_2}{N}\right) \end{align} $\blacksquare$

問題6.4.12:スコアにもとづく共通対数リスク比の推定値

スコアにもとづく推定値の公式 ${\hat{\theta}}_0=\frac{U \left(\theta_0\right)}{I \left(\theta_0\right)}$ より、 \begin{align} \hat{\beta}= \left(a-c \cdot \frac{m_1}{m_2}\right) \left(\frac{m_1}{m_2} \cdot \frac{n_1n_2}{N}\right) \end{align} また、スコアにもとづく推定値の分散の公式 $V \left({\hat{\theta}}_0\right)=\frac{1}{I \left(\theta_0\right)}$ より、 \begin{align} V \left(\hat{\beta}\right)=\frac{m_1}{m_2} \cdot \frac{n_1n_2}{N} \end{align} $\blacksquare$

参考文献

  • ジョン・ラチン 著, 宮岡 悦良 監訳, 遠藤 輝, 黒沢 健, 下川 朝有, 寒水 孝司 訳. 医薬データのための統計解析. 共立出版, 2020, p.297-298

関連記事

自己紹介

自分の写真

yama

大学時代に読書の面白さに気づいて以来、読書や勉強を通じて、興味をもったことや新しいことを学ぶことが生きる原動力。そんな人間が、その時々に学んだことを備忘録兼人生の軌跡として記録しているブログです。

このブログを検索

ブログ アーカイブ

QooQ