グラム・シュミットの正規直交化法（Gram-Schmidt orthonormalization）は、与えられた一次独立なベクトルの組から正規直交基底を構成する手法です。線形代数の理論と実践の両面で非常に重要な役割を果たします。
正規直交基底は、フーリエ解析、QR分解、最小二乗法など、多くの分野で基礎となる概念です。グラム・シュミットの正規直交化法はその構成方法として最も基本的なアルゴリズムです。
本記事の内容
- 直交射影の復習
- グラム・シュミットの正規直交化法のアルゴリズム
- 具体的な計算例
- QR分解との関係
- Pythonでの実装（古典版と修正版）
前提知識
この記事を読む前に、以下の概念を理解しておくと理解が深まります。
- 内積とノルムの定義
- 一次独立と基底の概念
正規直交基底とは
ベクトルの組 $\{\bm{e}_1, \bm{e}_2, \dots, \bm{e}_n\}$ が正規直交系であるとは、
$$ \bm{e}_i \cdot \bm{e}_j = \delta_{ij} = \begin{cases} 1 & (i = j) \\ 0 & (i \neq j) \end{cases} $$
を満たすことをいいます。つまり、各ベクトルの大きさが1（正規）で、異なるベクトル同士が直交（直交）しています。
正規直交基底を使うと、ベクトルの成分表示が内積で簡単に求まります。
$$ \bm{v} = \sum_{i=1}^{n} (\bm{v} \cdot \bm{e}_i) \bm{e}_i $$
直交射影
グラム・シュミットの正規直交化法の核となるのが直交射影です。
ベクトル $\bm{v}$ のベクトル $\bm{u}$（$\bm{u} \neq \bm{0}$）への直交射影は、
$$ \mathrm{proj}_{\bm{u}} \bm{v} = \frac{\bm{v} \cdot \bm{u}}{\bm{u} \cdot \bm{u}} \bm{u} $$
で定義されます。この射影ベクトルは $\bm{u}$ 方向の成分を表し、$\bm{v} – \mathrm{proj}_{\bm{u}} \bm{v}$ は $\bm{u}$ に直交する成分です。
グラム・シュミットのアルゴリズム
手順
一次独立なベクトルの組 $\{\bm{v}_1, \bm{v}_2, \dots, \bm{v}_n\}$ から正規直交系 $\{\bm{e}_1, \bm{e}_2, \dots, \bm{e}_n\}$ を構成します。
ステップ1: 最初のベクトルを正規化
$$ \bm{e}_1 = \frac{\bm{v}_1}{\|\bm{v}_1\|} $$
ステップ2: 2番目のベクトルを直交化して正規化
$$ \bm{u}_2 = \bm{v}_2 – (\bm{v}_2 \cdot \bm{e}_1) \bm{e}_1 $$
$$ \bm{e}_2 = \frac{\bm{u}_2}{\|\bm{u}_2\|} $$
ステップk: $k$ 番目のベクトルを直交化して正規化
$$ \bm{u}_k = \bm{v}_k – \sum_{j=1}^{k-1} (\bm{v}_k \cdot \bm{e}_j) \bm{e}_j $$
$$ \bm{e}_k = \frac{\bm{u}_k}{\|\bm{u}_k\|} $$
アルゴリズムの意味
各ステップで行っていることは、「新しいベクトル $\bm{v}_k$ から、既に得られた正規直交ベクトルの方向の成分をすべて引き去る」ことです。残った成分 $\bm{u}_k$ は、それまでのベクトルすべてに直交します。
具体例
$\bm{v}_1 = (1, 1, 0)^T$, $\bm{v}_2 = (1, 0, 1)^T$, $\bm{v}_3 = (0, 1, 1)^T$ を正規直交化します。
ステップ1:
$$ \|\bm{v}_1\| = \sqrt{1^2 + 1^2 + 0^2} = \sqrt{2} $$
$$ \bm{e}_1 = \frac{1}{\sqrt{2}}(1, 1, 0)^T $$
ステップ2:
$$ \bm{v}_2 \cdot \bm{e}_1 = \frac{1}{\sqrt{2}}(1 \cdot 1 + 0 \cdot 1 + 1 \cdot 0) = \frac{1}{\sqrt{2}} $$
$$ \bm{u}_2 = \bm{v}_2 – \frac{1}{\sqrt{2}} \bm{e}_1 = (1, 0, 1)^T – \frac{1}{2}(1, 1, 0)^T = \left(\frac{1}{2}, -\frac{1}{2}, 1\right)^T $$
$$ \|\bm{u}_2\| = \sqrt{\frac{1}{4} + \frac{1}{4} + 1} = \sqrt{\frac{3}{2}} = \frac{\sqrt{6}}{2} $$
$$ \bm{e}_2 = \frac{2}{\sqrt{6}}\left(\frac{1}{2}, -\frac{1}{2}, 1\right)^T = \frac{1}{\sqrt{6}}(1, -1, 2)^T $$
ステップ3:
$$ \bm{v}_3 \cdot \bm{e}_1 = \frac{1}{\sqrt{2}}(0 + 1 + 0) = \frac{1}{\sqrt{2}} $$
$$ \bm{v}_3 \cdot \bm{e}_2 = \frac{1}{\sqrt{6}}(0 – 1 + 2) = \frac{1}{\sqrt{6}} $$
$$ \bm{u}_3 = (0, 1, 1)^T – \frac{1}{\sqrt{2}} \cdot \frac{1}{\sqrt{2}}(1, 1, 0)^T – \frac{1}{\sqrt{6}} \cdot \frac{1}{\sqrt{6}}(1, -1, 2)^T $$
$$ = (0, 1, 1)^T – \frac{1}{2}(1, 1, 0)^T – \frac{1}{6}(1, -1, 2)^T = \left(-\frac{2}{3}, \frac{2}{3}, \frac{2}{3}\right)^T $$
$$ \|\bm{u}_3\| = \frac{2}{3}\sqrt{3} = \frac{2\sqrt{3}}{3} $$
$$ \bm{e}_3 = \frac{1}{\sqrt{3}}(-1, 1, 1)^T $$
QR分解との関係
グラム・シュミットの正規直交化法は、QR分解と密接に関係しています。
行列 $\bm{A}$ の列ベクトルを $\bm{v}_1, \dots, \bm{v}_n$ とし、グラム・シュミット法で得られた正規直交ベクトルを $\bm{e}_1, \dots, \bm{e}_n$ とすると、
$$ \bm{A} = \bm{Q}\bm{R} $$
と分解できます。ここで $\bm{Q} = (\bm{e}_1, \dots, \bm{e}_n)$ は直交行列、$\bm{R}$ は上三角行列で、$r_{jk} = \bm{v}_k \cdot \bm{e}_j$（$j \leq k$）です。
Pythonでの実装
古典的グラム・シュミット法
import numpy as np
def gram_schmidt_classical(V):
"""古典的グラム・シュミット正規直交化法
V: 各列がベクトルの行列 (n x m)
"""
n, m = V.shape
Q = np.zeros((n, m))
R = np.zeros((m, m))
for k in range(m):
# 直交化
u = V[:, k].copy()
for j in range(k):
R[j, k] = np.dot(Q[:, j], V[:, k])
u -= R[j, k] * Q[:, j]
# 正規化
R[k, k] = np.linalg.norm(u)
if R[k, k] < 1e-10:
raise ValueError("ベクトルが一次従属です")
Q[:, k] = u / R[k, k]
return Q, R
# テスト
V = np.array([[1, 1, 0],
[1, 0, 1],
[0, 1, 1]], dtype=float).T # 列ベクトルとして格納
# 転置して各列がベクトルになるようにする
A = np.array([[1, 1, 0],
[1, 0, 1],
[0, 1, 1]], dtype=float)
Q, R = gram_schmidt_classical(A)
print("Q (正規直交行列):")
print(np.round(Q, 4))
print("\nR (上三角行列):")
print(np.round(R, 4))
print("\n直交性の確認 Q^T Q:")
print(np.round(Q.T @ Q, 10))
print("\nA = QR の確認:")
print(np.round(Q @ R, 4))
print("元の行列 A:")
print(A)
修正グラム・シュミット法
古典的な方法は数値的に不安定になることがあります。修正版（Modified Gram-Schmidt）はより安定です。
import numpy as np
def gram_schmidt_modified(V):
"""修正グラム・シュミット正規直交化法"""
n, m = V.shape
Q = np.zeros((n, m))
R = np.zeros((m, m))
W = V.copy().astype(float)
for k in range(m):
R[k, k] = np.linalg.norm(W[:, k])
if R[k, k] < 1e-10:
raise ValueError("ベクトルが一次従属です")
Q[:, k] = W[:, k] / R[k, k]
for j in range(k + 1, m):
R[k, j] = np.dot(Q[:, k], W[:, j])
W[:, j] -= R[k, j] * Q[:, k]
return Q, R
# 数値安定性の比較
np.random.seed(42)
n = 50
# 条件数の悪い行列を生成
U, _ = np.linalg.qr(np.random.randn(n, n))
S = np.diag(np.logspace(0, -10, n))
A_ill = U @ S @ U.T
A_test = A_ill[:, :10]
Q_classical, R_classical = gram_schmidt_classical(A_test)
Q_modified, R_modified = gram_schmidt_modified(A_test)
Q_numpy, R_numpy = np.linalg.qr(A_test)
# 直交性の誤差を比較
err_classical = np.linalg.norm(Q_classical.T @ Q_classical - np.eye(10))
err_modified = np.linalg.norm(Q_modified.T @ Q_modified - np.eye(10))
err_numpy = np.linalg.norm(Q_numpy.T @ Q_numpy - np.eye(10))
print(f"直交性誤差（古典版）: {err_classical:.2e}")
print(f"直交性誤差（修正版）: {err_modified:.2e}")
print(f"直交性誤差（NumPy）: {err_numpy:.2e}")
可視化: 3次元での正規直交化
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
v1 = np.array([1, 1, 0])
v2 = np.array([1, 0, 1])
v3 = np.array([0, 1, 1])
A = np.column_stack([v1, v2, v3])
Q, R = gram_schmidt_classical(A)
e1, e2, e3 = Q[:, 0], Q[:, 1], Q[:, 2]
fig = plt.figure(figsize=(14, 6))
# 元のベクトル
ax1 = fig.add_subplot(121, projection='3d')
origin = np.zeros(3)
for v, label, color in zip([v1, v2, v3], ['$v_1$', '$v_2$', '$v_3$'], ['r', 'g', 'b']):
ax1.quiver(*origin, *v, color=color, arrow_length_ratio=0.1, linewidth=2)
ax1.text(*(v * 1.1), label, fontsize=14)
ax1.set_xlim([-0.5, 1.5])
ax1.set_ylim([-0.5, 1.5])
ax1.set_zlim([-0.5, 1.5])
ax1.set_xlabel('$x$')
ax1.set_ylabel('$y$')
ax1.set_zlabel('$z$')
ax1.set_title('Original vectors')
# 正規直交化後のベクトル
ax2 = fig.add_subplot(122, projection='3d')
for e, label, color in zip([e1, e2, e3], ['$e_1$', '$e_2$', '$e_3$'], ['r', 'g', 'b']):
ax2.quiver(*origin, *e, color=color, arrow_length_ratio=0.1, linewidth=2)
ax2.text(*(e * 1.15), label, fontsize=14)
ax2.set_xlim([-1, 1])
ax2.set_ylim([-1, 1])
ax2.set_zlim([-1, 1])
ax2.set_xlabel('$x$')
ax2.set_ylabel('$y$')
ax2.set_zlabel('$z$')
ax2.set_title('Orthonormalized vectors')
plt.tight_layout()
plt.savefig('gram_schmidt_3d.png', dpi=150, bbox_inches='tight')
plt.show()
まとめ
本記事では、グラム・シュミットの正規直交化法を解説し、Pythonで実装しました。
- 正規直交基底: 互いに直交し、各ベクトルのノルムが1
- アルゴリズム: 各ベクトルから既存の正規直交ベクトルの成分を引き去って直交化し、正規化する
- QR分解: グラム・シュミット法は $\bm{A} = \bm{Q}\bm{R}$ という分解を与える
- 修正版: 古典版よりも数値的に安定
次のステップとして、以下の記事も参考にしてください。