Lecture 11. 행렬공간, 랭크 1 행렬, 그래프 — 파이썬 실습 - 선형대수

L11 서술 파트에서 랭크 1 행렬과 행렬공간, 그리고 인접행렬을 다루었다. 이 노트북에서는 그것들을 직접 만들어 보도록 하자.

서술 파트의 내용	여기서 확인하는 방법
행렬도 벡터공간이다	행렬을 벡터로 펴서 랭크를 재면 차원이 나온다
랭크 1은 공간을 직선으로 뭉갠다	무작위 점을 통과시켜 회전해 본다
랭크 $r$ $=$ 랭크 1 조각 $r$ 개	쪼갠 조각을 하나씩 더하며 랭크를 센다
함수도 벡터다	다항식 공간에서 미분을 행렬로 쓴다
$(A^k)_{ij}$ 는 $k$ 걸음 경로의 개수	작은 그래프와 사회 연결망으로 확인한다

0. 준비¶

import numpy as np
import plotly.graph_objects as go
from scipy.linalg import null_space

from linalg_viz import COLORS, arrow, spin_figure, show_matrix

np.set_printoptions(precision=3, suppress=True)
rng = np.random.default_rng(0)
print("numpy", np.__version__)

numpy 2.5.2

1. 행렬도 벡터공간이다¶

행렬을 한 줄로 펴면 그냥 벡터가 된다. ravel() 이 그 일을 한다. 그러면 지금까지 쓰던 도구를 행렬의 집합에도 그대로 쓸 수 있다.

def 부분공간차원(행렬들):
    """행렬들을 벡터로 펴서 랭크를 재면 그들이 생성하는 공간의 차원이 나온다."""
    return np.linalg.matrix_rank(np.column_stack([M.ravel() for M in 행렬들]))


def E(i, j, n=3):
    """(i, j) 자리만 1인 행렬."""
    M = np.zeros((n, n))
    M[i, j] = 1.0
    return M

전체기저 = [E(i, j) for i in range(3) for j in range(3)]
대칭기저 = [E(i, i) for i in range(3)] + \
          [E(i, j) + E(j, i) for i in range(3) for j in range(i + 1, 3)]
위삼각기저 = [E(i, j) for i in range(3) for j in range(i, 3)]
대각기저 = [E(i, i) for i in range(3)]

for 이름, 기저 in [("3x3 행렬 전체", 전체기저), ("대칭", 대칭기저),
                  ("위삼각", 위삼각기저), ("대각", 대각기저)]:
    print(f"{이름:<14} 생성 벡터 {len(기저):>2}개,  차원 {부분공간차원(기저)}")

3x3 행렬 전체      생성 벡터  9개,  차원 9
대칭             생성 벡터  6개,  차원 6
위삼각            생성 벡터  6개,  차원 6
대각             생성 벡터  3개,  차원 3

서술 파트의 표와 같다. 합의 차원도 확인해 보자. 대칭 기저와 위삼각 기저를 모두 모아 랭크를 재면 된다.

합차원 = 부분공간차원(대칭기저 + 위삼각기저)
print(f"dim(대칭) + dim(위삼각) - dim(교집합) = 6 + 6 - 3 = {6 + 6 - 3}")
print(f"실제로 잰 dim(합)                     = {합차원}")
print(f"3x3 행렬 전체의 차원                  = {부분공간차원(전체기저)}")

dim(대칭) + dim(위삼각) - dim(교집합) = 6 + 6 - 3 = 9
실제로 잰 dim(합)                     = 9
3x3 행렬 전체의 차원                  = 9

합이 전체와 같으므로 어떤 행렬이든 대칭행렬과 위삼각행렬의 합으로 쓸 수 있어야 한다. 서술 파트의 방법대로 직접 만들어 보자. 대각선 아래만 맞추면 된다.

M = rng.standard_normal((3, 3))

S = np.zeros((3, 3))
for i in range(3):
    for j in range(i):                 # 대각선 아래
        S[i, j] = S[j, i] = M[i, j]
U = M - S

print(show_matrix(M, "M ="))
print(show_matrix(S, "S (대칭) ="))
print(show_matrix(U, "U (위삼각) ="))
print()
print("S 가 대칭인가        :", np.allclose(S, S.T))
print("U 가 위삼각인가      :", np.allclose(U, np.triu(U)))
print("S + U 가 M 인가      :", np.allclose(S + U, M))

M =
[   0.126   -0.132     0.64 ]
[   0.105   -0.536    0.362 ]
[     1.3    0.947   -0.704 ]
S (대칭) =
[      0   0.105     1.3 ]
[  0.105       0   0.947 ]
[    1.3   0.947       0 ]
U (위삼각) =
[   0.126   -0.237   -0.664 ]
[       0   -0.536   -0.585 ]
[       0        0   -0.704 ]

S 가 대칭인가        : True
U 가 위삼각인가      : True
S + U 가 M 인가      : True

$S$ 의 대각 성분을 아무 값으로 바꿔도 여전히 성립한다. 쪼개는 방법이 하나가 아니고, 그 자유도가 교집합의 차원 3과 같다.

S2 = S.copy()
S2[np.diag_indices(3)] = [7.0, -2.0, 0.5]      # 대각을 마음대로 바꾼다
U2 = M - S2

print("S2 가 대칭인가   :", np.allclose(S2, S2.T))
print("U2 가 위삼각인가 :", np.allclose(U2, np.triu(U2)))
print("S2 + U2 = M 인가 :", np.allclose(S2 + U2, M))

S2 가 대칭인가   : True
U2 가 위삼각인가 : True
S2 + U2 = M 인가 : True

2. 랭크 1 행렬¶

$\vv{u}\vv{v}^{\mathsf{T}}$ 를 만들어 보자. np.outer 가 이 계산을 한다.

u = np.array([1.0, 2.0, 1.0])
v = np.array([1.0, -1.0, 2.0])

A1 = np.outer(u, v)
print(show_matrix(A1, "A = u v^T ="))
print("rank :", np.linalg.matrix_rank(A1))
print()
print("u v^T 의 크기 :", A1.shape, " (행렬)")
print("u^T v 의 값   :", u @ v, " (숫자 하나)")

A = u v^T =
[   1   -1    2 ]
[   2   -2    4 ]
[   1   -1    2 ]
rank : 1

u v^T 의 크기 : (3, 3)  (행렬)
u^T v 의 값   : 1.0  (숫자 하나)

모든 열이 $\vv{u}$ 의 배수이고, 그 배수가 $\vv{v}$ 의 성분이다.

for j in range(3):
    print(f"열 {j + 1} = {A1[:, j]}  =  {v[j]:g} * u  =  {v[j] * u}")
print()
for i in range(3):
    print(f"행 {i + 1} = {A1[i]}  =  {u[i]:g} * v^T")

열 1 = [1. 2. 1.]  =  1 * u  =  [1. 2. 1.]
열 2 = [-1. -2. -1.]  =  -1 * u  =  [-1. -2. -1.]
열 3 = [2. 4. 2.]  =  2 * u  =  [2. 4. 2.]

행 1 = [ 1. -1.  2.]  =  1 * v^T
행 2 = [ 2. -2.  4.]  =  2 * v^T
행 3 = [ 1. -1.  2.]  =  1 * v^T

$A\vv{x}$ 는 언제나 $\vv{u}$ 의 배수이고, 그 배수는 $\vv{v} \cdot \vv{x}$ 이다.

for _ in range(4):
    x = rng.standard_normal(3)
    좌 = A1 @ x
    우 = (v @ x) * u
    print(f"A x = {좌}   (v . x) u = {우}   같은가 {np.allclose(좌, 우)}")

A x = [-0.559 -1.119 -0.559]   (v . x) u = [-0.559 -1.119 -0.559]   같은가 True
A x = [-4.598 -9.196 -4.598]   (v . x) u = [-4.598 -9.196 -4.598]   같은가 True
A x = [-0.821 -1.641 -0.821]   (v . x) u = [-0.821 -1.641 -0.821]   같은가 True
A x = [-0.888 -1.776 -0.888]   (v . x) u = [-0.888 -1.776 -0.888]   같은가 True

N = null_space(A1)
print("영공간의 차원 :", N.shape[1], "  (n - r = 3 - 1 = 2)")
print("영공간 기저가 v 와 수직인가 :", np.allclose(v @ N, 0))

영공간의 차원 : 2   (n - r = 3 - 1 = 2)
영공간 기저가 v 와 수직인가 : True

무작위 점 200개를 통과시켜 어디로 가는지 보자.

표본 = rng.standard_normal((3, 200))
상 = A1 @ 표본

d = u / np.linalg.norm(u)
traces = [
    go.Scatter3d(x=상[0], y=상[1], z=상[2], mode="markers",
                 marker=dict(size=3, color=COLORS["output"]), name="A x"),
    go.Scatter3d(x=[-9 * d[0], 9 * d[0]], y=[-9 * d[1], 9 * d[1]],
                 z=[-9 * d[2], 9 * d[2]], mode="lines",
                 line=dict(color=COLORS["nullspace"], width=4, dash="dash"),
                 name="u 방향 직선"),
]
traces += arrow([0, 0, 0], 3 * u, color=COLORS["third"], name="u")

spin_figure(traces, title="랭크 1 을 지나면 전부 한 직선 위에", extent=9)

3. 랭크 $r$ 행렬은 랭크 1 조각 $r$ 개의 합¶

서술 파트의 방법대로 쪼개 보자. 열공간의 기저를 $U$ 로 두고 계수 $C$ 를 구하면 $A = UC$ 이고, L3의 네 번째 관점으로 쪼개면 조각이 $r$ 개 나온다.

def 랭크1로_쪼개기(A):
    """A 를 랭크 1 조각들의 합으로 쪼갠다. 조각의 개수는 랭크와 같다."""
    A = np.asarray(A, dtype=float)
    r = np.linalg.matrix_rank(A)

    # 열공간의 기저를 왼쪽부터 독립인 열을 골라 만든다
    골라낸 = []
    for j in range(A.shape[1]):
        후보 = 골라낸 + [j]
        if np.linalg.matrix_rank(A[:, 후보]) == len(후보):
            골라낸 = 후보
        if len(골라낸) == r:
            break

    U = A[:, 골라낸]                      # m x r
    C, *_ = np.linalg.lstsq(U, A, rcond=None)   # r x n,  A = U C
    return [np.outer(U[:, k], C[k, :]) for k in range(r)]

A = rng.standard_normal((5, 3)) @ rng.standard_normal((3, 6))    # 랭크 3
r = np.linalg.matrix_rank(A)
조각들 = 랭크1로_쪼개기(A)

print("A 의 크기 :", A.shape, "  rank :", r)
print("조각 개수 :", len(조각들))
print("각 조각의 랭크 :", [int(np.linalg.matrix_rank(P)) for P in 조각들])
print("조각을 모두 더하면 A 인가 :", np.allclose(sum(조각들), A))

A 의 크기 : (5, 6)   rank : 3
조각 개수 : 3
각 조각의 랭크 : [1, 1, 1]
조각을 모두 더하면 A 인가 : True

조각을 하나씩 더해 가며 랭크를 세어 보자. 하나 더할 때마다 1씩 오른다.

누적 = np.zeros_like(A)
for k, P in enumerate(조각들, 1):
    누적 = 누적 + P
    남은오차 = np.linalg.norm(A - 누적)
    print(f"조각 {k}개까지 더한 랭크 : {np.linalg.matrix_rank(누적)}"
          f"   |A - 누적| = {남은오차:.3e}")

조각 1개까지 더한 랭크 : 1   |A - 누적| = 2.463e+01
조각 2개까지 더한 랭크 : 2   |A - 누적| = 2.817e+01
조각 3개까지 더한 랭크 : 3   |A - 누적| = 1.382e-14

남은 오차를 보면 이상한 점이 있다. 조각을 하나 더했을 때 오차가 줄기는커녕 늘어나는 경우가 생긴다. 조각을 앞에서부터 몇 개만 남겨도 원래에 가까워지는 것이 아니라는 뜻이다.

# 열 순서를 섞으면 다른 조각이 나온다
섞기 = rng.permutation(A.shape[1])
다른조각들 = 랭크1로_쪼개기(A[:, 섞기])

print("조각 개수는 같은가 :", len(다른조각들) == len(조각들))
print("첫 조각이 같은가   :", np.allclose(다른조각들[0], 조각들[0]))
print("합은 여전히 맞는가 :", np.allclose(sum(다른조각들), A[:, 섞기]))

조각 개수는 같은가 : True
첫 조각이 같은가   : False
합은 여전히 맞는가 : True

4. 함수도 벡터다¶

차수가 2 이하인 다항식 $a_0 + a_1x + a_2x^2$ 을 좌표 $(a_0, a_1, a_2)$ 로 나타내면 $\R^3$ 처럼 다룰 수 있다. 그러면 미분도 행렬이 된다.

\frac{d}{dx}(a_0 + a_1x + a_2x^2) = a_1 + 2a_2 x

(1)

이므로 계수는 $(a_0, a_1, a_2) \mapsto (a_1, 2a_2, 0)$ 으로 간다.

D = np.array([[0.0, 1, 0],
              [0.0, 0, 2],
              [0.0, 0, 0]])

print(show_matrix(D, "미분 연산자의 행렬 D ="))

계수 = np.array([5.0, 3.0, 4.0])          # 5 + 3x + 4x^2
print("p  의 계수 :", 계수)
print("p' 의 계수 :", D @ 계수, "  -> 3 + 8x")

미분 연산자의 행렬 D =
[  0   1   0 ]
[  0   0   2 ]
[  0   0   0 ]
p  의 계수 : [5. 3. 4.]
p' 의 계수 : [3. 8. 0.]   -> 3 + 8x

이 행렬에도 네 부분공간이 있다. 영공간은 미분해서 0이 되는 다항식, 즉 상수함수이다.

print("rank(D)      :", np.linalg.matrix_rank(D))
print("dim N(D)     :", null_space(D).shape[1])
print("합           :", np.linalg.matrix_rank(D) + null_space(D).shape[1], " (n = 3)")
print()
print("N(D) 의 기저 :", null_space(D)[:, 0], " -> 상수함수 (a0 만 남는다)")

rank(D)      : 2
dim N(D)     : 1
합           : 3  (n = 3)

N(D) 의 기저 : [1. 0. 0.]  -> 상수함수 (a0 만 남는다)

미분방정식 $y'' + y = 0$ 의 해공간도 확인해 보자. $\cos x$ 와 $\sin x$ 를 여러 점에서 값을 재어 벡터로 보면 서로 독립이다.

x = np.linspace(0, 2 * np.pi, 200)
F = np.column_stack([np.cos(x), np.sin(x)])

print("표본점에서 본 cos 와 sin 의 랭크 :", np.linalg.matrix_rank(F), " -> 독립, 차원 2")
print()

# y = c1 cos + c2 sin 이 정말 해인지 수치 미분으로 확인
c1, c2 = 2.0, -3.0
h = 1e-4
y = lambda t: c1 * np.cos(t) + c2 * np.sin(t)
t = np.linspace(0.5, 5.5, 50)
y2차 = (y(t + h) - 2 * y(t) + y(t - h)) / h ** 2

print("y'' + y 의 최대 절댓값 :", np.abs(y2차 + y(t)).max(), " (0 이어야 한다)")

표본점에서 본 cos 와 sin 의 랭크 : 2  -> 독립, 차원 2

y'' + y 의 최대 절댓값 : 1.3522990771619448e-07  (0 이어야 한다)

5. 그래프와 인접행렬¶

서술 파트의 다섯 점짜리 그래프를 그대로 만들어 보자.

변 = [(0, 1), (1, 2), (2, 3), (3, 4), (4, 0), (1, 4)]
Adj = np.zeros((5, 5))
for a, b in 변:
    Adj[a, b] = Adj[b, a] = 1.0

print(show_matrix(Adj, "인접행렬 A ="))
print(show_matrix(Adj @ Adj, "A^2 ="))

인접행렬 A =
[  0   1   0   0   1 ]
[  1   0   1   0   1 ]
[  0   1   0   1   0 ]
[  0   0   1   0   1 ]
[  1   1   0   1   0 ]
A^2 =
[  2   1   1   1   1 ]
[  1   3   0   2   1 ]
[  1   0   2   0   2 ]
[  1   2   0   2   0 ]
[  1   1   2   0   3 ]

print("(A^2)[1,1] =", (Adj @ Adj)[1, 1], " -> 2번 점의 이웃 수")
print("2번 점의 이웃 :", [i + 1 for i in range(5) if Adj[1, i] == 1])
print()
print("(A^2)[0,2] =", (Adj @ Adj)[0, 2], " -> 1번에서 3번으로 가는 두 걸음 경로 수")
print("  1 -> 2 -> 3 :", Adj[0, 1] * Adj[1, 2] == 1)
print("  1 -> 5 -> 3 :", Adj[0, 4] * Adj[4, 2] == 1)

(A^2)[1,1] = 3.0  -> 2번 점의 이웃 수
2번 점의 이웃 : [1, 3, 5]

(A^2)[0,2] = 1.0  -> 1번에서 3번으로 가는 두 걸음 경로 수
  1 -> 2 -> 3 : True
  1 -> 5 -> 3 : False

여섯 다리 건너¶

사회 연결망을 흉내 낸 그래프를 만들어, 몇 걸음이면 모든 쌍이 이어지는지 세어 보자. networkx 의 작은 세상 그래프를 쓴다. 이웃 몇 개만 연결된 고리에서 변 일부를 무작위로 바꿔 끼운 것이다.

import networkx as nx


def 도달비율(G, 최대걸음=8):
    """k 걸음 이내로 이어지는 점 쌍의 비율을 센다."""
    A = nx.to_numpy_array(G)
    누적 = np.eye(len(A))
    P = np.eye(len(A))
    for k in range(1, 최대걸음 + 1):
        P = P @ A
        누적 = 누적 + P
        yield k, np.count_nonzero(누적) / 누적.size

고리 = nx.watts_strogatz_graph(200, 6, 0.0, seed=0)      # 바꿔 끼운 변 없음
작은세상 = nx.watts_strogatz_graph(200, 6, 0.1, seed=0)   # 10% 를 바꿔 끼움

print(f"{'걸음':>4}{'고리만':>12}{'작은 세상':>12}")
for (k, a), (_, b) in zip(도달비율(고리), 도달비율(작은세상)):
    print(f"{k:>4}{a:>12.3f}{b:>12.3f}")

  걸음         고리만       작은 세상
   1       0.035       0.035
   2       0.065       0.097
   3       0.095       0.240
   4       0.125       0.501
   5       0.155       0.802
   6       0.185       0.968
   7       0.215       1.000
   8       0.245       1.000

변의 개수는 거의 같은데, 10퍼센트만 바꿔 끼워도 훨씬 적은 걸음으로 모두가 이어진다. 이것이 작은 세상이라는 말의 뜻이다.

여기서 쓴 계산은 전부 $(A^k)_{ij}$ 가 $k$ 걸음 경로의 개수라는 한 가지 사실에 기대고 있다.

for 이름, G in [("고리만", 고리), ("작은 세상", 작은세상)]:
    print(f"{이름:<10} 변 {G.number_of_edges()}개,"
          f"  지름(가장 먼 두 점 사이 거리) {nx.diameter(G)}")

고리만        변 600개,  지름(가장 먼 두 점 사이 거리) 34

작은 세상      변 600개,  지름(가장 먼 두 점 사이 거리) 8

마치며...¶

서술 파트의 내용	이 노트북의 코드
행렬도 벡터다	`M.ravel()` 후 랭크 재기
$\dim(S+U) = \dim S + \dim U - \dim(S \cap U)$	1절의 $6 + 6 - 3 = 9$
랭크 1 $=$ $\vv{u}\vv{v}^{\mathsf{T}}$	`np.outer(u, v)`
$A\vv{x} = (\vv{v} \cdot \vv{x})\vv{u}$	2절의 비교
랭크 $r$ $=$ 조각 $r$ 개	`랭크1로_쪼개기` 와 누적 랭크
미분도 행렬이다	3x3 행렬 `D`
$(A^k)_{ij}$ 는 경로의 개수	5절의 `A @ A`

더 해 볼 것¶

랭크1로_쪼개기 가 고르는 열을 오른쪽부터 고르도록 바꿔 보자. 조각이 달라지는가? 누적 랭크는 여전히 1, 2, 3 으로 오르는가?
차수 3 이하인 다항식으로 미분 행렬 D 를 $4 \times 4$ 로 만들어 보자. $D^2$ , $D^3$ , $D^4$ 는 무엇이 되는가? 랭크는 어떻게 변하는가?
5절에서 p 를 0에서 1까지 바꿔 가며 지름이 어떻게 변하는지 그려 보자. 어느 지점에서 가장 급격히 줄어드는가?
대칭행렬 부분공간과 반대칭행렬( $M^{\mathsf{T}} = -M$ ) 부분공간의 차원을 각각 세어 보자. 둘의 합은 무엇이 되는가? 교집합은?

다음 강의에서는 그래프를 행렬로 다루는 이야기를 이어 간다. 전기 회로에서 네 부분공간이 각각 물리적인 뜻을 갖는 것을 보게 된다.