問題設定
ロジスティック回帰は「ある事象が起こる確率」を予測する線形統計モデルである。
今回、ベルヌーイ分布で記述することで、二値分類問題を解く統計モデルについて説明する。
その際の出力変数 $y$は
y\in{0,1}
また、説明変数（特徴量）は $x = (x_1, x_2, ・・・, x_n)$と仮定する。
線形結合とシグモイド関数
線形モデル説明変数 x に対して、まずは線形結合を示す線形予測子$z$は
z=β_0+β_1x_1+β_2x_2+⋯+β_nx_n
しかし、$z$をそのまま使うと $−∞∼+∞$までの値を取るため、確率$（0～1の範囲）$に変換するため、シグモイド関数を利用する。シグモイド関数の値域は$(0,1)$であり確率表現に好適である。
p=\frac{1}{1+e^{-z}}
なお、$y=1とy=0$のとき、各々下式の通り「事象が起きる確率」が得られる。
P(y=1∣x)=\frac{1}{1+e^{−(β_0+β_1x_1+⋯+β_nx_n)}}
P(y=0∣x)=1−p(x)
データの生成モデルを考える
さて、次に与えられた目的変数が最も起こり得るβを推定できる最尤推定法を実行する。そこで、まずは出力が2値を前提とした確率分布を記述しなければならない。ちょうどベルヌーイ分布というピッタリのモデルがある。
出力 $y\in0, 1$ として、入力$xに対してy=1$のときの確率を$p(x)$とすると、
P(y∣p)=p^y(x)(1−p(x))^{1−y}, y∈{0,1}
上記式は以下のように出力$y=1またはy=0$の確率を一つの式で表現できる(確率質量関数)。
P(y=1∣p)=p(x)
P(y=0∣p)=1−p(x)
複数サンプルの同時確率
最尤推定法で最大化したいのは目的変数全データをよく再現できる$β$である。そこで、複数サンプルの同時確率を考える。ここでデータ(行数)を$m$として、各サンプル$(x^{(i)}, y^{(i)})$が独立同分布だと仮定すると、全体の同時確率$L(β)$は積で表すことができるため、
L(β)=\prod_{i=1}^mP(x^{(i)})^{y^{(i)}}(1-p(x^{(i)}))^{(1-y^{(i)})}
これがロジスティック回帰における尤度関数となる。
直観的理解
上記式はデータごとに「観測された$y^{(i)}$が実際に出る確率」を掛け合わせている。
$y^{(i)}=1$のとき$p(x)^{(i)}$が掛けられ、$y^{(i)}=0$のとき$1-p(x)^{(i)}$が掛けられる。
対数尤度とその理由
対数尤度にする理由は、そのままの尤度関数(積の形)は計算上扱いにくいからである。そこで、対数を取ると
\ell(β)=\sum_{i=1}^my^{(i)}logP(x^{(i)})+(1-y^{(i)})log(1-p(x^{(i)}))
これが「対数尤度関数」で、最尤推定ではこれを最大化することになる。
勾配の準備
1データ・1特徴量に対する勾配
最大化の問題であるため、微分して勾配を求めるというよくある手法を使っていく。その前に、シグモイド関数の重要な性質をおさらいする。
p(x)=\frac{1}{1+e^{-x}}のとき、\frac{\partial}{\partial x} p(x)=p(x)(1-p(x)),
さて、本題である対数尤度関数に話を戻したい。まず簡単のために$i番目(i行目)$のサンプルを考える。
\ell^{(i)}(β)=y^{(i)}logP(x^{(i)})+(1-y^{(i)})log(1-p(x^{(i)})
$\ell^{(i)}(β)をβ_j$で偏微分すると
\frac{\partial \ell^{(i)}(\beta)}{\partial \beta_j}
= \frac{y^{(i)}}{p(x^{(i)})} \cdot \frac{\partial p(x^{(i)})}{\partial \beta_j}
+ \frac{1-y^{(i)}}{1-p(x^{(i)})} \cdot \frac{\partial \bigl(1-p(x^{(i)})\bigr)}{\partial \beta_j}
ここで、
\frac{\partial p(x^{(i)})}{\partial \beta_j}
= p(x^{(i)}) \bigl(1 - p(x^{(i)})\bigr) x^{(i)}_j
\frac{\partial \bigl(1-p(x^{(i)})\bigr)}{\partial \beta_j}
= -p(x^{(i)}) \bigl(1 - p(x^{(i)})\bigr) x^{(i)}_j
よって、1サンプル(1行データ)の$j$番目の特徴量についての勾配は以下の通りとなる。
\frac{\partial \ell^{(i)}(\beta)}{\partial \beta_j}
= \bigl(y^{(i)}-p(x^{(i)})\bigr)x_j^{(i)}
全データに対する勾配
$j$番目の特徴量の勾配を全サンプル(行方向)まで考慮すると次の通りである。
\frac{\partial \ell(\beta)}{\partial \beta_j}
= \sum_{i=1}^m\bigl(y^{(i)}-p(x^{(i)})\bigr)x_j^{(i)}
$j$番目の特徴量だけではなく、上の式について全特徴量まで拡張し、全サンプルの全特徴(すべての$m行n$列)を考慮した勾配は次のベクトルになる
∇\ell(\beta)
= \sum_{i=1}^m\bigl(y^{(i)}-p(x^{(i)})\bigr)x^{(i)}= X^\top (y - p)
- $X \in \mathbb{R}^{m \times n}$
- $y, p \in \mathbb{R}^m$
- 結果は $n$ 次元ベクトル(特徴量の数)
念のためだが、勾配の形式は $\beta$ に関する偏微分をまとめたベクトルである：
\nabla \ell(\beta) =
\begin{bmatrix}
\frac{\partial \ell}{\partial \beta_1} \\
\frac{\partial \ell}{\partial \beta_2} \\
\vdots \\
\frac{\partial \ell}{\partial \beta_n}
\end{bmatrix}
\in \mathbb{R}^n
直感的解釈
$y^{(i)}−p(x^{(i)})$ は「実際のラベル − 予測確率」であり、つまり「誤差」と理解できる。それに入力 $x^{(i)}$を掛けて勾配を作っている。結果的に$y-p(x)$の項が出てくるのは、典型的な誤差関数を連想させる。元々は最尤推定法で$\beta$を求めていたはずだが、このような式になる理由は指数型分布族の一般的な性質に由来する。
多変数関数の極値判定法 (ヘッセ行列)
多変数関数の極値判定法は「勾配=0」であれば、極大または極小となる。極大または極小を判定するための方法は、二階偏微分作用を持つヘッセ行列$H$において、$H<0$であれば極大、$H>0$であれば極小となる。
ヘッセ行列とは次の通りである。
n変数関数p(x_1,x_2,・・・,x_n)に対して、jk成分が\frac{\partial^2 p(x)}{\partial x_j\partial x_k}であるようなn×n行列
$j行目,k列目$を代表例として、$(j,k)$成分を偏微分形式で表示すると
\frac{\partial^2 \ell}{\partial \beta_j\partial \beta_k}= -\sum_{i=1}^m\bigl(p_i(1-p_i)\bigr)x_{ij}x_{ik}
行列形式では
H(\beta) = ∇^2\ell(\beta) = -X^TWX
ここで$W は対角で各対角要素は p_i(1−p_i)∈(0,0.25]$。以上の行列形式で定義されたヘッセ行列$H(\beta)$の正負は$W$に依存し、$Wは(0. 0.25]$の値であるため常に負。よって極大解と確認できる。詳細な考え方は以下の通り。
任意のベクトル $( z \in \mathbb{R}^n )$ に対して考える。
$$
z^\top (X^\top W X) z = (X z)^\top W (X z)
$$
ここで $u = Xz \in \mathbb{R}^m$ とすると、内積は
$$
(X z)^\top W (X z) = u^\top W u = \sum_{i=1}^m w_i u_i^2
$$
各 $w_i = p_i (1 - p_i) \geq 0$ なので
$$
\sum_{i=1}^m w_i u_i^2 \geq 0
$$
したがって、
$$
X^\top W X \succeq 0 \quad （半正定値）
$$
これが「 $W \geq 0$ だから $X^\top W X \geq 0$ 」と表現される理由である。
ニュートン法による勾配0になる点の探索
ニュートン法の更新式は次の通りである。
x^{(t+1)} = x^{t} - \frac{f(x^{(t)})}{f'(x^{(t)})}
ニュートン法の導出: ニュートン法の導出は以下の通りである。下の図のように$x_1, x_2, f(x_1)$で作られる直角三角形に注目する。$⊿x, ⊿y$はそれぞれ直角三角形の辺の長さを表す。
$f'(x_1)$ は $f(x_1)$ における接線なので、傾きの定義式の関係より、$x_2$ を求めることができる。
$$
f'(x_1) = \frac{\Delta y}{\Delta x}
$$
$$
= \frac{f(x_1)}{x_1 - x_2}
$$
$$
\Longleftrightarrow \quad x_1 - x_2 = \frac{f(x_1)}{f'(x_1)}
$$
$$
\therefore \quad x_2 = x_1 - \frac{f(x_1)}{f'(x_1)}
$$
対数尤度関数の勾配(下式)0になる条件の探索は上記ニュートン法などで行う。勾配は
∇\ell(\beta)
= \sum_{i=1}^m\bigl(y^{(i)}-p(x^{(i)})\bigr)x^{(i)}
であり、ニュートン法に当てはめると$f(x)$が勾配$∇\ell(\beta)$(1階偏微分)となるため、ニュートン法の分母の導関数は二階偏微分となる。この二階偏微分は行列形式であればヘッセ行列$H$である。また、行列形式の場合、分母に相当する計算は逆行列をかけることになる。よって、
\beta^{(t+1)} = \beta^{(t)} - H(\beta^{(t)})^{-1}∇\ell(\beta^{(t)})=β^{(t)}+(X^⊤WX)^{−1}X^⊤(y−p)
という式で$\beta$を数値解析的に求めることになる。なお、$\beta$は特徴量の次元nを持つベクトルである。