アラフォーからの機械学習

アラフォーです。機械学習や統計学の解説記事を書いていきます

ロジスティック回帰をはじめから

ロジスティック回帰について

ロジスティック回帰は機械学習の分類問題において、おそらく基本的な手法に位置するものであると思います。この手法は線形回帰モデルをその目的変数について、カテゴリカルな量に対応できるよう拡張させたものと考えることができます。カテゴリカルな量としてはたとえば、合格と不合格、生と死などそれらが連続的ではないものを指します。たとえば、合格と不合格かを判定するロジックを考えてみます。直感的には、それを説明する変数（たとえば勉強時間や難易度など）を考慮して合格する確率はxx%だというようにすると良さそうです。予測した確率に対してしきい値を設定して、たとえば50%とすればそれを境に合格と不合格を判定することができます。ということを考えると、目的変数としてはあるカテゴリになる確率を設定すれば良さそうです。つまり、ロジスティック回帰は確率の回帰を行うモデルです。しかしながら、線形回帰モデルで出力される値の範囲は-∞から∞です。一方で、確率は0から1です。そのため、線形結合した値について0から1の範囲に収まるように変換が必要になります。この変換にロジスティック分布を用いたものをロジスティック回帰モデル（ロジットモデル）と言います。

ロジスティック回帰モデルを導出する

予測の対象が合格か不合格かなど2値を取るような場合、確率モデルとしてベルヌーイ分布を考えることができます。たとえば合格する確率をpとすれば、以下のように書けます。

\displaystyle{ P(y = 合格) = p }

一般化するために、合格を1で表し、複数回分の結果があるとすれば、

\displaystyle{ P(y_i = 1) = p_i }

となり、それぞれが独立であるとすれば、

\displaystyle{ y_i \sim Ber(p_i) }

と表現できます。合格する確率がp_iのベルヌーイ分布に従います。では、勉強時間などの要因と考えられる説明変数とその回帰係数の線形結合\boldsymbol{x}^\mathrm{T}_i\boldsymbol{\beta}の結果から確率p_iを見たいとします。つまり、線形結合\boldsymbol{x}^\mathrm{T}_i\boldsymbol{\beta}が増加するにつれて、このp_iも増加するように線形結合の値を変換することを考えます。なお、線形結合の範囲は[-\infty, \infty]ですが確率の範囲は[0, 1]です。よって、変換するための関数の要件は以下になります。

  • 線形結合が増加するにつれて、変換後の値も増加する（単調増加）
  • [-\infty, \infty]を[0, 1]にする

これには確率分布の累積分布関数（分布関数）が使えます。累積分布関数は単調増加関数であり、

\displaystyle{ \lim_{x \to -\infty}F(x) = 0, \lim_{x \to +\infty}F(x) = 1 }

という性質を持っています。この数式の意味するところは、xを-\inftyに持っていくと結果が0に近づき、xを\inftyに持っていくと結果が0に近づくということです。この性質は要件にピッタリです。この性質は要件にピッタリです。ただし、どんな確率分布の分布関数でも良いわけではなく、たとえば正規分布のように確率変数の範囲が[-\infty, \infty]である必要はあります。下図は標準正規分布の確率密度関数（左）と分布関数（右）です。

f:id:matakanobu:20200703172330p:plain
標準正規分布の確率密度関数と分布関数

分布関数をFとして表現すると、

\displaystyle{ p_i = F(\boldsymbol{x}^\mathrm{T}_i\boldsymbol{\beta}) \tag{1} }

となります。線形結合を分布関数を用いて確率p_iに変換しているという意味です。この分布関数に、標準正規分布の分布関数である

\displaystyle{ F(x) = \int^x_{-\infty}\frac{1}{\sqrt{2 \pi}} \exp{\left(-\frac{x^2}{2}\right)} \mathrm{d}x }

を用いたモデルをプロビットモデルと言います。しかし、このモデルは分布関数に積分があるのでちょっと抵抗があります。分布関数に積分が無いものが扱いやすいですが、そのような確率分布はあるのでしょうか？それがロジスティック分布です。（標準）ロジスティック分布の分布関数はこちらです。

\displaystyle{ F(x) = \frac{1}{1 + e^{-x}} }

f:id:matakanobu:20200703174248p:plain
ロジスティック分布の確率密度関数と分布関数

こちらは標準正規分布と違って積分が無いので扱いやすいです。分母と分子にe^xを掛けると、

\displaystyle{ F(x) = \frac{e^{x}}{1 + e^{x}} \tag{2} }

と変形ができ、(1)式と統合すると、

\displaystyle{ p_i = \frac{\exp(\boldsymbol{x}^\mathrm{T}_i\boldsymbol{\beta})}{1 + \exp (\boldsymbol{x}^\mathrm{T}_i\boldsymbol{\beta})} \tag{3} }

となります。これがロジスティック回帰モデル（ロジットモデル）です。

ロジスティック回帰モデルの数式を考察する

次にロジスティック回帰モデルの(3)式について、線形結合\boldsymbol{x}^\mathrm{T}_i\boldsymbol{\beta}=の形にしてみます。これはロジスティック分布の分布関数の逆関数で下式になります。

\displaystyle{ \log \frac{p_i}{1-p_i} = \boldsymbol{x}^\mathrm{T}_i\boldsymbol{\beta} \tag{4} }

左辺を対数オッズ（またはロジット）と言います。オッズは

\displaystyle{ オッズ = \frac{p}{1-p} }

で表され、ある事象が発生する確率を発生しない確率で割った値です。なので、

  • odds > 1 : 発生する確率のほうが大きい
  • odds < 1 : 発生しない確率のほうが大きい
  • odds = 1 : 発生する確率としない確率が等しい

です。対数を取る理由は発生する確率が大きくなる（p_iが大きくなる）と非常に大きな値になるから、そのスケールを抑えるためです。対数オッズを計算することはロジット変換と呼ぶこともあります。

\displaystyle{ 対数オッズ：\mathrm{logit}(p) = \log \frac{p}{1-p} \tag{5} }

これは確率（[0, 1]）を[-\infty, \infty]の範囲のデータに変換することを意味します。ロジット変換は確率や比率に対する特徴量エンジニアリングとしても使うことができ、これによって線形回帰モデルを推定することも可能になります。

話をロジスティック回帰モデルに戻すと、(4)式からロジスティック回帰モデルは対数オッズの形で表されるモデルであるということがわかります。

ロジスティック回帰モデルの解釈

たとえば、勉強時間から合格する確率を説明（予測）するとして、ロジスティック回帰モデルを用いることにします。得られたデータから、勉強時間の回帰係数は以下のように推定されたものとします。

  • 回帰係数:0.05
  • 標準誤差:0.01

勉強時間の単位は1時間であるとします。このとき、勉強時間が増えた場合に合格確率はどう変化するかを考えてみます。線形回帰モデルのときはダイレクトに回帰係数の値を解釈すればよかったですが、ロジスティック回帰モデルでは変換が行われているため、そう単純でもなさそうです。そのために、まず(4)式を変形しますが、左辺の\logを外すため両辺に指数を取ります。

\displaystyle{ \frac{p_i}{1-p_i} = \exp(\boldsymbol{x}^\mathrm{T}_i\boldsymbol{\beta}) \tag{6} }

これは上の話で出てきたオッズの形です。ベクトル\boldsymbol{x}_iの1要因（変数）に着目するので、

\displaystyle{ \frac{p_i}{1-p_i} = \exp(x_{i1}\beta_1 + \beta_0) }

ここで\beta_0は切片項です。x_{i1}が1単位増えた場合とそうでない場合を比較したいので、それらの比を取ります。

\displaystyle{ \frac{\exp\{(x_{i1}+1)\beta_1 + \beta_0)\}}{\exp(x_{i1}\beta_1 + \beta_0)} = \exp(\beta_1) \tag{7} }

これはオッズの比を取っているので、オッズ比と言います。さて、この式が意味するところは、x_{i1}が1単位増えたときのオッズは\exp(\beta_1)倍になるということです。仮に1単位増えたときに変化がないとした場合（これは回帰係数が0の場合）は分子と分母が同じになり、1になります。また、回帰係数が負の場合はオッズ比は1より小さくなります。これらの関係をまとめると以下のようになります。

  • オッズ比 > 1 : 可能性の増加（回帰係数が正）
  • オッズ比 < 1 : 可能性の減少（回帰係数が負）
  • オッズ比 = 1 : 変化なし（回帰係数が0）

実際に値を入れてみると理解できると思いますので、先の勉強時間と合格の例で計算してみます。

\displaystyle{ \exp(0.05) - 1 = 0.051(5.1\%) }

よって、1時間勉強時間を増やすと合格する確率は5.1%上昇するという解釈になります。標準誤差を用いれば95%信頼区間も構成できて、

\displaystyle{ \exp(0.05 \pm 2 \times 0.01) = (1.030, 1.073) }

となるので、およそ3〜7%の上昇であると見積もれます。

回帰係数の推定

ロジスティック回帰モデルの回帰係数は最尤推定を用いて推定するのが一般的です。y_iをラベルとして、y_i=1のときの確率をp_iとすればy_i=0のときの確率は1-p_iと書けます。尤度関数は確率の積になりますので、

\displaystyle{ L(\boldsymbol{\beta}) = \prod^n_{i=1} p^{y_i}_i(1-p_i)^{1-y_i} }

これはベルヌーイ分布の尤度関数です。なお、p_iはロジスティック回帰モデル(3)式であることに注意してください。これに対数を取った対数尤度関数は

\displaystyle{ l(\boldsymbol{\beta}) = \sum^n_{i=1} \left\{ y_i \log p_i + (1 - y_i) \log (1-p_i) \right\} \tag{8} }

となります。この対数尤度が最大となる\boldsymbol{\beta}が解になりますので、\boldsymbol{\beta}で偏微分して0として解きます。しかしながら、これは解析的に解くことはできず数値的に求める必要があります。その手法としてニュートン・ラフソンアルゴリズムがありますが、今回の説明はここまでといたします。

参考

  • 現代数理統計学の基礎
  • 統計的学習の基礎