.. _sec_eigendecompositions:
固有分解
========
固有値は、線形代数を学ぶ際にしばしば最も有用な概念の一つである。
しかし、初心者にとってはその重要性を見落としやすいものでもある。
以下では、固有分解を導入し、
それがなぜこれほど重要なのかを少しでも伝えたいと思いる。
次の成分をもつ行列 :math:`A` があるとする。
.. math::
\mathbf{A} = \begin{bmatrix}
2 & 0 \\
0 & -1
\end{bmatrix}.
:math:`A` を任意のベクトル :math:`\mathbf{v} = [x, y]^\top`
に適用すると、 :math:`\mathbf{A}\mathbf{v} = [2x, -y]^\top`
というベクトルが得られる。 これは直感的には、 ベクトルを :math:`x`
方向に2倍に引き伸ばし、 その後 :math:`y`
方向に反転させる、と解釈できる。
しかし、\ *ある種の*\ ベクトルでは、何かが変わらないまま残る。
すなわち、\ :math:`[1, 0]^\top` は :math:`[2, 0]^\top` に写され、
:math:`[0, 1]^\top` は :math:`[0, -1]^\top` に写される。
これらのベクトルは依然として同じ直線上にあり、
変化するのは行列がそれらをそれぞれ :math:`2` と :math:`-1`
の倍率で伸縮することだけである。
このようなベクトルを\ *固有ベクトル*\ 、
その伸縮の倍率を\ *固有値*\ と呼ぶ。
一般に、ある数 :math:`\lambda` とベクトル :math:`\mathbf{v}` が存在して
.. math::
\mathbf{A}\mathbf{v} = \lambda \mathbf{v}.
となるとき、\ :math:`\mathbf{v}` は :math:`A`
の固有ベクトルであり、\ :math:`\lambda` は固有値であると言う。
固有値の求め方
--------------
それらをどのように求めるかを考えてみよう。両辺から
:math:`\lambda \mathbf{v}` を引き、
さらにベクトルをくくり出すと、上式は次と同値であることがわかる。
.. math:: (\mathbf{A} - \lambda \mathbf{I})\mathbf{v} = 0.
:label: eq_eigvalue_der
:eq:`eq_eigvalue_der`
が成り立つためには、\ :math:`(\mathbf{A} - \lambda \mathbf{I})`
がある方向をゼロまで押しつぶす必要がある。
したがってこれは可逆ではなく、ゆえに行列式はゼロである。
したがって、\ :math:`\det(\mathbf{A}-\lambda \mathbf{I}) = 0` となる
:math:`\lambda` を求めることで、 *固有値*\ を見つけることができる。
固有値が求まれば、 :math:`\mathbf{A}\mathbf{v} = \lambda \mathbf{v}`
を解くことで、対応する\ *固有ベクトル*\ を求められる。
例
~~
もう少し難しい行列で見てみよう。
.. math::
\mathbf{A} = \begin{bmatrix}
2 & 1\\
2 & 3
\end{bmatrix}.
:math:`\det(\mathbf{A}-\lambda \mathbf{I}) = 0` を考えると、
これは多項式方程式
:math:`0 = (2-\lambda)(3-\lambda)-2 = (4-\lambda)(1-\lambda)`
と同値であることがわかる。 したがって、固有値は :math:`4` と :math:`1`
の2つである。 対応するベクトルを求めるには、次を解く必要がある。
.. math::
\begin{bmatrix}
2 & 1\\
2 & 3
\end{bmatrix}\begin{bmatrix}x \\ y\end{bmatrix} = \begin{bmatrix}x \\ y\end{bmatrix} \; \textrm{and} \;
\begin{bmatrix}
2 & 1\\
2 & 3
\end{bmatrix}\begin{bmatrix}x \\ y\end{bmatrix} = \begin{bmatrix}4x \\ 4y\end{bmatrix} .
それぞれ、ベクトル :math:`[1, -1]^\top` と :math:`[1, 2]^\top`
で解ける。
これを組み込みの ``numpy.linalg.eig``
ルーチンを使ってコードで確認できる。
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
%matplotlib inline
from d2l import torch as d2l
from IPython import display
import torch
torch.linalg.eig(torch.tensor([[2, 1], [2, 3]], dtype=torch.float64))
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
torch.return_types.linalg_eig(
eigenvalues=tensor([1.+0.j, 4.+0.j], dtype=torch.complex128),
eigenvectors=tensor([[-0.7071+0.j, -0.4472+0.j],
[ 0.7071+0.j, -0.8944+0.j]], dtype=torch.complex128))
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
%matplotlib inline
from d2l import mxnet as d2l
from IPython import display
import numpy as np
np.linalg.eig(np.array([[2, 1], [2, 3]]))
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
(array([1., 4.]),
array([[-0.70710678, -0.4472136 ],
[ 0.70710678, -0.89442719]]))
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
%matplotlib inline
from d2l import tensorflow as d2l
from IPython import display
import tensorflow as tf
tf.linalg.eig(tf.constant([[2, 1], [2, 3]], dtype=tf.float64))
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
(,
)
.. raw:: html
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
A = torch.tensor([[1.0, 0.1, 0.1, 0.1],
[0.1, 3.0, 0.2, 0.3],
[0.1, 0.2, 5.0, 0.5],
[0.1, 0.3, 0.5, 9.0]])
v, _ = torch.linalg.eig(A)
v
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
tensor([0.9923+0.j, 9.0803+0.j, 4.9539+0.j, 2.9734+0.j])
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
A = np.array([[1.0, 0.1, 0.1, 0.1],
[0.1, 3.0, 0.2, 0.3],
[0.1, 0.2, 5.0, 0.5],
[0.1, 0.3, 0.5, 9.0]])
v, _ = np.linalg.eig(A)
v
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
array([9.08033648, 0.99228545, 4.95394089, 2.97343718])
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
A = tf.constant([[1.0, 0.1, 0.1, 0.1],
[0.1, 3.0, 0.2, 0.3],
[0.1, 0.2, 5.0, 0.5],
[0.1, 0.3, 0.5, 9.0]])
v, _ = tf.linalg.eigh(A)
v
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
.. raw:: html
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
torch.manual_seed(42)
k = 5
A = torch.randn(k, k, dtype=torch.float64)
A
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
tensor([[ 0.2996, 0.2424, 0.2832, -0.2329, 0.6712],
[ 0.7818, -1.7903, -1.7484, 0.1735, -0.1182],
[-1.7446, -0.4695, 0.4573, 0.5177, -0.2771],
[-0.6641, 0.6551, 0.2616, -1.5265, -0.3311],
[-0.6378, 0.1072, 0.7096, 0.3009, -0.2869]], dtype=torch.float64)
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
np.random.seed(8675309)
k = 5
A = np.random.randn(k, k)
A
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
array([[ 0.58902366, 0.73311856, -1.1621888 , -0.55681601, -0.77248843],
[-0.16822143, -0.41650391, -1.37843129, 0.74925588, 0.17888446],
[ 0.69401121, -1.9780535 , -0.83381434, 0.56437344, 0.31201299],
[-0.87334496, 0.15601291, -0.38710108, -0.23920821, 0.88850104],
[ 1.29385371, -0.76774106, 0.20131613, 0.91800842, 0.38974115]])
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
k = 5
A = tf.random.normal((k, k), dtype=tf.float64)
A
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
.. raw:: html
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Calculate the sequence of norms after repeatedly applying `A`
v_in = torch.randn(k, 1, dtype=torch.float64)
norm_list = [torch.norm(v_in).item()]
for i in range(1, 100):
v_in = A @ v_in
norm_list.append(torch.norm(v_in).item())
d2l.plot(torch.arange(0, 100), norm_list, 'Iteration', 'Value')
.. figure:: output_eigendecomposition_89942e_39_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Calculate the sequence of norms after repeatedly applying `A`
v_in = np.random.randn(k, 1)
norm_list = [np.linalg.norm(v_in)]
for i in range(1, 100):
v_in = A.dot(v_in)
norm_list.append(np.linalg.norm(v_in))
d2l.plot(np.arange(0, 100), norm_list, 'Iteration', 'Value')
.. figure:: output_eigendecomposition_89942e_42_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Calculate the sequence of norms after repeatedly applying `A`
v_in = tf.random.normal((k, 1), dtype=tf.float64)
norm_list = [tf.norm(v_in).numpy()]
for i in range(1, 100):
v_in = tf.matmul(A, v_in)
norm_list.append(tf.norm(v_in).numpy())
d2l.plot(tf.range(0, 100), norm_list, 'Iteration', 'Value')
.. figure:: output_eigendecomposition_89942e_45_0.svg
.. raw:: html
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Compute the scaling factor of the norms
norm_ratio_list = []
for i in range(1, 100):
norm_ratio_list.append(norm_list[i]/norm_list[i - 1])
d2l.plot(torch.arange(1, 100), norm_ratio_list, 'Iteration', 'Ratio')
.. figure:: output_eigendecomposition_89942e_51_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Compute the scaling factor of the norms
norm_ratio_list = []
for i in range(1, 100):
norm_ratio_list.append(norm_list[i]/norm_list[i - 1])
d2l.plot(np.arange(1, 100), norm_ratio_list, 'Iteration', 'Ratio')
.. figure:: output_eigendecomposition_89942e_54_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Compute the scaling factor of the norms
norm_ratio_list = []
for i in range(1, 100):
norm_ratio_list.append(norm_list[i]/norm_list[i - 1])
d2l.plot(tf.range(1, 100), norm_ratio_list, 'Iteration', 'Ratio')
.. figure:: output_eigendecomposition_89942e_57_0.svg
.. raw:: html
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Compute the eigenvalues
eigs = torch.linalg.eig(A).eigenvalues.tolist()
norm_eigs = [torch.abs(torch.tensor(x)) for x in eigs]
norm_eigs.sort()
print(f'norms of eigenvalues: {norm_eigs}')
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
norms of eigenvalues: [tensor(0.3490), tensor(1.1296), tensor(1.1296), tensor(1.1828), tensor(2.4532)]
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Compute the eigenvalues
eigs = np.linalg.eigvals(A).tolist()
norm_eigs = [np.absolute(x) for x in eigs]
norm_eigs.sort()
print(f'norms of eigenvalues: {norm_eigs}')
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
norms of eigenvalues: [0.8786205280381857, 1.2757952665062624, 1.4983381517710659, 1.4983381517710659, 1.974459321485074]
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Compute the eigenvalues
eigs = tf.linalg.eigh(A)[0].numpy().tolist()
norm_eigs = [tf.abs(tf.constant(x, dtype=tf.float64)) for x in eigs]
norm_eigs.sort()
print(f'norms of eigenvalues: {norm_eigs}')
.. raw:: latex
\diilbookstyleoutputcell
.. parsed-literal::
:class: output
norms of eigenvalues: [, , , , ]
.. raw:: html
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Rescale the matrix `A`
A /= norm_eigs[-1]
# Do the same experiment again
v_in = torch.randn(k, 1, dtype=torch.float64)
norm_list = [torch.norm(v_in).item()]
for i in range(1, 100):
v_in = A @ v_in
norm_list.append(torch.norm(v_in).item())
d2l.plot(torch.arange(0, 100), norm_list, 'Iteration', 'Value')
.. figure:: output_eigendecomposition_89942e_75_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Rescale the matrix `A`
A /= norm_eigs[-1]
# Do the same experiment again
v_in = np.random.randn(k, 1)
norm_list = [np.linalg.norm(v_in)]
for i in range(1, 100):
v_in = A.dot(v_in)
norm_list.append(np.linalg.norm(v_in))
d2l.plot(np.arange(0, 100), norm_list, 'Iteration', 'Value')
.. figure:: output_eigendecomposition_89942e_78_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Rescale the matrix `A`
A /= norm_eigs[-1]
# Do the same experiment again
v_in = tf.random.normal((k, 1), dtype=tf.float64)
norm_list = [tf.norm(v_in).numpy()]
for i in range(1, 100):
v_in = tf.matmul(A, v_in)
norm_list.append(tf.norm(v_in).numpy())
d2l.plot(tf.range(0, 100), norm_list, 'Iteration', 'Value')
.. figure:: output_eigendecomposition_89942e_81_0.svg
.. raw:: html
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Also plot the ratio
norm_ratio_list = []
for i in range(1, 100):
norm_ratio_list.append(norm_list[i]/norm_list[i-1])
d2l.plot(torch.arange(1, 100), norm_ratio_list, 'Iteration', 'Ratio')
.. figure:: output_eigendecomposition_89942e_87_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Also plot the ratio
norm_ratio_list = []
for i in range(1, 100):
norm_ratio_list.append(norm_list[i]/norm_list[i-1])
d2l.plot(np.arange(1, 100), norm_ratio_list, 'Iteration', 'Ratio')
.. figure:: output_eigendecomposition_89942e_90_0.svg
.. raw:: html
.. raw:: html
.. raw:: latex
\diilbookstyleinputcell
.. code:: python
# Also plot the ratio
norm_ratio_list = []
for i in range(1, 100):
norm_ratio_list.append(norm_list[i]/norm_list[i-1])
d2l.plot(tf.range(1, 100), norm_ratio_list, 'Iteration', 'Ratio')
.. figure:: output_eigendecomposition_89942e_93_0.svg
.. raw:: html
.. raw:: html