.. _sec_linear_regression:
線形回帰
========
*回帰* 問題は、数値を予測したいときに現れる。
代表的な例としては、価格(住宅、株式など)の予測、
入院患者の滞在期間の予測、 需要予測(小売売上)など、数多くある。
すべての予測問題が古典的な回帰であるわけではない。
後ほど、カテゴリ集合の中での所属を予測することを目的とする
分類問題を紹介する。
例として、家の面積(平方フィート)と築年数(年)に基づいて、
住宅価格(ドル)を推定したいとする。
住宅価格を予測するモデルを構築するには、
各住宅について売却価格、面積、築年数を含むデータを
手に入れる必要がある。 機械学習の用語では、このデータセットは
*訓練データセット* または *訓練集合* と呼ばれ、
各行(1件の売買に対応するデータを含む)は *例*\ (または
*データ例*\ 、\ *インスタンス*\ 、\ *サンプル*\ )と呼ばれる。
私たちが予測しようとしているもの(価格)は *ラベル*\ (または
*ターゲット*\ )と呼ばれる。 予測の基礎となる変数(築年数と面積)は
*特徴量*\ (または *共変量*\ )と呼ばれる。
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
%matplotlib inline
from d2l import torch as d2l
import math
import torch
import numpy as np
import time
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
%matplotlib inline
from d2l import mxnet as d2l
import math
from mxnet import np
import time
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
%matplotlib inline
from d2l import jax as d2l
from jax import numpy as jnp
import math
import time
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
No GPU/TPU found, falling back to CPU. (Set TF_CPP_MIN_LOG_LEVEL=0 and rerun for more info.)
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
%matplotlib inline
from d2l import tensorflow as d2l
import math
import tensorflow as tf
import numpy as np
import time
.. raw:: html
.. raw:: html
基礎
----
*線形回帰* は、回帰問題に取り組むための標準的な手法の中で
最も単純であり、かつ最も広く使われている。
19世紀初頭にさかのぼる線形回帰 :cite:`Legendre.1805,Gauss.1809` は、
いくつかの単純な仮定から導かれる。 まず、特徴量 :math:`\mathbf{x}`
とターゲット :math:`y` の関係が おおむね線形である、
すなわち条件付き平均 :math:`E[Y \mid X=\mathbf{x}]` が 特徴量
:math:`\mathbf{x}` の重み付き和として表せると仮定する。
この設定では、観測ノイズのために
ターゲット値が期待値からずれることは許される。 次に、そのようなノイズが
ガウス分布に従う、よく振る舞うものだと仮定できる。 通常、\ :math:`n`
をデータセット中の例の数を表すのに用いる。
サンプルとターゲットの列挙には上付き添字を、
座標の添字には下付き添字を用いる。 より具体的には、
:math:`\mathbf{x}^{(i)}` は :math:`i` 番目のサンプルを表し、
:math:`x_j^{(i)}` はその :math:`j` 番目の座標を表す。
.. _subsec_linear_model:
モデル
~~~~~~
あらゆる解法の中心には、
特徴量をターゲットの推定値へと変換する方法を記述するモデルがある。
線形性の仮定とは、ターゲット(価格)の期待値が
特徴量(面積と築年数)の重み付き和として表せることを意味する。
.. math:: \textrm{price} = w_{\textrm{area}} \cdot \textrm{area} + w_{\textrm{age}} \cdot \textrm{age} + b.
:label: eq_price-area
ここで :math:`w_{\textrm{area}}` と :math:`w_{\textrm{age}}` は *重み*
と呼ばれ、\ :math:`b` は *バイアス* (または
*オフセット*\ 、\ *切片*\ )と呼ばれる。
重みは、各特徴量が予測に与える影響を決定する。
バイアスは、すべての特徴量が 0 のときの推定値を決定する。
新築住宅で面積が厳密に 0 のものは決して見られないが、
それでもバイアスが必要なのは、
特徴量のすべての線形関数を表現できるようにするためです
(原点を通る直線に限定されないようにするためです)。
厳密には、:eq:`eq_price-area` は入力特徴量の *アフィン変換*
であり、 特徴量の *線形変換*\ (重み付き和)と、
加えられたバイアスによる *平行移動* から成る。
データセットが与えられたとき、私たちの目標は 重み :math:`\mathbf{w}`
とバイアス :math:`b` を選び、
平均的に見て、モデルの予測がデータ中で観測された真の価格に
できるだけよく一致するようにすることである。
特徴量が少数のデータセットに注目することが一般的な分野では、
:eq:`eq_price-area`
のようにモデルを明示的に長い形で表すことがよくある。
機械学習では、通常は高次元データセットを扱うため、
コンパクトな線形代数の記法を用いる方が便利である。 入力が :math:`d`
個の特徴量からなるとき、 それぞれに(1 から :math:`d`
までの)添字を割り当て、 予測 :math:`\hat{y}` を
(一般に「ハット」記号は推定値を表す)次のように表せる。
.. math:: \hat{y} = w_1 x_1 + \cdots + w_d x_d + b.
すべての特徴量をベクトル :math:`\mathbf{x} \in \mathbb{R}^d` に、
すべての重みをベクトル :math:`\mathbf{w} \in \mathbb{R}^d`
にまとめると、 モデルは :math:`\mathbf{w}` と :math:`\mathbf{x}`
の内積を用いて簡潔に表せる。
.. math:: \hat{y} = \mathbf{w}^\top \mathbf{x} + b.
:label: eq_linreg-y
:eq:`eq_linreg-y` では、ベクトル :math:`\mathbf{x}` は
1つの例の特徴量に対応する。 データセット全体が :math:`n`
個の例からなるときには、 *設計行列*
:math:`\mathbf{X} \in \mathbb{R}^{n \times d}` を用いて
特徴量を表すと便利なことがよくある。 ここで、\ :math:`\mathbf{X}`
は各例に対して1行、 各特徴量に対して1列を持つ。 特徴量の集合
:math:`\mathbf{X}` に対して、 予測
:math:`\hat{\mathbf{y}} \in \mathbb{R}^n` は 行列–ベクトル積で表せる。
.. math:: {\hat{\mathbf{y}}} = \mathbf{X} \mathbf{w} + b,
:label: eq_linreg-y-vec
ここでは、和を取る際にブロードキャスト(:numref:`subsec_broadcasting`\ )が適用される。
訓練データセットの特徴量 :math:`\mathbf{X}` と 対応する(既知の)ラベル
:math:`\mathbf{y}` が与えられたとき、 線形回帰の目標は、
:math:`\mathbf{X}`
と同じ分布からサンプリングされた新しいデータ例の特徴量が与えられたときに、
その新しい例のラベルを(期待値として)最小誤差で予測できるような
重みベクトル :math:`\mathbf{w}` とバイアス項 :math:`b`
を見つけることである。
たとえ、\ :math:`\mathbf{x}` が与えられたときの :math:`y`
を予測する最良のモデルが 線形であると信じていたとしても、 現実世界の
:math:`n` 個の例からなるデータセットで すべての :math:`1 \leq i \leq n`
に対して :math:`y^{(i)} = \mathbf{w}^\top \mathbf{x}^{(i)}+b`
が厳密に成り立つとは期待しない。 たとえば、\ :math:`\mathbf{X}` と
:math:`\mathbf{y}` を観測するためにどのような機器を使っても、
わずかな測定誤差があるかもしれない。
したがって、基礎となる関係が線形であると確信している場合でも、
そのような誤差を考慮するためにノイズ項を組み込みる。
最良の *パラメータ* (または *モデルパラメータ*\ )\ :math:`\mathbf{w}`
と :math:`b` を探し始める前に、 さらに2つ必要である。 (i)
与えられたモデルの品質を測る尺度、 (ii)
モデルの品質を改善するためにモデルを更新する手順である。
.. _subsec_linear-regression-loss-function:
損失関数
~~~~~~~~
当然ながら、モデルをデータに適合させるには、 *適合度*\ (あるいは同値に
*不適合度*\ )の尺度について合意する必要がある。 *損失関数*
は、ターゲットの *真の* 値と *予測* 値の間の距離を定量化する。
損失は通常、非負の数であり、 小さいほど良く、完全な予測では損失 0
になる。 回帰問題で最も一般的な損失関数は二乗誤差である。 例 :math:`i`
に対する予測が :math:`\hat{y}^{(i)}`\ 、 対応する真のラベルが
:math:`y^{(i)}` のとき、 *二乗誤差* は次のように与えられる。
.. math:: l^{(i)}(\mathbf{w}, b) = \frac{1}{2} \left(\hat{y}^{(i)} - y^{(i)}\right)^2.
:label: eq_mse
定数 :math:`\frac{1}{2}` は本質的には違いを生まないが、
損失の微分を取るときに打ち消し合うので、記法上便利である。
訓練データセットは与えられており、したがって私たちの制御外にあるため、
経験誤差はモデルパラメータの関数にすぎない。 :numref:`fig_fit_linreg`
では、1次元入力の問題における 線形回帰モデルの適合を可視化している。
.. _fig_fit_linreg:
.. figure:: ../img/fit-linreg.svg
1次元データに線形回帰モデルを適合させる。
推定値 :math:`\hat{y}^{(i)}` とターゲット :math:`y^{(i)}`
の差が大きいと、 その二次形式のために損失への寄与はさらに大きくなる
(この二次性は諸刃の剣である。大きな誤差を避けるようモデルを促す一方で、
異常データに対して過度に敏感になることもある)。 :math:`n`
個の例からなるデータセット全体に対するモデルの品質を測るには、
訓練集合上の損失を単純に平均(あるいは同値に総和)すればよいである。
.. math:: L(\mathbf{w}, b) =\frac{1}{n}\sum_{i=1}^n l^{(i)}(\mathbf{w}, b) =\frac{1}{n} \sum_{i=1}^n \frac{1}{2}\left(\mathbf{w}^\top \mathbf{x}^{(i)} + b - y^{(i)}\right)^2.
モデルを訓練するとき、私たちは すべての訓練例にわたる総損失を最小にする
パラメータ (:math:`\mathbf{w}^*, b^*`) を求める。
.. math:: \mathbf{w}^*, b^* = \operatorname*{argmin}_{\mathbf{w}, b}\ L(\mathbf{w}, b).
解析解
~~~~~~
私たちが扱うほとんどのモデルとは異なり、
線形回帰は驚くほど簡単な最適化問題を与えてくれる。
特に、単純な公式を適用することで、
(訓練データ上で評価した)最適パラメータを解析的に求めることができる。
まず、設計行列にすべて 1 からなる列を追加することで、 バイアス :math:`b`
をパラメータ :math:`\mathbf{w}` に組み込める。 すると、予測問題は
:math:`\|\mathbf{y} - \mathbf{X}\mathbf{w}\|^2` を最小化することになる。
設計行列 :math:`\mathbf{X}` がフルランク
(どの特徴量も他の特徴量の線形従属でない)である限り、
損失面にはただ1つの臨界点しかなく、
それは定義域全体にわたる損失の最小値に対応する。 損失を
:math:`\mathbf{w}` で微分して 0 に等置すると、次が得られる。
.. math::
\begin{aligned}
\partial_{\mathbf{w}} \|\mathbf{y} - \mathbf{X}\mathbf{w}\|^2 =
2 \mathbf{X}^\top (\mathbf{X} \mathbf{w} - \mathbf{y}) = 0
\textrm{ and hence }
\mathbf{X}^\top \mathbf{y} = \mathbf{X}^\top \mathbf{X} \mathbf{w}.
\end{aligned}
:math:`\mathbf{w}` について解けば、最適化問題の最適解が得られる。
この解は
.. math:: \mathbf{w}^* = (\mathbf X^\top \mathbf X)^{-1}\mathbf X^\top \mathbf{y}
が一意であるのは、 行列 :math:`\mathbf X^\top \mathbf X` が可逆な場合、
すなわち設計行列の列が線形独立である場合に限られます
:cite:`Golub.Van-Loan.1996`\ 。
線形回帰のような単純な問題では 解析解が得られることがあるが、
そのような幸運に慣れてはいけない。
解析解は美しい数学的解析を可能にするが、
解析解を要求することは非常に制約が強く、
深層学習のほとんどすべての魅力的な側面を排除してしまいる。
ミニバッチ確率的勾配降下法
~~~~~~~~~~~~~~~~~~~~~~~~~~
幸いなことに、モデルを解析的に解けない場合でも、
実際にはしばしば効果的にモデルを訓練できる。
さらに、多くのタスクでは、そのような最適化の難しいモデルの方が
はるかに優れているため、訓練方法を見つける手間は
十分に見合うものになる。
ほぼすべての深層学習モデルを最適化するための鍵となる技法であり、
本書全体を通して用いることになるのが、
パラメータを少しずつ損失関数を下げる方向に更新することで
誤差を反復的に減らす方法である。 このアルゴリズムは *勾配降下法*
と呼ばれる。
勾配降下法の最も素朴な適用は、
データセット中のすべての例で計算された損失の平均である
損失関数の微分を取ることである。
実際には、これは非常に遅くなることがある。 更新が非常に強力であっても、
1回更新する前にデータセット全体を1周しなければならないからです
:cite:`Liu.Nocedal.1989`\ 。
さらに悪いことに、訓練データに冗長性が多い場合、
フル更新の利点は限られる。
もう一方の極端は、一度に1つの例だけを考え、
1つの観測に基づいて更新することである。
その結果得られるアルゴリズムである *確率的勾配降下法*\ (SGD)は、
大規模データセットに対しても有効な戦略になりえます
:cite:`Bottou.2010`\ 。 残念ながら、SGD
には計算上および統計上の欠点がある。 1つの問題は、プロセッサが
主記憶からプロセッサキャッシュへデータを移動するよりも、
数値の乗算や加算をはるかに高速に行えるという事実から生じる。
行列–ベクトル積を実行する方が、
対応する数のベクトル–ベクトル演算を行うよりも、
最大で1桁程度効率的である。 これは、1サンプルずつ処理するのに
フルバッチよりもはるかに時間がかかりうることを意味する。
2つ目の問題は、バッチ正規化(:numref:`sec_batch_norm`
で説明する)のような一部の層は、
一度に複数の観測にアクセスできる場合にのみ うまく機能することである。
これら2つの問題に対する解決策は、中間的な戦略を選ぶことである。
フルバッチでも1サンプルずつでもなく、 観測の *ミニバッチ* を取ります
:cite:`Li.Zhang.Chen.ea.2014`\ 。
このミニバッチのサイズの具体的な選択は、
メモリ量、アクセラレータの数、層の選択、データセット全体のサイズなど、
多くの要因に依存する。 それでも、32 から 256 の間の数、 できれば 2
の大きな冪の倍数がよい出発点である。 これが *ミニバッチ確率的勾配降下法*
につながりる。
最も基本的な形では、各反復 :math:`t` において、 まず固定個数
:math:`|\mathcal{B}|` の訓練例からなるミニバッチ :math:`\mathcal{B}_t`
をランダムにサンプリングする。
次に、ミニバッチ上の平均損失をモデルパラメータに関して微分(勾配計算)する。
最後に、その勾配に *学習率* と呼ばれるあらかじめ定めた小さな正の値
:math:`\eta` を掛け、 得られた項を現在のパラメータ値から引きる。
更新は次のように表せる。
.. math:: (\mathbf{w},b) \leftarrow (\mathbf{w},b) - \frac{\eta}{|\mathcal{B}|} \sum_{i \in \mathcal{B}_t} \partial_{(\mathbf{w},b)} l^{(i)}(\mathbf{w},b).
要するに、ミニバッチ SGD は次のように進みる。 (i)
通常はランダムに、モデルパラメータの値を初期化する。 (ii)
データからランダムなミニバッチを反復的にサンプリングし、
負の勾配の方向にパラメータを更新する。
二次損失とアフィン変換に対しては、これは閉形式で展開できる。
.. math:: \begin{aligned} \mathbf{w} & \leftarrow \mathbf{w} - \frac{\eta}{|\mathcal{B}|} \sum_{i \in \mathcal{B}_t} \partial_{\mathbf{w}} l^{(i)}(\mathbf{w}, b) && = \mathbf{w} - \frac{\eta}{|\mathcal{B}|} \sum_{i \in \mathcal{B}_t} \mathbf{x}^{(i)} \left(\mathbf{w}^\top \mathbf{x}^{(i)} + b - y^{(i)}\right)\\ b &\leftarrow b - \frac{\eta}{|\mathcal{B}|} \sum_{i \in \mathcal{B}_t} \partial_b l^{(i)}(\mathbf{w}, b) && = b - \frac{\eta}{|\mathcal{B}|} \sum_{i \in \mathcal{B}_t} \left(\mathbf{w}^\top \mathbf{x}^{(i)} + b - y^{(i)}\right). \end{aligned}
:label: eq_linreg_batch_update
ミニバッチ :math:`\mathcal{B}` を選ぶので、 そのサイズ
:math:`|\mathcal{B}|` で正規化する必要がある。
しばしば、ミニバッチサイズと学習率はユーザが定義する。
訓練ループの中で更新されないこのような調整可能なパラメータは
*ハイパーパラメータ* と呼ばれる。
それらはベイズ最適化などのいくつかの手法で自動的に調整できる
:cite:`Frazier.2018`\ 。 最終的に、解の品質は通常、別の
*検証データセット*\ (または *検証集合*\ )で評価される。
あらかじめ定めた回数だけ反復して訓練した後
(あるいは別の停止条件が満たされた後)、 推定されたモデルパラメータ
:math:`\hat{\mathbf{w}}, \hat{b}` を記録する。
たとえ関数が本当に線形でノイズがないとしても、
これらのパラメータは損失の厳密な最小化解にはならず、
ましてや決定論的でもない。
アルゴリズムは最小化解に向かってゆっくり収束するが、
有限回のステップでそれらを厳密に見つけることは通常ない。
さらに、パラメータ更新に用いるミニバッチ :math:`\mathcal{B}` は
ランダムに選ばれる。 これにより決定論性は失われる。
線形回帰は、たまたま大域的最小値を持つ学習問題です
(\ :math:`\mathbf{X}` がフルランクであるとき、あるいは同値に
:math:`\mathbf{X}^\top \mathbf{X}` が可逆であるとき)。
しかし、深層ネットワークの損失面には多くの鞍点と最小値がある。
幸いなことに、私たちは通常、
厳密なパラメータ集合を見つけることには関心がなく、
正確な予測(したがって低い損失)を与える
何らかのパラメータ集合が得られれば十分である。
実際には、深層学習の実践者は
訓練集合上で損失を最小化するパラメータを見つけるのに
苦労することはめったにない
:cite:`Izmailov.Podoprikhin.Garipov.ea.2018,Frankle.Carbin.2018`\ 。
より困難なのは、これまで見たことのないデータに対して
正確な予測を与えるパラメータを見つけることである。 この課題は *汎化*
と呼ばれる。 これらの話題には本書を通して何度も立ち返りる。
予測
~~~~
モデル :math:`\hat{\mathbf{w}}^\top \mathbf{x} + \hat{b}`
が与えられれば、 新しい例に対して *予測* を行える。 たとえば、面積
:math:`x_1` と築年数 :math:`x_2` が与えられたときに、
これまで見たことのない住宅の売却価格を予測できる。
深層学習の実践者は予測段階を *推論* と呼ぶようになっているが、
これは少し不正確である。\ *推論* は一般に、
証拠に基づいて到達したあらゆる結論を指し、
パラメータの値と、見えないインスタンスのラベルの可能性の両方を含む。
むしろ統計学の文献では、 *推論* はパラメータ推論を指すことの方が多く、
この用語の多義性は、深層学習の実践者が統計学者と話すときに
不要な混乱を生みる。 以下では、可能な限り *予測* を用いる。
高速化のためのベクトル化
------------------------
モデルを訓練するとき、通常は
例のミニバッチ全体を同時に処理したいと考える。
これを効率よく行うには、計算をベクトル化し、 Python で高価な for
ループを書くのではなく、 高速な線形代数ライブラリを活用する必要がある。
これがなぜそれほど重要なのかを見るために、
ベクトルを加算する2つの方法を考えてみよう。 まず、すべて 1 からなる
10,000 次元ベクトルを2つ用意する。 1つ目の方法では、Python の for
ループでベクトルを走査する。 2つ目では、\ ``+``
を1回呼び出すだけに頼りる。
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
n = 10000
a = d2l.ones(n)
b = d2l.ones(n)
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
[06:59:48] ../src/storage/storage.cc:196: Using Pooled (Naive) StorageManager for CPU
これで処理時間をベンチマークできる。 まず、for
ループを使って、1座標ずつ加算する。
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
c = d2l.zeros(n)
t = time.time()
for i in range(n):
c[i] = a[i] + b[i]
f'{time.time() - t:.5f} sec'
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
'0.12482 sec'
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
c = d2l.zeros(n)
t = time.time()
for i in range(n):
c[i] = a[i] + b[i]
f'{time.time() - t:.5f} sec'
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
'4.90001 sec'
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# JAX arrays are immutable, meaning that once created their contents
# cannot be changed. For updating individual elements, JAX provides
# an indexed update syntax that returns an updated copy
c = d2l.zeros(n)
t = time.time()
for i in range(n):
c = c.at[i].set(a[i] + b[i])
f'{time.time() - t:.5f} sec'
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
'24.00628 sec'
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
c = tf.Variable(d2l.zeros(n))
t = time.time()
for i in range(n):
c[i].assign(a[i] + b[i])
f'{time.time() - t:.5f} sec'
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
'12.79235 sec'
.. raw:: html
.. raw:: html
あるいは、再ロードされた ``+`` 演算子を使って要素ごとの和を計算する。
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
t = time.time()
d = a + b
f'{time.time() - t:.5f} sec'
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
'0.00029 sec'
2つ目の方法は1つ目よりも劇的に高速である。
コードをベクトル化すると、しばしば桁違いの高速化が得られる。
さらに、数学の多くをライブラリ側に任せられるので、
自分で計算を書く量が減り、
誤りの可能性が下がり、コードの移植性も高まりる。
.. _subsec_normal_distribution_and_squared_loss:
正規分布と二乗損失
------------------
ここまで、二乗損失目的関数についてかなり実用的な動機づけを与えてきた。
すなわち、基礎となるパターンが本当に線形であるなら、
最適パラメータは条件付き期待値 :math:`E[Y\mid X]` を返し、
またこの損失は外れ値に大きな罰則を与える。
さらに、ノイズの分布について確率的な仮定を置くことで、
二乗損失目的関数に対してより形式的な動機づけも与えられる。
線形回帰は19世紀の変わり目に発明された。
ガウスとルジャンドルのどちらが最初にこのアイデアを思いついたかは長く議論されてきたが、
正規分布(\ *ガウス分布* とも呼ばれます)を発見したのもガウスでした。
正規分布と二乗損失を伴う線形回帰は、
単に同じ祖先を持つという以上の深い関係を共有していることがわかる。
まず、平均 :math:`\mu`\ 、分散 :math:`\sigma^2`\ (標準偏差
:math:`\sigma`\ )の正規分布は次で与えられる。
.. math:: p(x) = \frac{1}{\sqrt{2 \pi \sigma^2}} \exp\left(-\frac{1}{2 \sigma^2} (x - \mu)^2\right).
以下では 正規分布を計算する関数を定義する。
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
def normal(x, mu, sigma):
p = 1 / math.sqrt(2 * math.pi * sigma**2)
if tab.selected('jax'):
return p * jnp.exp(-0.5 * (x - mu)**2 / sigma**2)
if tab.selected('pytorch', 'mxnet', 'tensorflow'):
return p * np.exp(-0.5 * (x - mu)**2 / sigma**2)
これで 正規分布を可視化 できる。
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
if tab.selected('jax'):
# Use JAX NumPy for visualization
x = jnp.arange(-7, 7, 0.01)
if tab.selected('pytorch', 'mxnet', 'tensorflow'):
# Use NumPy again for visualization
x = np.arange(-7, 7, 0.01)
# Mean and standard deviation pairs
params = [(0, 1), (0, 2), (3, 1)]
d2l.plot(x, [normal(x, mu, sigma) for mu, sigma in params], xlabel='x',
ylabel='p(x)', figsize=(4.5, 2.5),
legend=[f'mean {mu}, std {sigma}' for mu, sigma in params])
.. figure:: output_linear-regression_433d43_39_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Use NumPy again for visualization
x = np.arange(-7, 7, 0.01)
# Mean and standard deviation pairs
params = [(0, 1), (0, 2), (3, 1)]
d2l.plot(x.asnumpy(), [normal(x, mu, sigma).asnumpy() for mu, sigma in params], xlabel='x',
ylabel='p(x)', figsize=(4.5, 2.5),
legend=[f'mean {mu}, std {sigma}' for mu, sigma in params])
.. figure:: output_linear-regression_433d43_42_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
if tab.selected('jax'):
# Use JAX NumPy for visualization
x = jnp.arange(-7, 7, 0.01)
if tab.selected('pytorch', 'mxnet', 'tensorflow'):
# Use NumPy again for visualization
x = np.arange(-7, 7, 0.01)
# Mean and standard deviation pairs
params = [(0, 1), (0, 2), (3, 1)]
d2l.plot(x, [normal(x, mu, sigma) for mu, sigma in params], xlabel='x',
ylabel='p(x)', figsize=(4.5, 2.5),
legend=[f'mean {mu}, std {sigma}' for mu, sigma in params])
.. figure:: output_linear-regression_433d43_45_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
if tab.selected('jax'):
# Use JAX NumPy for visualization
x = jnp.arange(-7, 7, 0.01)
if tab.selected('pytorch', 'mxnet', 'tensorflow'):
# Use NumPy again for visualization
x = np.arange(-7, 7, 0.01)
# Mean and standard deviation pairs
params = [(0, 1), (0, 2), (3, 1)]
d2l.plot(x, [normal(x, mu, sigma) for mu, sigma in params], xlabel='x',
ylabel='p(x)', figsize=(4.5, 2.5),
legend=[f'mean {mu}, std {sigma}' for mu, sigma in params])
.. figure:: output_linear-regression_433d43_48_0.svg
.. raw:: html
.. raw:: html
平均を変えると :math:`x` 軸に沿った平行移動に対応し、
分散を大きくすると分布は広がり、 ピークは低くなる。
二乗損失を伴う線形回帰を動機づける1つの方法は、
観測がノイズ付き測定から生じると仮定し、 そのノイズ :math:`\epsilon`
が正規分布 :math:`\mathcal{N}(0, \sigma^2)` に従うとすることである。
.. math:: y = \mathbf{w}^\top \mathbf{x} + b + \epsilon \textrm{ where } \epsilon \sim \mathcal{N}(0, \sigma^2).
したがって、与えられた :math:`\mathbf{x}` に対して特定の :math:`y`
が観測される *尤度* は次のように書ける。
.. math:: P(y \mid \mathbf{x}) = \frac{1}{\sqrt{2 \pi \sigma^2}} \exp\left(-\frac{1}{2 \sigma^2} (y - \mathbf{w}^\top \mathbf{x} - b)^2\right).
このように、尤度は因数分解される。 *最尤原理* によれば、パラメータ
:math:`\mathbf{w}` と :math:`b` の最良の値は、 データセット全体の *尤度*
を最大化するものである。
.. math:: P(\mathbf y \mid \mathbf X) = \prod_{i=1}^{n} p(y^{(i)} \mid \mathbf{x}^{(i)}).
この等式は、すべての組 :math:`(\mathbf{x}^{(i)}, y^{(i)})` が
互いに独立に抽出されたために成り立ちる。
最尤原理に従って選ばれる推定量は *最尤推定量* と呼ばれる。
多くの指数関数の積を最大化するのは難しそうに見えるが、
目的関数を変えずに、尤度の対数を最大化することで 大幅に簡単にできる。
歴史的な理由から、最適化は最大化よりも
最小化として表現されることが多いである。 そこで、何も変えずに、
*負の対数尤度* を *最小化* できる。これは次のように表せる。
.. math:: -\log P(\mathbf y \mid \mathbf X) = \sum_{i=1}^n \frac{1}{2} \log(2 \pi \sigma^2) + \frac{1}{2 \sigma^2} \left(y^{(i)} - \mathbf{w}^\top \mathbf{x}^{(i)} - b\right)^2.
:math:`\sigma` が固定であると仮定すれば、 第1項は :math:`\mathbf{w}` や
:math:`b` に依存しないので無視できる。
第2項は、前に導入した二乗誤差損失と同じであり、 乗法定数
:math:`\frac{1}{\sigma^2}` が付いているだけである。 幸いなことに、解は
:math:`\\sigma` にも依存しない。
したがって、平均二乗誤差を最小化することは、
加法的ガウスノイズを仮定した線形モデルの最尤推定と等価である。
ニューラルネットワークとしての線形回帰
--------------------------------------
線形モデルは、本書でこれから紹介する
複雑なネットワークを表現するには十分に豊かではないが、
(人工)ニューラルネットワークは、
すべての特徴量が入力ニューロンで表され、
それらがすべて出力に直接接続されたネットワークとして
線形モデルを包含できるほど豊かである。
:numref:`fig_single_neuron` は
線形回帰をニューラルネットワークとして描いている。
この図は、各入力が出力にどのように接続されているかといった
接続パターンを強調しているが、
重みやバイアスが取る具体的な値は示していない。
.. _fig_single_neuron:
.. figure:: ../img/singleneuron.svg
線形回帰は1層ニューラルネットワークである。
入力は :math:`x_1, \ldots, x_d` である。 :math:`d` を入力層における
*入力数* または *特徴次元* と呼ぶ。 ネットワークの出力は :math:`o_1`
である。 私たちは単一の数値を予測したいだけなので、
出力ニューロンは1つしかない。 入力値はすべて *与えられた*
ものであることに注意されたい。 *計算される* ニューロンは1つだけである。
要するに、線形回帰は 1層の全結合ニューラルネットワークとみなせる。
後の章では、はるかに多くの層を持つネットワークに出会いる。
生物学
~~~~~~
線形回帰は計算神経科学よりも前に登場したため、
線形回帰をニューラルネットワークの観点で説明するのは
時代錯誤に見えるかもしれない。
それでも、サイバネティクス研究者で神経生理学者の
ウォーレン・マカロックとウォルター・ピッツが
人工ニューロンのモデルを開発し始めたとき、 それは自然な出発点でした。
:numref:`fig_Neuron` に示す生物学的ニューロンの
漫画的な図を考えてみよう。 そこには、\ *樹状突起*\ (入力端子)、
*細胞核*\ (CPU)、\ *軸索*\ (出力線)、 および
*軸索終末*\ (出力端子)があり、 *シナプス*
を介して他のニューロンとの接続を可能にしている。
.. _fig_Neuron:
.. figure:: ../img/neuron.svg
実際のニューロン(出典: 米国国立がん研究所 Surveillance, Epidemiology
and End Results (SEER) Program の “Anatomy and Physiology”)。
他のニューロン(あるいは環境センサー)から到着する情報 :math:`x_i`
は樹状突起で受け取られる。 特に、その情報は *シナプス重み* :math:`w_i`
によって重み付けされ、 積 :math:`x_i w_i`
による活性化や抑制など、入力の効果を決定する。
複数のソースから到着した重み付き入力は 細胞核で重み付き和
:math:`y = \sum_i x_i w_i + b` として集約され、 場合によっては関数
:math:`\sigma(y)` による非線形な後処理を受ける。
この情報はその後、軸索を通って軸索終末へ送られ、 そこで目的地
(たとえば筋肉のようなアクチュエータ)に到達するか、
あるいは樹状突起を介して別のニューロンに入力される。
もちろん、正しい接続性と学習アルゴリズムを備えていれば、
このような多数の単位を組み合わせることで、
単一のニューロンだけでは表現できない
はるかに興味深く複雑な振る舞いを生み出せるという 高レベルの考え方は、
実際の生物学的神経系の研究から生まれたものである。
同時に、今日の深層学習研究の大部分は より広い源泉から着想を得ている。
:cite:t:`Russell.Norvig.2016` を引く。 彼らは、飛行機が鳥に
*着想を得た* ものであっても、
鳥類学がここ数世紀の航空工学革新の主たる原動力ではなかったと指摘した。
同様に、今日の深層学習における着想は、
数学、言語学、心理学、統計学、計算機科学、
そして他の多くの分野から、同程度かそれ以上に得られている。
要約
----
この節では、伝統的な線形回帰を導入した。
そこでは、線形関数のパラメータを
訓練集合上の二乗損失を最小化するように選ぶ。
また、この目的関数の選択を、 いくつかの実用的な考慮と、
線形性およびガウスノイズの仮定の下での
線形回帰の最尤推定としての解釈の両方から動機づけた。
計算上の考慮と統計との関係の両方を議論した後、 このような線形モデルを、
入力が出力に直接配線された単純なニューラルネットワークとして
どのように表せるかを示した。
まもなく線形モデルからは完全に先へ進みるが、
それらは、すべてのモデルに必要なほとんどの構成要素
――パラメトリックな形、微分可能な目的関数、
ミニバッチ確率的勾配降下法による最適化、
そして最終的には、これまで見たことのないデータでの評価――
を導入するのに十分である。
演習
----
1. いくつかのデータ :math:`x_1, \ldots, x_n \in \mathbb{R}`
があると仮定する。私たちの目標は、\ :math:`\sum_i (x_i - b)^2`
が最小になる定数 :math:`b` を見つけることである。
1. :math:`b` の最適値の解析解を求めなさい。
2. この問題とその解は正規分布とどのように関係しているか?
3. 損失を :math:`\sum_i (x_i - b)^2` から :math:`\sum_i |x_i-b|`
に変えたらどうなるか? :math:`b` の最適解を見つけられますか?
2. :math:`\mathbf{x}^\top \mathbf{w} + b`
で表せるアフィン関数が、\ :math:`(\mathbf{x}, 1)`
上の線形関数と等価であることを証明しなさい。
3. :math:`\mathbf{x}` の二次関数、すなわち
:math:`f(\mathbf{x}) = b + \sum_i w_i x_i + \sum_{j \leq i} w_{ij} x_{i} x_{j}`
を求めたいとする。これを深層ネットワークでどのように定式化するか?
4. 線形回帰問題が解ける条件の1つは、設計行列
:math:`\mathbf{X}^\top \mathbf{X}`
がフルランクであることだったことを思い出しよう。
1. そうでない場合はどうなるか?
2. どうすれば修正できるか? :math:`\mathbf{X}`
のすべての要素に、座標ごとに独立なガウスノイズを少量加えるとどうなるか?
3. この場合、設計行列 :math:`\mathbf{X}^\top \mathbf{X}`
の期待値はいくらですか?
4. :math:`\mathbf{X}^\top \mathbf{X}`
がフルランクでないとき、確率的勾配降下法では何が起こりますか?
5. 加法ノイズ :math:`\epsilon`
を支配するノイズモデルが指数分布であると仮定する。すなわち、\ :math:`p(\epsilon) = \frac{1}{2} \exp(-|\epsilon|)`
である。
1. このモデルの下でのデータの負の対数尤度
:math:`-\log P(\mathbf y \mid \mathbf X)` を書き下しなさい。
2. 閉形式解を見つけられますか?
3. この問題を解くためのミニバッチ確率的勾配降下法アルゴリズムを提案しなさい。何がうまくいかなくなる可能性があるか(ヒント:
パラメータを更新し続けると、停留点の近くで何が起こりますか)?
それを修正できるか?
6. 2つの線形層を合成して2層のニューラルネットワークを設計したいとする。すなわち、最初の層の出力が2番目の層の入力になる。なぜこのような素朴な合成はうまくいかないのだろうか?
7. 現実的な住宅価格や株価の推定に回帰を使いたいとすると、何が起こりますか?
1. 加法的ガウスノイズの仮定が適切でないことを示しなさい。ヒント:
価格が負になることはあるか? 変動はどうだろうか?
2. 価格の対数に回帰する方が、つまり :math:`y = \log \textrm{price}`
とする方が、なぜずっと良いのだろうか?
3. ペニー株、すなわち非常に低価格の株を扱うとき、何に注意する必要があるか?
ヒント: すべての可能な価格で取引できるか?
なぜ安い株の方がこの問題は大きいのだろうか?
詳しくは、オプション価格付けの有名な Black–Scholes
モデルを復習せよ :cite:`Black.Scholes.1973`\ 。
8. 食料品店で売れたリンゴの *個数*
を推定するために回帰を使いたいとする。
1. ガウス加法ノイズモデルの問題点は何ですか? ヒント:
売っているのは石油ではなくリンゴである。
2. `ポアソン分布