SVD và xấp xỉ hạng thấp

~4-5 giờ · cần xong bài 8 trước · BÀI KẾT CỦA MÔN 2

Đích của bài này: hiểu vì sao mọi ma trận đều tách được thành xoay → co giãn theo trục → xoay, và vì sao vứt đi những trục co giãn yếu nhất là cách nén tốt nhất về mặt toán học — không phải một heuristic hay, mà là tối ưu có chứng minh.

🧭 Làm theo thứ tự này

1. Không có video nào cho bài này (3Blue1Brown không có tập về SVD).
2. Làm ô ⚡ bên dưới.
3. Đọc §1-2 và tính tay SVD theo ba bước.
4. Đọc §3-7 (§4: giá trị suy biến khác eigenvalue; §7: Eckart-Young).
5. Đọc MML ch.4 §4.5 và §4.6 (nguồn chính của bài này), đọc kỹ.
6. Code §8: tự cài SVD, không dùng np.linalg.svd. Tự viết xong mới mở code mẫu trong bài.
7. Làm 3 bài tay ở §12 (SVD tay của [[3,0],[4,5]], nén một ảnh thật, dựng ma trận phổ phẳng và phổ dốc).
8. Tự kiểm §11. Trả lời được cả 8 câu mới coi là xong môn 2.
9. Sau cùng: xem lại tập 13 và 14. Tập 16 là phần thưởng cuối môn.

Xong bài khi: nói được "SVD nén theo nghĩa nào" và giữ k giá trị lớn nhất thì mất gì, giữ gì. Câu 4 của phần tự kiểm là câu sâu nhất của cả môn 2.

⚡ 30 phút đầu — làm trước, đọc sau

import numpy as np
A = np.array([[2., 1.], [1., 2.]])
U, S, Vt = np.linalg.svd(A)
print(S)                                # đích: [3. 1.]
A1 = S[0] * np.outer(U[:, 0], Vt[0])    # giữ ĐÚNG 1 giá trị suy biến
print(A1)                               # đích: [[1.5 1.5] [1.5 1.5]]
print(np.linalg.norm(A - A1))           # đích: 1.0   ← chính là sigma_2

Ba đích: [3. 1.], [[1.5 1.5],[1.5 1.5]], 1.0.

Con số cuối là chỗ đáng dừng: bạn vừa vứt một nửa ma trận đi, và sai số bạn phải trả bằng đúng giá trị suy biến mà bạn vứt. Không xấp xỉ, không "khoảng chừng" — bằng đúng. Đó là toàn bộ bài này, cô đặc trong một dòng.

1. Ý trung tâm: mọi ma trận đều là "xoay → co giãn → xoay"

Bài 2 dạy: ma trận là một hành động lên không gian. Nhưng hành động đó có thể là bất cứ thứ gì — quay, trượt, bóp dẹp, lật. Nghe như một mớ hỗn độn.

SVD nói: không hề hỗn độn.

Mọi ma trận A kích thước m×n đều viết được thành

A = U · Σ · Vᵀ

với V và U là ma trận trực giao (chỉ xoay/lật, không làm méo), và Σ là ma trận đường chéo với các số không âm σ₁ ≥ σ₂ ≥ ... ≥ 0.

Đọc từ phải sang trái (bài 2: A·x nghĩa là làm Vᵀ trước):

Đó là tất cả những gì một phép biến đổi tuyến tính có thể làm. Ba bước. Không có bước thứ tư.

Kiểm bằng hình: cho A ăn toàn bộ đường tròn đơn vị (mọi vector độ dài 1). Ảnh ra luôn là một elip. Nửa trục dài của elip = σ₁, nửa trục ngắn = σ₂. Với A = [[2,1],[1,2]], quét 200.001 điểm trên đường tròn:

max ||A·x|| = 3.0000000000000004     ← σ₁
min ||A·x|| = 1.0                    ← σ₂

σ₁ không phải một số trừu tượng. Nó là hệ số phóng đại lớn nhất mà ma trận có thể gây ra cho một vector. σ_min là hệ số nhỏ nhất. Trong ngôn ngữ backend: σ₁ là worst-case latency của phép biến đổi, σ_min là best-case. Tỉ số σ₁/σ_min là số điều kiện — đúng con số mà thư viện số học dùng để cảnh báo "ma trận này giải ra sẽ mất chữ số".

Vì sao đây là chuyện lớn — so với bài 7

Bài 7 cho bạn A = P·D·P⁻¹ (chéo hoá bằng eigenvector). Nhưng nó có điều kiện đầu vào:

Chéo hoá là một API có thể ném exception. SVD là API không bao giờ ném. Đó là lý do trong thực hành số học người ta dùng SVD chứ không dùng eigen, kể cả khi cả hai đều chạy.

2. Tính tay — công thức ba bước

Cách tìm SVD bằng tay dùng lại nguyên xi máy móc của bài 7. Chỉ cần một mẹo:

A = U Σ Vᵀ   ⟹   AᵀA = (U Σ Vᵀ)ᵀ (U Σ Vᵀ) = V Σᵀ Uᵀ U Σ Vᵀ = V Σ² Vᵀ
                                              └── Uᵀ U = I vì U trực giao

AᵀA là ma trận đối xứng, và V Σ² Vᵀ chính là chéo hoá của nó. Suy ra:

Công thức 3 bước:

1. Tính AᵀA. Tìm eigenvalue λᵢ của nó (máy móc bài 7). Đặt σᵢ = √λᵢ.
2. Eigenvector của AᵀA (chuẩn hoá) = các cột của V.
3. uᵢ = A·vᵢ / σᵢ — cho từng σᵢ > 0.

Ví dụ chính: A = [[2,1],[1,2]]

Bước 1.

AᵀA = [ 2  1 ] [ 2  1 ]  =  [ 5  4 ]
      [ 1  2 ] [ 1  2 ]     [ 4  5 ]

Ma trận [[5,4],[4,5]] có eigenvector (1,1) và (1,−1) — nhìn ra ngay vì nó đối xứng và hai hàng là hoán vị của nhau:

[ 5  4 ][ 1 ]   [ 9 ]                    [ 5  4 ][  1 ]   [  1 ]
[ 4  5 ][ 1 ] = [ 9 ] = 9·(1, 1)         [ 4  5 ][ −1 ] = [ −1 ] = 1·(1, −1)

→ λ₁ = 9, λ₂ = 1 → σ₁ = 3, σ₂ = 1.

Bước 2. v₁ = (1,1)/√2, v₂ = (1,−1)/√2.

Bước 3.

u₁ = A·v₁/σ₁ = (3,3)/√2 / 3 = (1,1)/√2
u₂ = A·v₂/σ₂ = (1,−1)/√2 / 1 = (1,−1)/√2

Ở đây U = V vì A đối xứng — trường hợp đặc biệt, đừng tưởng luôn vậy.

Ba phép kiểm tra chéo — làm luôn, đừng bỏ

Đẳng thức thứ hai là linh hồn của phép nén ở §5. σᵢ² là năng lượng trục i gánh. Tổng năng lượng cố định. Giữ trục lớn = giữ được nhiều năng lượng.

Bốn ma trận mốc — thuộc lòng chúng

🪤 BẪY: dấu của U và V

Chạy np.linalg.svd([[2,1],[1,2]]) bạn sẽ nhận:

U  = [[−0.7071, −0.7071], [−0.7071,  0.7071]]
Vᵀ = [[−0.7071, −0.7071], [−0.7071,  0.7071]]

Cột đầu của U là (−0.7071, −0.7071) chứ không phải (0.7071, 0.7071) như bạn tính tay. Không phải bạn sai. Nếu đổi dấu cả uᵢ và vᵢ thì σᵢ·uᵢvᵢᵀ không đổi — hai dấu âm triệt tiêu.

Không bao giờ so sánh U hay V từng phần tử với numpy. Chỉ so σ và so kết quả tái tạo U Σ Vᵀ. Đây là nguồn lỗi số 1 khi tự cài SVD, và nó không phải bug.

Bài tập tính tay bắt buộc

C = [[3,0],[4,5]]. Tự làm ba bước rồi đối chiếu:

CᵀC = [[25, 20], [20, 25]]      λ = 45, 5      σ = 3√5 ≈ 6.7082,  √5 ≈ 2.2361
det C = 15                       kiểm: 3√5 · √5 = 3·5 = 15 ✓

3. Ma trận hạng 1 — chỗ SVD gặp lại bài 3

B = [[1,2],[2,4]], ma trận bạn đã gặp ở bài 3 (hàng 2 = 2× hàng 1).

BᵀB = [ 5  10]        λ = 25, 0        σ = 5, 0
      [10  20]

Một giá trị suy biến khác 0. Và:

v₁ = u₁ = (1,2)/√5
σ₁·u₁v₁ᵀ = 5 · (1/5)·[[1,2],[2,4]] = [[1., 2.], [2., 4.]] = B  ← khớp CHÍNH XÁC

Nói cách khác B = (1,2)ᵀ·(1,2) — một tích ngoài (outer product) của một vector với chính nó. Bốn số của B thật ra chỉ chứa hai số thông tin.

Kiểm thêm: ‖B‖_F = 5 = σ₁. Chỉ có một trục, nên nó gánh toàn bộ năng lượng.

hạng = số giá trị suy biến khác 0. Bài 3 định nghĩa hạng bằng "số cột độc lập tuyến tính" — đúng nhưng khó đo trên máy: cột "gần phụ thuộc" thì tính sao? SVD trả lời: đo bằng σ. σ nhỏ nhưng khác 0 = "gần phụ thuộc". Hạng thôi là câu hỏi có/không; phổ σ là câu hỏi bao nhiêu.

Và đây là chỗ trả nốt BẪY 3 của bài 3 — "hạng theo cột luôn bằng hạng theo hàng", cái định lý bài 3 bảo bạn "chấp nhận và dùng, sẽ thấy hiển nhiên sau khi học SVD". Giờ nó hiển nhiên thật:

A = U Σ Vᵀ        ⟹        Aᵀ = V Σ Uᵀ

A và Aᵀ dùng chung một Σ — chỉ đổi chỗ U với V. Cùng bộ σ ⟹ cùng số σ khác 0 ⟹ cùng hạng. Không cần chứng minh gì thêm: hạng hàng và hạng cột không phải hai đại lượng tình cờ bằng nhau, chúng là một đại lượng đọc từ hai phía.

Hạng thấp không có nghĩa là ma trận nhỏ

Ma trận 50×50 với G[i,j] = i + j. 2500 số. Hạng bằng bao nhiêu?

hạng(G) = 2
σ = [2646.7067, 196.7067, 0, 0, 0, ...]      σ₃ lớn nhất còn lại = 2.5e-13 (nhiễu làm tròn)

Đúng 2. Vì G = i·1ᵀ + 1·jᵀ — tổng của hai tích ngoài. 2500 số, nhưng chỉ 2 lớp thông tin.

Đây là toàn bộ trực giác về nén: kích thước ma trận nói lên khả năng chứa, σ nói lên nội dung thật sự đang chứa. Giống bảng SQL 200 cột nhưng 198 cột suy ra được từ 2 cột kia.

4. Giá trị suy biến ≠ eigenvalue — chỗ eigenvalue nói dối

Đây là mục ngắn nhưng đắt nhất của bài. Lấy phép trượt ngang từ bài 2:

S = [ 1  1 ]
    [ 0  1 ]

Eigenvalue của S: 1 và 1. Đọc theo bài 7: "không co giãn gì cả, hệ số 1".

Nhưng nhìn hình đi — phép trượt rõ ràng có kéo giãn không gian. Eigenvalue nói dối.

Singular value của S:

SᵀS = [[1, 1], [1, 2]]     λ = 0.382, 2.618     σ = 1.618, 0.618

1.618... chính là tỉ số vàng (1+√5)/2, và 0.618... = (√5−1)/2 là nghịch đảo của nó. Kiểm: σ₁·σ₂ = 1.0 = |det S| ✓.

→ Phép trượt kéo dài 1.618 lần theo một hướng và nén còn 0.618 lần theo hướng vuông góc. Hướng vào nghiêng 58.28°, hướng ra nghiêng 31.72° — nó xoay hẳn.

Tệ hơn nữa: ma trận eigenvector của S có hạng 1 — hai eigenvector trùng nhau, S không chéo hoá được. Bài 7 bó tay hoàn toàn. SVD vẫn cho ra ba ma trận sạch sẽ.

Và ma trận quay R = [[0,−1],[1,0]]:

eigenvalue = 0+1j, 0−1j     ← SỐ PHỨC
σ          = 1, 1           ← thực, dương, đọc được ngay: "chỉ xoay, không co giãn"

Câu chốt: eigenvalue đo "cái gì đứng yên". Singular value đo "cái gì bị bóp méo nhiều nhất". Khi mục tiêu là nén, bạn cần cái thứ hai. Đó là lý do bài 9 không phải là bài 7 nói lại. (Ngoại lệ đáng nhớ: với ma trận đối xứng nửa xác định dương, hai thứ trùng nhau — và đó chính là chỗ PCA sống, xem §9.)

5. Cắt = nén. Truncated SVD

Viết lại A = U Σ Vᵀ thành tổng:

A = σ₁·u₁v₁ᵀ + σ₂·u₂v₂ᵀ + ... + σ_r·u_r v_rᵀ
    └───────┘   └───────┘         └─────────┘
     hạng 1      hạng 1             hạng 1

Mọi ma trận là một chồng các lớp hạng 1, đã được xếp sẵn theo thứ tự quan trọng giảm dần. SVD không "tính" ra điều gì mới — nó sắp xếp thứ đã có sẵn.

Giữ k lớp đầu, vứt phần còn lại:

A_k = σ₁·u₁v₁ᵀ + ... + σ_k·u_k v_kᵀ           ← hạng đúng bằng k

Chi phí lưu trữ. Ma trận m×n gốc cần m·n số. A_k cần:

k vector u (mỗi cái m số) + k vector v (mỗi cái n số) + k số σ  =  k·(m + n + 1)

Với m = n = 50: k·101 so với 2500.

🪤 BẪY LỚN NHẤT CỦA CẢ BÀI. Nén chỉ có lãi khi k·(m+n+1) < m·n, tức k < mn/(m+n+1). Với 50×50 điểm hoà vốn là k = 24.75 — giữ 25 thành phần là đã lỗ (101% dung lượng). Với 1000×1000 điểm hoà vốn là k ≈ 499.

Đây là bài toán Big-O quen thuộc: chi phí gốc O(mn), chi phí nén O(k(m+n)). Nén chỉ thắng khi k ≪ min(m,n). Ma trận càng to, thắng càng đậm. Ma trận vuông nhỏ, gần như không đáng.

6. Nén thật — bảng số đo

Ma trận thí nghiệm, 50×50, có cấu trúc và tự sinh được:

i, j = np.meshgrid(np.arange(50), np.arange(50), indexing="ij")
M = np.exp(-np.abs(i - j) / 8.0)

Đọc là: "50 mục xếp theo thứ tự; hai mục càng gần nhau càng giống nhau, độ giống giảm theo hàm mũ." Đây đúng là hình dạng của một ma trận điểm attention có thiên lệch cục bộ, hoặc một ma trận tương quan chuỗi thời gian.

Hạng đầy đủ = 50 (không có cột nào thừa!). ‖M‖_F = 19.241815.

Mười giá trị suy biến đầu:

13.9562  9.8643  6.4255  4.2099  2.8686  2.0439  1.5167  1.1650  0.9211  0.7460

Đọc lại dòng k = 5: giữ 5 trên 50 thành phần — 10% — và 97.048% năng lượng vẫn còn. Đó là câu trả lời cho "SVD nén theo nghĩa nào?".

Nhìn theo chiều ngược lại — cần bao nhiêu σ để đạt một mức năng lượng:

Đối chứng — nhiễu ngẫu nhiên KHÔNG nén được

Cùng kích thước 50×50, nhưng mỗi phần tử là một số ngẫu nhiên độc lập:

N = np.random.default_rng(7).normal(size=(50, 50))     # ghim seed để bạn ra đúng số dưới đây

Tám giá trị suy biến đầu:

13.410  13.127  12.465  11.828  11.418  11.109  10.615  10.439  ...   σ nhỏ nhất = 0.2062

Không hề dốc xuống — phẳng lì.

Để đạt 90% năng lượng cần k = 26 → 105.04% dung lượng. Nén xong to hơn bản gốc.

Đây là điều quan trọng nhất trong toàn bài 9. SVD không nén được gì cả. DỮ LIỆU mới là thứ nén được — SVD chỉ tìm ra chỗ nén và đo nó.

Ẩn dụ backend: gzip không làm file nhỏ đi. Sự lặp lại trong file làm file nhỏ đi; gzip chỉ phát hiện ra. Đưa gzip một file /dev/urandom thì nó phình ra. Phổ σ chính là bản báo cáo "dữ liệu của bạn nén được tới đâu" — và bạn đọc được nó trước khi nén.

7. Eckart-Young — vì sao "giữ k cái lớn nhất" là TỐT NHẤT

Cắt k thành phần đầu là cách hiển nhiên. Nhưng hiển nhiên ≠ tối ưu. Định lý Eckart-Young nói nó đúng là tối ưu:

Định lý Eckart-Young (–Mirsky). Trong tất cả các ma trận B có hạng(B) ≤ k, ma trận A_k (giữ k giá trị suy biến lớn nhất) làm ‖A − B‖ nhỏ nhất — đúng với cả chuẩn Frobenius lẫn chuẩn phổ.

Và sai số biết trước, không cần tìm kiếm:

‖A − A_k‖_F  = √(σ²_{k+1} + σ²_{k+2} + ... )
‖A − A_k‖₂   = σ_{k+1}

Kiểm trên ma trận M ở §6 với k = 5:

‖M − M₅‖_F = 3.30610858        √(Σ_{i>5} σᵢ²) = 3.30610858   ✓
‖M − M₅‖₂  = 2.04393807        σ₆             = 2.04393807   ✓

Ba đối thủ, ba lần thua

Đối thủ 1 — chọn nhầm bộ thành phần. Vẫn giữ đúng 5 lớp, chỉ đổi chọn lớp nào:

Đối thủ 2 — 20.000 không gian con hạng 5 ngẫu nhiên. Với mỗi không gian, lấy xấp xỉ tốt nhất bên trong nó (phép chiếu trực giao):

tốt nhất trong 20.000 lần  = 15.836391
tệ nhất                    = 19.040407
SVD                        =  3.306109      ← tốt hơn kẻ may mắn nhất 4.79 lần

Đối thủ 3 — tối ưu hoá thật. ALS (bình phương tối thiểu xen kẽ), 400 vòng lặp, 30 điểm khởi đầu ngẫu nhiên khác nhau — một thuật toán thực sự đi tìm lời giải:

min = 3.306109      max = 3.306109      có lần nào thắng SVD không? False

Cả 30 lần đều hội tụ về đúng con số của SVD, không lần nào vượt qua. Không thể vượt được.

Vì sao đây là chuyện lạ. Trong hầu hết bài toán tối ưu, chiến lược tham lam (greedy) cho lời giải tạm được, còn tối ưu thật phải tìm kiếm. Ở đây tham lam CHÍNH LÀ tối ưu, có chứng minh.

Với dân backend: giống có một query planner mà bạn chứng minh được nó luôn ra kế hoạch rẻ nhất, thay vì một đống heuristic. Loại đảm bảo đó gần như không tồn tại. Ở đây nó tồn tại — và đó là lý do low-rank có mặt ở khắp mọi nơi trong ML.

Hệ quả thực dụng: σ_{k+1} là trần lỗi đọc được trước. Muốn biết cắt còn k chiều thì hỏng bao nhiêu? Không cần thử. Nhìn vào phổ. Đây là khác biệt giữa đo và đoán.

8. Code — tự cài, KHÔNG dùng np.linalg.svd

Cách cài dùng lại đúng power iteration của bài 8 (PageRank), cộng thêm một bước gọi là bóc tách (deflation): tìm lớp mạnh nhất, trừ nó ra khỏi ma trận, lặp lại.

import numpy as np


def gia_tri_suy_bien_lon_nhat(R, so_vong=1000, seed=0):
    """Power iteration trên RᵀR — đúng kỹ thuật lặp của bài 7/8."""
    rng = np.random.default_rng(seed)
    v = rng.normal(size=R.shape[1])
    v /= np.linalg.norm(v)
    for _ in range(so_vong):
        w = R.T @ (R @ v)               # KHÔNG dựng RᵀR tường minh — đắt và mất chính xác
        chuan = np.linalg.norm(w)
        if chuan < 1e-300:
            return 0.0, np.zeros(R.shape[0]), v
        v = w / chuan
    Rv = R @ v
    sigma = np.linalg.norm(Rv)
    u = Rv / sigma if sigma > 1e-300 else np.zeros(R.shape[0])
    return sigma, u, v


def svd_tu_cai(Ain, k):
    """Lấy k bộ (sigma, u, v) đầu: power iteration + BÓC TÁCH."""
    R = np.array(Ain, dtype=float)
    sigmas, us, vs = [], [], []
    for t in range(k):
        s, u, v = gia_tri_suy_bien_lon_nhat(R, seed=t)
        sigmas.append(s); us.append(u); vs.append(v)
        R = R - s * np.outer(u, v)      # bóc lớp vừa tìm được ra khỏi ma trận
    return np.array(sigmas), np.array(us).T, np.array(vs)


# (a) đối chiếu trên ma trận tính tay ở §2
A = np.array([[2., 1.], [1., 2.]])
s_toi, U_toi, Vt_toi = svd_tu_cai(A, 2)
assert np.allclose(s_toi, [3., 1.], atol=1e-9)
assert np.allclose(U_toi @ np.diag(s_toi) @ Vt_toi, A, atol=1e-9)
print("2x2 tính tay: sigma =", np.round(s_toi, 8), " khớp [3, 1] ✓")

# (b) đối chiếu trên ma trận 50x50 thật
n = 50
i, j = np.meshgrid(np.arange(n), np.arange(n), indexing="ij")
M = np.exp(-np.abs(i - j) / 8.0)
s_toi, U_toi, Vt_toi = svd_tu_cai(M, 10)
s_np = np.linalg.svd(M)[1][:10]
assert np.allclose(s_toi, s_np, atol=1e-9)
print("50x50: lệch sigma so với numpy =", np.abs(s_toi - s_np).max())

# (c) BẢY đẳng thức phải đúng — đây là phần kiểm tra hiểu, không phải kiểm tra code
U, S, Vt = np.linalg.svd(A)
assert np.allclose(U @ np.diag(S) @ Vt, A)                     # tái tạo
assert np.allclose(U.T @ U, np.eye(2))                         # U trực giao
assert np.allclose(Vt @ Vt.T, np.eye(2))                       # V trực giao
assert np.isclose(np.linalg.norm(A), np.sqrt((S**2).sum()))    # Frobenius = sqrt(sum sigma^2)
assert np.isclose(abs(np.linalg.det(A)), np.prod(S))           # |det| = tích sigma
assert np.linalg.matrix_rank(A) == int((S > 1e-12).sum())      # hạng = số sigma khác 0
assert np.isclose(np.linalg.cond(A), S[0] / S[-1])             # số điều kiện
print("các đẳng thức khớp numpy ✓")

# (d) ECKART-YOUNG: không ai thắng được phép cắt
k = 5
U, S, Vt = np.linalg.svd(M)
Mk = (U[:, :k] * S[:k]) @ Vt[:k]
sai_so_svd = np.linalg.norm(M - Mk)
assert np.isclose(sai_so_svd, np.sqrt((S[k:]**2).sum()))       # Frobenius
assert np.isclose(np.linalg.svd(M - Mk)[1][0], S[k])           # chuẩn phổ = sigma_{k+1}


def als(Min, k, so_vong=400, seed=0):
    """Đối thủ: bình phương tối thiểu xen kẽ, khởi đầu ngẫu nhiên."""
    rng = np.random.default_rng(seed)
    W = rng.normal(size=(Min.shape[0], k))
    for _ in range(so_vong):
        H = np.linalg.lstsq(W, Min, rcond=None)[0]
        W = np.linalg.lstsq(H.T, Min.T, rcond=None)[0].T
    return W @ H


sai_so_als = [np.linalg.norm(M - als(M, k, seed=s)) for s in range(30)]
assert min(sai_so_als) >= sai_so_svd - 1e-9, "có kẻ thắng SVD?! kiểm tra lại"
print(f"SVD = {sai_so_svd:.6f}   ALS tốt nhất sau 30 lần = {min(sai_so_als):.6f}  → không ai thắng ✓")

Kết quả phải ra:

2x2 tính tay: sigma = [3. 1.]  khớp [3, 1] ✓
50x50: lệch sigma so với numpy = 4.440892098500626e-15
các đẳng thức khớp numpy ✓
SVD = 3.306109   ALS tốt nhất sau 30 lần = 3.306109  → không ai thắng ✓

Quy tắc 1 đang luyện ở đây: khi svd_tu_cai lệch numpy, thứ tự nghi ngờ là (1) dấu u/v — không phải bug, đọc lại BẪY ở §2; (2) so_vong quá ít khi hai σ liền nhau gần bằng nhau (power iteration hội tụ theo tỉ số (σ_{i+1}/σᵢ)² — hai σ sát nhau thì bò rất chậm); (3) sai chỉ số khi bóc tách. Ghi vào nhật ký lỗi.

🪤 Bẫy số học: đừng bao giờ dựng AᵀA tường minh

Trong code trên có R.T @ (R @ v) chứ không phải (R.T @ R) @ v. Lý do: số điều kiện của AᵀA là bình phương số điều kiện của A. Với A = [[2,1],[1,2]], cond(A) = 3, cond(AᵀA) = 9. Với ma trận thật cond 10⁶ thì AᵀA có cond 10¹² — ăn mất 12 trong ~16 chữ số có nghĩa của float64, gần như không còn gì. Công thức 3 bước ở §2 rất tốt để tính tay, và rất tệ để cài vào máy.

Chi phí — đo thật

Tính SVD đầy đủ là O(mn·min(m,n)), tức O(n³) với ma trận vuông. Đo trên máy:

Gấp đôi cạnh → chậm khoảng 3–5 lần và tỉ số đang tiến dần về 8 (= 2³). Với ma trận 100k×100k thì SVD đầy đủ là không tưởng — đó là lý do tồn tại randomized SVD và truncated SVD (chỉ tìm k thành phần đầu, đúng như svd_tu_cai ở trên).

9. Nối vào AI

PCA = SVD trên dữ liệu đã trừ trung bình

300 điểm 2 chiều, tâm ở (50, −20):

rng = np.random.default_rng(42)
X = rng.normal(size=(300, 2)) @ np.array([[3., 1.], [0., 0.7]]) + np.array([50., -20.])

Xc = X - X.mean(axis=0)                     # BƯỚC SỐNG CÒN
U, S, Vt = np.linalg.svd(Xc, full_matrices=False)
print(S**2 / (len(Xc) - 1))    # -> [10.20903819  0.380446  ]   phương sai theo trục chính
print(Vt[0])                   # -> [-0.93406786 -0.35709555]   trục chính thứ nhất

PCA không phải một thuật toán riêng. Nó là "trừ trung bình, rồi lấy SVD, rồi gọi V là trục chính". Ba dòng. Môn 6 sẽ chứng minh vì sao trục có σ lớn nhất chính là trục phương sai lớn nhất — và mắt xích là đúng công thức §2: ma trận hiệp phương sai chính là XᵀX/(n−1), nên eigenvalue của nó = σ²/(n−1). PCA và SVD không "giống nhau"; chúng là một thứ, viết bằng hai bộ từ vựng.

🪤 Bẫy: quên trừ trung bình. Cùng dữ liệu, bỏ dòng X - X.mean():

v₁ đúng (đã trừ)      = (−0.9341, −0.3571)
v₁ khi QUÊN trừ       = (−0.9287,  0.3709)     ← lệch hẳn hướng
hướng trung bình      = ( 0.9280, −0.3725)     ← v₁ bám theo đúng cái này

Không trừ trung bình thì σ₁ chỉ đo dữ liệu nằm xa gốc toạ độ bao nhiêu, không đo phương sai. Một lỗi một dòng, phá sạch kết quả, và không báo lỗi.

LoRA — vì sao thêm dưới 1% tham số mà đổi được hành vi model

Fine-tune đầy đủ nghĩa là học cả ΔW kích thước d×d. LoRA ép ΔW = B·A với B là d×r, A là r×d. Theo bài 3, hạng(BA) ≤ r — luôn luôn, dù B, A có là gì.

(Bảng chi phí này bạn đã gặp ở bài 5 §8 — ở đó nó chỉ là hệ quả của hạng(BA) ≤ r. Lặp lại ở đây vì bây giờ bạn đọc được vế còn thiếu: r chính là k của §5, và Eckart-Young ở §7 nói cắt ở đâu thì mất bao nhiêu.)

Kiểm bằng code: rank(B @ A) với d=4096, r=8 cho ra đúng 8. B@A có 16.777.216 phần tử, nhưng chỉ 65.536 tham số tự do.

Giả thiết ngầm của LoRA, phát biểu bằng ngôn ngữ bài 9: "phần cập nhật cần thiết ΔW vốn dĩ có phổ σ dốc — vài σ đầu đã gánh gần hết năng lượng." Nếu giả thiết đó sai — nếu ΔW có phổ phẳng như ma trận nhiễu ở §6 — LoRA sẽ không hoạt động. Nó hoạt động là bằng chứng thực nghiệm rằng phổ dốc. Và r chính là k của §5: quá nhỏ thì thiếu, quá lớn thì hết rẻ.

MLA — nén KV cache

Số học minh hoạ theo cấu hình kiểu DeepSeek-V2: 60 lớp, 128 đầu, d_head=128, latent d_c=512 cộng phần RoPE 64.

MHA thường : 2 · n_head · d_head · L = 2 · 128 · 128 · 60 = 1,966,080 số / token
MLA        : (d_c + d_rope) · L      = (512 + 64) · 60    =    34,560 số / token
                                                     tỉ lệ = 1.758%  → giảm 98.24%

Với ngữ cảnh 100k token, lưu fp16: 393.22 GB → 6.912 GB. Đó là khác biệt giữa "không chạy nổi" và "chạy được trên một máy".

Cơ chế: thay vì cache đầy đủ K và V, cache một vector tiềm ẩn c nhỏ, rồi dựng lại K, V từ c bằng hai ma trận chiếu ngược. Chính là A ≈ U_k Σ_k V_kᵀ — chỉ khác ở chỗ ma trận nén được học thay vì được tính bằng SVD.

Nén embedding

1.000.000 vector × 768 chiều × 4 byte  = 3.072 GB
cắt còn 128 chiều                      = 0.512 GB
+ ma trận cơ sở V (768×128)            = 0.393 MB      ← lưu MỘT lần cho cả kho
                                  tổng → 16.68% dung lượng gốc

Ma trận V chỉ tốn 0.393 MB nhưng cần cho mọi truy vấn — vì query cũng phải chiếu xuống cùng không gian 128 chiều trước khi so cosine (bài 1). Quên chiếu query = so hai vector ở hai hệ trục khác nhau = rác.

Vì sao attention "gần" hạng thấp

Với T = 64 token và d_head = 8:

hạng(QKᵀ/√d)      = 8      ← không thể lớn hơn d_head, theo bài 3
hạng(softmax(...)) = 64     ← softmax phá vỡ tính hạng thấp
σ của QKᵀ/√d:   28.028  26.143  24.807  22.163  21.382  18.249  16.238  13.537  0  0 ...
σ của softmax:   1.0286  0.5015  0.4434  0.4152  0.3624  0.3412  0.2768  0.2590  0.2207 ...

Ma trận sau softmax có hạng đầy đủ 64, nhưng 8 giá trị suy biến đầu vẫn giữ 85.99% năng lượng.

Đây là móng của cả một dòng nghiên cứu ("linear attention", "Performer", "Linformer"): ma trận attention không hạng thấp về mặt chặt chẽ, nhưng phổ của nó dốc đủ để xấp xỉ hạng thấp có lãi. Toàn bộ tranh luận đó diễn ra bằng đúng từ vựng của §6-§7 bài này.

Trước khi sang môn 3: σ₁ của tích Jacobian

Backprop qua L lớp là nhân L ma trận Jacobian. Chuẩn phổ có tính dưới nhân: ‖A·B‖₂ ≤ ‖A‖₂·‖B‖₂. Nhân 40 ma trận 20×20 ngẫu nhiên:

mỗi lớp co lại một chút → σ₁ của tích sau 40 lớp = 3.38e-04     ← gradient TAN BIẾN
mỗi lớp giãn ra một chút → σ₁ của tích sau 40 lớp = 2.43e+03     ← gradient BÙNG NỔ

Gradient vanishing/exploding không phải hiện tượng bí ẩn của deep learning. Nó là câu "σ nhân với nhau qua các lớp" — một câu thuần đại số tuyến tính, viết được bằng chữ của bài 9.

10. Nhìn lại môn 2 — chín bài đã dựng nên cái gì

Kiểm tra thật xem môn 2 đã xong chưa: quay lại mon-02-dai-so-tuyen-tinh/README.md §6 và trả lời cả 7 câu "tận cùng" ở đó. Câu 6 ("SVD nén theo nghĩa nào?") giờ phải trả lời được bằng ba tầng: cắt σ nhỏ (§5), đo bằng năng lượng (§6), và tối ưu theo Eckart-Young (§7).

Và làm nốt capstone: PageRank tự viết khớp nx.pagerank sai số < 1e-6 trên ít nhất 3 đồ thị.

11. Tự kiểm — trả lời được mới coi là xong môn 2

1. "Mọi ma trận đều là xoay → co giãn → xoay." Ba chữ đó ứng với U, Σ, V thế nào, và vì sao đọc từ phải sang trái?
2. Vì sao σᵢ = √(eigenvalue của AᵀA)? Dẫn lại trong 2 dòng, và chỉ rõ chỗ nào dùng tới tính trực giao của U.
3. Phép trượt [[1,1],[0,1]] có eigenvalue 1, 1 nhưng singular value 1.618, 0.618. Con số nào mô tả đúng thứ ma trận đang làm với không gian, và vì sao con số kia nói dối?
4. Ma trận nhiễu ngẫu nhiên 50×50 cần k = 26 mới đạt 90% năng lượng, tức 105% dung lượng gốc. Điều đó nói gì về câu "SVD là thuật toán nén"? Phát biểu lại cho đúng.
5. Eckart-Young nói A_k là xấp xỉ hạng k tốt nhất. Nêu hai hệ quả thực dụng: một về khỏi phải tìm kiếm gì, một về biết trước sai số.
6. LoRA chỉ thêm 0.3906% tham số mà đổi được hành vi model. Giả thiết ngầm về ΔW là gì, và nếu giả thiết đó sai thì LoRA hỏng ở đâu?
7. PCA quên trừ trung bình thì v₁ bám theo hướng trung bình thay vì hướng phương sai lớn nhất. Giải thích bằng ngôn ngữ SVD vì sao đúng như vậy — đừng dẫn công thức PCA.
8. Ma trận 50×50 với G[i,j] = i + j có hạng đúng bằng 2. Chứng minh trong một dòng bằng tích ngoài, không cần tính σ.

Câu 4 là câu sâu nhất của cả môn 2. Nó tách "biết gọi np.linalg.svd" khỏi "hiểu SVD là gì". SVD không tạo ra khả năng nén; nó phát hiện và đo khả năng nén vốn đã nằm sẵn trong dữ liệu. Ai chưa nói được câu đó thì chưa qua môn.

Câu 3 và câu 6 là hai câu đáng thứ nhì — một cho phần toán, một cho phần ứng dụng.

12. Bài tập về nhà — trước khi sang môn 3

🚨 3Blue1Brown KHÔNG có tập nào về SVD

Nói thẳng để bạn khỏi mất thời gian đi tìm: series Essence of Linear Algebra dừng ở tập 16 và không hề chạm tới SVD. Không có tập nào để xem cho bài này. Đừng tìm nữa.

Thay vào đó:

MML ch.4 §4.7 Matrix Phylogeny là một trang sơ đồ xếp mọi loại ma trận vào một cây phả hệ. Xem nó sau khi đọc xong §4.5-4.6, coi như bản đồ tổng kết môn 2.

Ba bài tay bắt buộc

1. Tính SVD bằng tay cho C = [[3,0],[4,5]] — đủ ba bước ở §2. Đối chiếu: σ = 3√5, √5, det = 15. Rồi kiểm bằng np.linalg.svd, chỉ so σ và bản tái tạo, đừng so dấu.
2. Nén một tấm ảnh thật. Lấy một ảnh xám bất kỳ dạng ma trận m×n, dựng bảng như §6: k | sai số tương đối | % năng lượng | % dung lượng, với k = 1, 5, 20, 50, 100. Trả lời hai câu: k nhỏ nhất để mắt thường không phân biệt được là bao nhiêu? Và điểm hoà vốn k = mn/(m+n+1) của ảnh đó nằm ở đâu?
3. Tự dựng một ma trận có phổ phẳng (nén không được) và một ma trận có phổ dốc (nén tốt), rồi giải thích bằng lời cấu trúc nào trong ma trận thứ hai khiến phổ dốc xuống.

Teaser môn 3 — Giải tích & Tối ưu

Ở §9 bạn vừa thấy: nhân 40 ma trận Jacobian, mỗi cái co lại một chút, σ₁ của tích tụt xuống 3.38e-04. Mỗi cái giãn ra một chút, σ₁ vọt lên 2.43e+03.

Môn 3 sẽ dạy chain rule, và sự thật lớn nhất của nó là: "chain rule qua nhiều lớp chính là nhân một chuỗi ma trận Jacobian". Nghĩa là backpropagation không phải thuật toán mới — nó là phép nhân ma trận của bài 2, lặp lại. Và gradient vanishing chỉ là bài toán σ của bài 9, mặc áo giải tích.

Bạn đã có sẵn nửa câu trả lời trước khi môn 3 bắt đầu.