一般ロジスティック回帰モデル

本稿では、ロジスティック・モデルを共変量が複数の場合に拡張した一般ロジスティック回帰モデルを紹介しています。一般ロジスティック回帰モデルは、関心のあるイベントの発生確率の予測や交絡因子の影響の除去などに使用されます。

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

一般ロジスティック回帰モデル

定義

マッチングなしのコホート研究（積二項モデル）において、発症確率が、$p$ 変量ベクトルと対応する $p+1$ 個の成分からなるパラメータベクトル $$\begin{array}{r}\mathit{X}=\left\{\begin{array}{c}{X}_1\\ X_2\\ ⋮\\ X_p\end{array}\right\}\mathit{\theta }=\left\{\begin{array}{c}\alpha \\ \mathit{\beta }\end{array}\right\}\mathit{\beta }=\left\{\begin{array}{c}{\beta }_1\\ {\beta }_2\\ ⋮\\ {\beta }_p\end{array}\right\}\end{array}$$ の線形関数として表現できると仮定する。 このとき、$i$ 番目の被験者に対し、その被験者の共変量ベクトル $$\begin{array}{r}{\mathit{x}}_{\mathit{i}}=\left\{\begin{array}{c}{x}_{i1}\\ x_{i2}\\ ⋮\\ x_{ip}\end{array}\right\}\end{array}$$ の関数として その確率を $$\begin{array}{r}{\pi }_i({\mathit{x}}_{\mathit{i}},\mathit{\theta })=P(y_i=1|{\mathit{x}}_{\mathit{i}},\mathit{\theta })\end{array}$$ と書く。 このとき、これらの $p$ 個の共変量によって表現された確率のロジットは、パラメータベクトル $\mathit{\theta }$ の線形関数 $$\begin{array}{c}\log \frac{{\pi }_i}{1-{\pi }_i}=\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }=\alpha +\sum _{k=1}^p{x}_{ik}\cdot {\beta }_k\\ {\pi }_i={\pi }_i({\mathit{x}}_{\mathit{i}},\mathit{\theta })\end{array}$$ のように表現することができるものとする 統計モデルを一般ロジスティック回帰モデル
General Logistic Regression Models と呼ぶ。

$i$ 番目の被験者の発症確率は、ロジスティック関数を用いて $$\begin{array}{c}{\pi }_i=\frac{e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}{1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}1-{\pi }_i=\frac{1}{1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}\end{array}$$ で与えられる。

パラメータの最尤推定

観測値はベルヌーイ変数 $\left\{y_i\right\}$ なので、確率 $\mathit{\pi }=\{{\pi }_1,{\pi }_2,\cdots ,{\pi }_N\}$ に関する尤度は $$\begin{array}{r}L\left(\mathit{\pi }\right)=\prod _{i=1}^N{\pi }_i^{y_i}{(1-{\pi }_i)}^{1-y_i}\end{array}$$ 共変量とパラメータのロジスティック関数として確率を表すと、 $$\begin{array}{r}L\left(\mathit{\theta }\right)=\prod _{i=1}^N{\left(\frac{e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}{1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}\right)}^{y_i}{\left(\frac{1}{1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}\right)}^{1-y_i}\end{array}$$ 対数尤度関数 $l(\theta ,\mathit{x})=\log L(\theta ,\mathit{x})$ は、 $$\begin{array}{c}l\left(\mathit{\pi }\right)=\sum _{i=1}^N\{y_i\log {\pi }_i+(1-y_i)\log (1-{\pi }_i)\}\\ l\left(\mathit{\theta }\right)=\sum _{i=1}^N\{y_i(\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta })-\log (1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }})\}\end{array}$$ ここで、パラメータに対するスコアベクトルを $$\begin{array}{r}U\left(\mathit{\theta }\right)=\left\{\begin{array}{c}{U}_{\alpha }\left(\mathit{\theta }\right)\\ U_{{\beta }_1}\left(\mathit{\theta }\right)\\ ⋮\\ U_{{\beta }_p}\left(\mathit{\theta }\right)\end{array}\right\}\end{array}$$ とおく。 〔a〕パラメータ $\alpha$ に関するスコア関数 $U\left(\theta \right)=\frac{\partial }{\partial \theta }l\left(\theta \right)$ は、 $$\begin{array}{rl}{U}_{\alpha }\left(\mathit{\theta }\right)& =\sum _{i=1}^N{y}_i-\sum _{i=1}^N\frac{e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}{1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}\\ & =\sum _{i=1}^N{y}_i-\sum _{i=1}^N{\pi }_i\\ & =m_1-\sum _{i=1}^N{\pi }_i\\ & =\sum _{i=1}^N\{y_i-E\left(y_i\right|{\mathit{x}}_{\mathit{i}})\}\\ & =m_1-E\left(m_1\right|\mathit{\theta })\end{array}$$ 〔b〕パラメータ ${\beta }_k$ に関するスコア関数 $U\left(\theta \right)=\frac{\partial }{\partial \theta }l\left(\theta \right)$ は、 $$\begin{array}{rl}{U}_{{\beta }_k}\left(\mathit{\theta }\right)& =\sum _{i=1}^N{y}_i\cdot x_{ik}-\sum _{i=1}^N\frac{e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}{1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}\cdot x_{ik}\\ & =\sum _{i=1}^N(y_i-{\pi }_i)x_{ik}\end{array}$$ ヘッセ行列の成分の定義より、 $$\begin{array}{rl}{H}_{\alpha }\left(\mathit{\theta }\right)& =\frac{\partial }{\partial \alpha }{U}_{\alpha }\left(\mathit{\theta }\right)\\ & =-\sum _{i=1}^N\frac{e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}{{(1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }})}^2}\\ & =-\sum _{i=1}^N{\pi }_i(1-{\pi }_i)\end{array}$$ $$\begin{array}{rl}{H}_{{\beta }_k}\left(\mathit{\theta }\right)& =\frac{\partial }{\partial {\beta }_k}{U}_{{\beta }_k}\left(\mathit{\theta }\right)\\ & =-\sum _{i=1}^N\frac{e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}{{(1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }})}^2}\cdot x_{ik}\cdot x_{ik}\\ & =-\sum _{i=1}^N{\pi }_i(1-{\pi }_i)\cdot x_{ik}^2\end{array}$$ $$\begin{array}{rl}{H}_{{\beta }_k{\beta }_l}\left(\mathit{\theta }\right)& =\frac{\partial }{\partial {\beta }_l}{U}_{{\beta }_k}\left(\mathit{\theta }\right)\\ & =-\sum _{i=1}^N\frac{e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}{{(1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }})}^2}\cdot x_{ik}\cdot x_{il}\\ & =-\sum _{i=1}^N{\pi }_i(1-{\pi }_i)\cdot x_{ik}{x}_{il}\end{array}$$ $$\begin{array}{rl}{H}_{\alpha {\beta }_k}\left(\mathit{\theta }\right)& =\frac{\partial }{\partial \alpha }{U}_{{\beta }_k}\left(\mathit{\theta }\right)\\ & =-\sum _{i=1}^N\frac{e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }}}{{(1+e^{\alpha +{\mathit{x}}_{\mathit{i}}^{\mathit{T}}\mathit{\beta }})}^2}\cdot x_{ik}\\ & =-\sum _{i=1}^N{\pi }_i(1-{\pi }_i)\cdot x_{ik}\end{array}$$ 二項分布における観測情報行列の定義式 $\mathit{I}\left(\mathit{\theta }\right)=\mathit{i}\left(\mathit{\theta }\right)=-\mathit{H}\left(\mathit{\theta }\right)$ より、 $$\begin{array}{c}{I}_{\alpha }\left(\mathit{\theta }\right)=i_{\alpha }\left(\mathit{\theta }\right)=\sum _{i=1}^N{\pi }_i(1-{\pi }_i)\\ I_{{\beta }_k}\left(\mathit{\theta }\right)=i_{{\beta }_k}\left(\mathit{\theta }\right)=\sum _{i=1}^N{\pi }_i(1-{\pi }_i)\cdot x_{ik}^2\\ I_{{\beta }_k{\beta }_l}\left(\mathit{\theta }\right)=i_{{\beta }_k{\beta }_l}\left(\mathit{\theta }\right)=\sum _{i=1}^N{\pi }_i(1-{\pi }_i)\cdot x_{ik}{x}_{il}\\ I_{\alpha {\beta }_k}\left(\mathit{\theta }\right)=i_{\alpha {\beta }_k}\left(\mathit{\theta }\right)=\sum _{i=1}^N{\pi }_i(1-{\pi }_i)\cdot x_{ik}\end{array}$$

各パラメータの最尤推定量は、方程式 $U\left(\mathit{\theta }\right)=0$ の解、つまり、$p+1$ 個のスコア方程式が0 となる解によって得られる。閉形式の解は存在しないため、解を得るためにはニュートン・ラフソンアルゴリズムなどの反復計算によって求める必要がある。

参考文献

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

関連記事

• リサーチ・クエスチョンの整理と評価
• 交絡
• コホート研究
• ロジスティック回帰分析
• 横断研究・コホート研究【有病率・発生割合】（マッチングなし・層化なし）
• ロジスティック・モデル
• ロジスティック・モデルによる対数オッズ比の最尤推定
• ロジスティック・モデルによる対数リスク比の最尤推定
• ロジスティック・モデルによる共通対数オッズ比の最尤推定
• 条件付きロジスティック・モデル