戰(zhàn))
1. 從“解方程”到“找空間”Ax0問題的本質(zhì)在工程計(jì)算、數(shù)據(jù)分析乃至機(jī)器學(xué)習(xí)模型訓(xùn)練中我們經(jīng)常會(huì)遇到形如Ax0的線性方程組。乍一看這似乎比常見的Axb要簡(jiǎn)單畢竟右邊全是零。但恰恰是這種“簡(jiǎn)單”的形式蘊(yùn)含著線性代數(shù)最核心、也最讓人著迷的概念之一——零空間。它不像求解具體數(shù)值解那樣直接而是要找出所有可能的解所構(gòu)成的一個(gè)“集合”或“空間”。很多朋友在初次接觸時(shí)會(huì)覺得這部分內(nèi)容抽象解題時(shí)無從下手。今天我就結(jié)合自己多年在算法開發(fā)和數(shù)值計(jì)算中處理這類問題的經(jīng)驗(yàn)把Ax0的求解方法掰開揉碎了講清楚不僅告訴你“怎么做”更重點(diǎn)解釋“為什么這么做”以及在實(shí)際應(yīng)用中如何理解和運(yùn)用這個(gè)“零空間”。簡(jiǎn)單來說求解Ax0就是在尋找矩陣A的零空間的一組基。這有什么用呢舉個(gè)例子在計(jì)算機(jī)視覺中相機(jī)標(biāo)定、三維重建等任務(wù)最終都會(huì)歸結(jié)為求解一個(gè)齊次線性方程組其非零解就對(duì)應(yīng)著我們想要的空間點(diǎn)坐標(biāo)或變換矩陣。在推薦系統(tǒng)里矩陣分解的潛在因子有時(shí)也需要滿足某種齊次約束。因此掌握Ax0的求解絕非紙上談兵而是打通后續(xù)許多高級(jí)應(yīng)用的關(guān)鍵基礎(chǔ)。2. 核心思路拆解為什么是“自由變量”和“基礎(chǔ)解系”面對(duì)一個(gè)m×n的矩陣Am個(gè)方程n個(gè)未知數(shù)求解Ax0的核心思路可以概括為化繁為簡(jiǎn)找出“自由”的未知數(shù)并用它們表示出所有“受約束”的未知數(shù)從而得到通解。2.1 秩Rank的決定性作用矩陣A的秩r是整個(gè)求解過程的“總開關(guān)”。它代表了矩陣中真正獨(dú)立的行或列的個(gè)數(shù)也即有效約束的數(shù)量。如果 r n這意味著約束數(shù)量有效方程等于未知數(shù)個(gè)數(shù)且所有約束都是獨(dú)立的。通常唯一解就是零向量x0我們稱之為平凡解。此時(shí)零空間只包含一個(gè)點(diǎn)零向量。如果 r n這是更常見、也更有趣的情況。約束數(shù)量少于未知數(shù)個(gè)數(shù)意味著存在n - r個(gè)“自由度”。這些自由度對(duì)應(yīng)的未知數(shù)可以自由取值我們稱之為自由變量。而剩下的r個(gè)未知數(shù)則被這些自由變量和方程所決定。此時(shí)零空間是一個(gè)維數(shù)為n - r的向量空間存在無窮多個(gè)非零解非平凡解。我們所有求解方法的最終目標(biāo)就是系統(tǒng)地找出這n - r個(gè)自由變量并構(gòu)造出n - r個(gè)線性無關(guān)的解向量即基礎(chǔ)解系使得零空間中的任何一個(gè)解都可以表示為這些基礎(chǔ)解系的線性組合。2.2 方法選型高斯消元法行最簡(jiǎn)形是基石為什么教材和實(shí)踐中都首選高斯消元法將矩陣化為行最簡(jiǎn)形因?yàn)樗宰钪庇^、機(jī)械化的方式同時(shí)完成了兩件關(guān)鍵事確定矩陣的秩 r行最簡(jiǎn)形中非零行的數(shù)量就是秩 r。顯式地標(biāo)識(shí)出主元列和自由列每個(gè)非零行的首個(gè)非零元主元所在的列是主元列其余列則是自由列。自由列對(duì)應(yīng)的變量自然就被選為自由變量。這是一種穩(wěn)定、普適的方法無論是手算還是編程實(shí)現(xiàn)都是最可靠的起點(diǎn)。其他更高級(jí)的方法如SVD分解通常用于數(shù)值穩(wěn)定性要求極高或矩陣性質(zhì)特殊如接近奇異的場(chǎng)合但理解行最簡(jiǎn)形法是理解所有方法的基礎(chǔ)。3. 手算實(shí)戰(zhàn)一步一步求解基礎(chǔ)解系理論說再多不如動(dòng)手算一遍。我們用一個(gè)具體例子貫穿整個(gè)手算過程。設(shè)矩陣A為A [ 1 2 2 1 ] [ 2 4 1 2 ] [ 3 6 0 3 ]求解Ax0。3.1 第一步化為行最簡(jiǎn)形RREF我們對(duì)增廣矩陣[A | 0]進(jìn)行行初等變換因?yàn)橛疫吺?所以只對(duì)A操作即可。R2 R2 - 2*R1,R3 R3 - 3*R1:[ 1 2 2 1 ] [ 0 0 -3 0 ] [ 0 0 -6 0 ]R3 R3 - 2*R2:[ 1 2 2 1 ] [ 0 0 -3 0 ] [ 0 0 0 0 ]R2 R2 / (-3)(將主元化為1)然后R1 R1 - 2*R2消去主元上方的元素[ 1 2 0 1 ] [ 0 0 1 0 ] [ 0 0 0 0 ]至此我們得到了行最簡(jiǎn)形。可以看到非零行有2行所以矩陣的秩r 2??偽粗獢?shù)n 4。自由度的數(shù)量為n - r 2。這意味著零空間是二維的基礎(chǔ)解系應(yīng)包含2個(gè)線性無關(guān)的解向量。3.2 第二步識(shí)別主元列與自由變量在行最簡(jiǎn)形[ 1 2 0 1; 0 0 1 0; 0 0 0 0 ]中主元列第1列主元為1和第3列主元為1。對(duì)應(yīng)的變量x1和x3是基本變量。自由列第2列和第4列。對(duì)應(yīng)的變量x2和x4被選為自由變量。注意自由變量的選擇不是唯一的你可以選擇自由列對(duì)應(yīng)的變量也可以有其他選法但選擇自由列對(duì)應(yīng)的變量是最直接、最不容易出錯(cuò)的方法。一旦選定后續(xù)步驟就要保持一致。3.3 第三步將基本變量用自由變量表示并賦值求解向量根據(jù)行最簡(jiǎn)形我們可以直接“讀”出方程x1 2*x2 x4 0x3 0將基本變量x1,x3用自由變量x2,x4表示x3 0x1 -2*x2 - x4現(xiàn)在我們通過給自由變量賦值來構(gòu)造基礎(chǔ)解系。為了得到線性無關(guān)的解向量我們每次只讓一個(gè)自由變量為1其余為0。令x2 1,x4 0則x1 -2*1 - 0 -2x3 0得到解向量v1 [-2, 1, 0, 0]^T(T表示轉(zhuǎn)置即列向量)。令x2 0,x4 1則x1 -2*0 - 1 -1x3 0得到解向量v2 [-1, 0, 0, 1]^T。3.4 第四步寫出通解形式矩陣A的零空間N(A)就是所有解向量的集合它可以由基礎(chǔ)解系{v1, v2}線性張成。因此方程Ax0的通解為x c1 * v1 c2 * v2 c1 * [-2, 1, 0, 0]^T c2 * [-1, 0, 0, 1]^T其中c1,c2是任意實(shí)數(shù)。實(shí)操心得檢查養(yǎng)成好習(xí)慣將得到的基礎(chǔ)解系向量代回原方程Ax0驗(yàn)證。例如計(jì)算A * v1看看結(jié)果是否為零向量。這是防止計(jì)算錯(cuò)誤的最有效手段。標(biāo)準(zhǔn)化雖然基礎(chǔ)解系不唯一給自由變量賦不同的值會(huì)得到不同的基但通過“每次一個(gè)自由變量為1”的方法得到的是最簡(jiǎn)潔、標(biāo)準(zhǔn)的一組基非常便于理解和后續(xù)計(jì)算。4. 編程實(shí)現(xiàn)用NumPy進(jìn)行數(shù)值求解在實(shí)際的科研或工程項(xiàng)目中我們幾乎不會(huì)手算而是借助數(shù)值計(jì)算庫(kù)。Python的NumPy和SciPy庫(kù)是首選。這里重點(diǎn)講NumPy的方法。4.1 使用np.linalg.svd進(jìn)行奇異值分解推薦奇異值分解是數(shù)值計(jì)算中求解零空間最穩(wěn)定、最通用的方法。對(duì)于矩陣A其SVD分解為A U * S * V^T。其中V^T是右奇異向量矩陣的轉(zhuǎn)置。零空間的一組標(biāo)準(zhǔn)正交基就藏在 V 矩陣的最后 n-r 列中。import numpy as np # 定義矩陣A A np.array([[1, 2, 2, 1], [2, 4, 1, 2], [3, 6, 0, 3]], dtypefloat) # 進(jìn)行奇異值分解 U, S, Vh np.linalg.svd(A) # Vh 即 V^T # 計(jì)算矩陣的秩通過奇異值閾值 tol 1e-10 # 一個(gè)很小的閾值用于判斷奇異值是否為0 r np.sum(S tol) print(f矩陣的秩 r {r}) n A.shape[1] # 列數(shù)即未知數(shù)個(gè)數(shù) # 零空間基向量是 Vh 的最后 (n - r) 行因?yàn)閂h是V的轉(zhuǎn)置 null_space_basis Vh[r:].T # 轉(zhuǎn)置回來使得每一列是一個(gè)基向量 print(零空間的一組標(biāo)準(zhǔn)正交基列向量形式:) print(null_space_basis)運(yùn)行這段代碼你會(huì)得到兩個(gè)列向量它們張成了零空間。你會(huì)發(fā)現(xiàn)它們可能與我們手算的[-2, 1, 0, 0]^T和[-1, 0, 0, 1]^T看起來不同但它們是同一空間的兩組不同的基且是正交歸一的。你可以驗(yàn)證np.dot(A, null_space_basis[:, i])是否接近零向量。為什么推薦SVD數(shù)值穩(wěn)定性即使矩陣A是病態(tài)的或秩接近虧損SVD也能穩(wěn)健地確定其秩和零空間。直接得到標(biāo)準(zhǔn)正交基得到的基向量是兩兩正交且長(zhǎng)度為1的這在很多后續(xù)計(jì)算中非常方便。通用性適用于任意形狀的矩陣包括行數(shù)不等于列數(shù)。4.2 使用scipy.linalg.null_spaceSciPy庫(kù)提供了一個(gè)更直接的封裝函數(shù)from scipy.linalg import null_space Z null_space(A) print(Z)這個(gè)函數(shù)內(nèi)部通常也是基于SVD實(shí)現(xiàn)的是最高效快捷的方式。4.3 利用sympy進(jìn)行符號(hào)計(jì)算如果你需要得到像手算那樣精確的、分?jǐn)?shù)形式的基礎(chǔ)解系可以使用SymPy庫(kù)進(jìn)行符號(hào)運(yùn)算。import sympy as sp A sp.Matrix([[1, 2, 2, 1], [2, 4, 1, 2], [3, 6, 0, 3]]) # 計(jì)算零空間返回一個(gè)列表其中每個(gè)元素是基礎(chǔ)解系的一個(gè)向量 nullspace A.nullspace() for i, vec in enumerate(nullspace): print(f基礎(chǔ)解系向量 v{i1}:) sp.pprint(vec) print()SymPy會(huì)輸出[-2, 1, 0, 0]和[-1, 0, 0, 1]與我們手算結(jié)果完全一致。編程注意事項(xiàng)浮點(diǎn)數(shù)誤差使用NumPy/SciPy進(jìn)行數(shù)值計(jì)算時(shí)由于浮點(diǎn)數(shù)精度所謂的“零向量”可能是一個(gè)范數(shù)極小的向量如1e-15量級(jí)。判斷時(shí)應(yīng)用范數(shù)np.linalg.norm(A v)并與一個(gè)容差如1e-10比較而不是直接判斷是否等于0。秩的判斷數(shù)值計(jì)算中矩陣的“秩”是一個(gè)模糊概念。像上面代碼中通過奇異值閾值tol來判斷是標(biāo)準(zhǔn)做法。閾值的選擇需要根據(jù)具體問題的尺度來調(diào)整。5. 深入理解零空間的幾何意義與重要性質(zhì)理解Ax0的解不能只停留在代數(shù)計(jì)算層面從幾何視角看會(huì)清晰得多。5.1 幾何解釋矩陣變換下的“壓縮”與“消失”將矩陣A視為一個(gè)線性變換。方程Ax0就是在問有哪些向量 x在經(jīng)過 A 變換后被壓縮到了原點(diǎn)這些向量x的集合就是零空間N(A)。零空間的維數(shù)n-r直觀反映了這個(gè)變換“丟失”了多少信息或者說有多少個(gè)獨(dú)立的方向被“壓扁”成了零維的點(diǎn)。例如一個(gè)將三維空間投影到二維平面的變換其零空間就是一條垂直于該平面的直線一維因?yàn)檫@條直線上的所有點(diǎn)都被投影到了原點(diǎn)。5.2 與列空間、行空間的關(guān)系秩-零度定理這是線性代數(shù)中最優(yōu)美的定理之一對(duì)于 m×n 矩陣 A有 n rank(A) nullity(A)。其中rank(A)是秩列空間的維數(shù)nullity(A)是零化度零空間的維數(shù)。n定義域的維度x所在空間的維度。rank(A)值域列空間的維度即變換后像空間的“有效”維度。nullity(A)被“壓縮掉”的維度。 這個(gè)定理定量地描述了定義域在變換下如何被分割為“有效部分”列空間的原像和“無效部分”零空間。5.3 在最小二乘問題中的應(yīng)用在求解超定方程組Ax ≈ b的最小二乘解時(shí)我們求解的是A^T A x A^T b。如果A的列線性相關(guān)即秩虧那么A^T A是奇異矩陣其零空間非零。這意味著最小二乘解不唯一會(huì)有無窮多解。其中范數(shù)最小的解最小范數(shù)解可以通過將通解投影到A^T A的行空間或零空間的正交補(bǔ)上得到。這時(shí)對(duì)零空間的理解就至關(guān)重要。6. 常見陷阱、疑難解答與擴(kuò)展6.1 為什么自由變量不能選主元列對(duì)應(yīng)的變量這是一個(gè)常見的概念混淆點(diǎn)。主元列對(duì)應(yīng)的變量基本變量已經(jīng)被方程嚴(yán)格約束了。如果我們強(qiáng)行指定一個(gè)基本變量比如x1為自由變量并賦值那么由于方程的存在其他基本變量和自由變量的值可能會(huì)產(chǎn)生矛盾導(dǎo)致無法構(gòu)造出一個(gè)有效的解向量。自由變量之所以“自由”正是因?yàn)樗鼈冊(cè)谛凶詈?jiǎn)形對(duì)應(yīng)的方程中沒有對(duì)應(yīng)的主元對(duì)其進(jìn)行直接約束。6.2 矩陣行數(shù)少于列數(shù)m n就一定有無窮多解嗎不一定但可能性極大。因?yàn)橹萺 ≤ min(m, n) m。如果r n這幾乎總是成立除非矩陣非常特殊且滿行秩那么n - r 0零空間維數(shù)大于零存在無窮多非零解。只有極其特殊的情況下r n這要求m ≥ n且列滿秩但當(dāng)m n時(shí)不可能列滿秩所以對(duì)于m n的矩陣只要 A 不是零矩陣其零空間一定至少是一維的。6.3 如何判斷求出的基礎(chǔ)解系是否正確線性無關(guān)性檢查你得到的幾個(gè)解向量是否線性無關(guān)。對(duì)于二維零空間兩個(gè)向量不應(yīng)成比例。代入驗(yàn)證這是黃金準(zhǔn)則。將每個(gè)基礎(chǔ)解系向量代入原方程Ax計(jì)算結(jié)果應(yīng)為零向量允許有微小的數(shù)值誤差。維數(shù)核對(duì)基礎(chǔ)解系中向量的個(gè)數(shù)應(yīng)等于n - r。6.4 與“代數(shù)余子式”和“特征值”的聯(lián)系代數(shù)余子式在求解行列式或某些特定結(jié)構(gòu)的齊次方程組時(shí)代數(shù)余子式可能會(huì)出現(xiàn)在克萊姆法則的推導(dǎo)中但對(duì)于一般的Ax0求解行最簡(jiǎn)形法是更系統(tǒng)的方法。特征值與特征向量方程(A - λI)x 0是特征值的定義式。當(dāng) λ0 時(shí)它就退化為我們討論的Ax0。因此零空間中的非零向量就是矩陣 A 對(duì)應(yīng)于特征值 λ0 的特征向量。從這個(gè)角度看求解Ax0就是在求矩陣的“零特征值”對(duì)應(yīng)的特征空間。6.5 處理數(shù)值計(jì)算中的秩虧問題在實(shí)際數(shù)據(jù)中矩陣可能不是嚴(yán)格的秩虧而是接近秩虧某些奇異值非常小。這時(shí)嚴(yán)格數(shù)學(xué)意義上的零空間可能只包含零向量但存在一個(gè)“近似零空間”其中的向量x使得||Ax||非常小。 處理方法是設(shè)置一個(gè)閾值。在SVD中將所有小于閾值的奇異值視為0其對(duì)應(yīng)的右奇異向量就張成了這個(gè)“數(shù)值零空間”。閾值的選擇需要根據(jù)具體應(yīng)用和數(shù)據(jù)的噪聲水平來決定通??梢匀∽畲笃娈愔档哪硞€(gè)比例如1e-6倍。我個(gè)人在處理大規(guī)模數(shù)據(jù)或病態(tài)矩陣時(shí)會(huì)優(yōu)先選擇SVD方法。它的穩(wěn)定性遠(yuǎn)超基于高斯消元的QR分解求零空間的方法。手算和理解概念時(shí)行最簡(jiǎn)形法無可替代但一旦進(jìn)入代碼實(shí)戰(zhàn)scipy.linalg.null_space()是我最常用的工具它簡(jiǎn)潔且足夠穩(wěn)健。理解Ax0的求解最終是為了讓你在遇到更復(fù)雜的模型約束、優(yōu)化問題或系統(tǒng)分析時(shí)能一眼看穿其中隱藏的“自由度”和“冗余度”這是從計(jì)算員邁向設(shè)計(jì)者的關(guān)鍵一步。