建、索引到運算與避坑)
學MATLAB繞不開矩陣就像學中文繞不開漢字。前幾天整理自己的學習筆記發(fā)現(xiàn)大部分基礎內容其實都和矩陣有關。這篇筆記主要是記錄我入門MATLAB時最常碰到的矩陣操作怎么創(chuàng)建、怎么索引、怎么運算以及一堆當時踩過的坑。如果你是剛開始接觸MATLAB或者想系統(tǒng)回顧矩陣這部分的基礎知識這份筆記應該能幫上忙。1. 為什么MATLAB里一切皆矩陣1.1 矩陣不是普通的“二維表”MATLAB的全稱是Matrix Laboratory翻譯過來就是“矩陣實驗室”所以矩陣這個概念幾乎刻在它骨子里。剛開始學習時我總習慣把矩陣理解成Excel表格后來才發(fā)現(xiàn)這個類比很危險。矩陣在MATLAB里不僅僅是一個表格而是一個有嚴格維度定義的數(shù)據(jù)對象。更關鍵的是MATLAB里的標量、向量、矩陣、多維數(shù)組全都統(tǒng)一在“數(shù)組”這個框架下一個數(shù)可以看成1×1矩陣一個行向量是1×n矩陣一個列向量是n×1矩陣二維矩陣自然就是m×n再往上還可以有頁、可以有更高維度。這個設計帶來的直接好處是你不需要像C語言那樣區(qū)分數(shù)組、指針、結構體幾乎所有數(shù)據(jù)都可以用矩陣來表達。比如你讀一張圖片它其實是一個三維矩陣長寬是圖片尺寸第三維是RGB通道你錄一段聲音它往往是一個很長的列向量一個表格數(shù)據(jù)導入后也常常變成一個二維矩陣。理解這一點后你再看MATLAB的很多操作就順了因為大部分內置函數(shù)都在和矩陣打交道。1.2 矩陣的維度、行數(shù)、列數(shù)和元素總數(shù)在入門階段有三個函數(shù)很容易混淆size、length和numel。size(A)返回的是A每個維度的大小比如一個3行4列的矩陣size(A)返回[3 4]numel(A)返回矩陣里總的元素個數(shù)也就是行數(shù)乘列數(shù)3×4的矩陣返回12length(A)則有點坑它返回的是最大維度的值對于3×4矩陣length返回4如果這個矩陣恰好是100×2length居然返回100。我當時就犯過一個大意要用一個矩陣的行數(shù)結果用length(A)去拿遇到列數(shù)比行數(shù)多的情況行數(shù)拿錯后續(xù)整個循環(huán)全亂套。后來我給自己定了規(guī)矩——要行數(shù)就用size(A,1)要列數(shù)就用size(A,2)要總數(shù)就用numel(A)盡量不要用length只有在明確知道矩陣是向量時才用它。這個習慣看起來小但能少踩很多坑。把矩陣想象成一棟房子行是樓層列是房間號。只知道“最長的那一排”并不能讓你知道房子到底有幾層。size才是那個告訴你完整結構的尺子。2. 矩陣創(chuàng)建手動輸入與快捷函數(shù)2.1 方括號、空格、逗號和分號的使用規(guī)則創(chuàng)建矩陣最直接的方式就是使用方括號。行元素之間用空格或逗號分隔行與行之間用分號分隔。例如A [1 2 3; 4 5 6; 7 8 9]這個命令生成一個3×3矩陣。如果你想創(chuàng)建行向量可以寫作v [1 2 3]創(chuàng)建列向量可以寫作v [1; 2; 3]或者寫成v [1 2 3].。我特別喜歡用分號創(chuàng)建因為一眼能看出矩陣形狀。這里有一個細節(jié)如果你在命令行運行命令時結尾沒有加分號MATLAB會立刻把結果打印在屏幕上。這在學習時很有用但如果你用腳本跑一個很大的矩陣忘記加分號會把整個命令窗口刷爆運行速度也會變慢。我的習慣是所有創(chuàng)建矩陣的語句后面都加分號需要查看時單獨運行變量名。2.2 冒號表達式與linspace生成等差向量MATLAB里最靈活的向量生成工具是冒號表達式基本格式是start:step:end。比如1:2:9生成[1 3 5 7 9]10:-2:2生成[10 8 6 4 2]如果中間不寫步長默認步長為1比如1:5生成[1 2 3 4 5]。這里有個容易忽略的地方如果start大于end而步長又是正數(shù)得到的是空矩陣。我見過不少初學者寫5:2:1結果得到一個空數(shù)組然后在后面索引時直接報錯。如果你需要生成一個固定數(shù)量的向量而不是靠步長計算可以用linspace。linspace(0,1,5)生成從0到1之間均勻分布的5個點0 0.25 0.5 0.75 1。這個方法特別適合畫圖時構造橫坐標比冒號表達式更不容易出錯。還有一個類似的logspace用于生成對數(shù)等分的向量比如做頻響分析時常用。我建議初學者優(yōu)先掌握linspace因為它不依賴步長的精度計算數(shù)量由你直接指定思路更直觀。2.3 zeros、ones、eye與隨機矩陣預分配和初始化創(chuàng)建特殊矩陣的函數(shù)非常多但最核心的四個是zeros、ones、eye和rand。zeros(m,n)生成m×n全零矩陣ones(m,n)生成全一矩陣eye(n)生成n×n單位矩陣即對角線為1、其余為0rand(m,n)生成0到1之間均勻分布的隨機數(shù)矩陣randn(m,n)生成服從標準正態(tài)分布的隨機數(shù)矩陣。為什么要單獨提zeros因為它在編程中經常用來“預分配”。如果你要寫一個循環(huán)把結果逐個存進矩陣最忌諱的做法是每次循環(huán)讓矩陣自動增長比如A(end1) x這樣每增加一個元素MATLAB都要重新分配內存循環(huán)次數(shù)一多運行時間會成倍增加。正確做法是先把矩陣尺寸定好用zeros初始化然后循環(huán)里直接賦值。比如要存100個計算結果就預先寫result zeros(100,1)這是MATLAB性能優(yōu)化的基本功。2.4 矩陣拼接和重復cat、repmat與repelem創(chuàng)建矩陣之后經常需要把兩個矩陣拼在一起。橫向拼接可以用[A B]縱向拼接用[A; B]但更規(guī)范的寫法是用cat函數(shù)C cat(1, A, B); % 按行方向拼接 D cat(2, A, B); % 按列方向拼接cat(1,...)相當于vertcatcat(2,...)相當于horzcat。對于三維數(shù)據(jù)還可以用cat(3,A,B)在第三維上拼接這在處理多通道圖像時很常見。有時候需要把一個小矩陣重復成大矩陣這時候用repmat和repelem。repmat(A, 2, 3)把A在行方向重復2次、列方向重復3次相當于整體復制。repelem(A, 2, 3)則是將每個元素各重復2次和3次。舉個例子A [1 2; 3 4]; repmat(A, 1, 2) % 結果[1 2 1 2; 3 4 3 4] repelem(A, 1, 2) % 結果[1 1 2 2; 3 3 4 4]兩者用途完全不同。需要把一組數(shù)據(jù)平鋪到更大區(qū)域時我會用repelem需要把整個模式復制多份時用repmat。3. 矩陣索引訪問元素與子區(qū)域3.1 下標索引與線性索引兩個角度訪問同一個數(shù)訪問矩陣元素有兩種方式下標索引和線性索引。下標索引就是指定行和列比如A(2,3)表示第2行第3列的元素。線性索引則只用一個數(shù)字表示按“列優(yōu)先”的順序數(shù)過去。所謂列優(yōu)先就是先數(shù)完第一列再數(shù)第二列依次往后。例如A [10 20 30; 40 50 60]; A(3) % 結果是20因為存儲順序是10, 40, 20, 50, 30, 60很多初學者在這里卡住覺得A(3)應該是第三個數(shù)也就是“第一行第三列”的30但實際是20。原因就是MATLAB在內存中把矩陣按列來排列。如果想從行列下標得到線性索引可以用sub2ind反向用ind2sublinearIndex sub2ind(size(A), 2, 3); % 求第2行第3列的線性索引 [row, col] ind2sub(size(A), 5); % 線性索引5對應的行列有人說用下標不就行了嗎為什么還要理解線性索引因為很多內置函數(shù)返回的位置是線性索引比如后面要說的find如果不理解列優(yōu)先拿到索引值后根本對應不上元素。3.2 冒號用法與end關鍵字提取行、列和子塊提取矩陣的某一行、某一列幾乎離不開冒號。A(:,1)提取第一列A(1,:)提取第一行A(1:2,2:3)提取左上角2×2子塊。這里的冒號表示“所有”。end是一個Magic Word它表示“最后一個索引位置”。比如A(end,:)取最后一行A(end-1, end)取倒數(shù)第二行最后一列。end在循環(huán)和動態(tài)編程里特別有用。我記得自己最早寫代碼時要往矩陣末尾追加一行總是先數(shù)矩陣有多少行再賦值。后來發(fā)現(xiàn)直接寫A(end1,:) newRow就行這個操作會在原矩陣下方新增一行。這種寫法不僅簡潔還避免了手動維護行數(shù)出錯。3.3 邏輯索引用條件直接篩選元素邏輯索引是我覺得MATLAB最優(yōu)雅的地方之一。你可以用一個真/假矩陣作為索引只有對應位置為true的元素才會被選中。例如A [1 5 3; 2 8 4]; logicalIndex A 4; result A(logicalIndex);result會返回一個列向量里面是所有滿足大于4的元素。這等價于兩步合并成一步A(A4)。這個寫法在處理數(shù)據(jù)清洗時極好用比如把矩陣里所有負數(shù)清零直接寫A(A0)0。需要注意邏輯索引的矩陣維度必須和原矩陣相同或者維度可以被隱式擴展。邏輯索引的選擇結果永遠是列向量形式因為MATLAB會把選出來的元素按線性順序依次排列。如果你希望保留行列位置可以先轉成下標也就是用find。3.4 修改矩陣賦值、刪除與擴展索引不但可以讀還可以寫。A(2,2)100修改單個元素A(:,2)[10;20;30]整列替換A(2,:)[]刪除第2行注意這里用空方括號[]作為賦值內容時MATLAB會直接刪掉相應行同時矩陣行數(shù)減一。刪除操作在循環(huán)里要格外小心。如果你在for i1:size(A,1)循環(huán)中刪除某一行循環(huán)次數(shù)不會動態(tài)變化那么跳過某些行或者索引越界幾乎是必然的。我的經驗是如果要做逐行篩選刪除最好先通過邏輯索引一次性選出要保留的行或者倒著循環(huán)。比如從后往前刪for i size(A,1):-1:1這樣索引不會亂。4. 矩陣運算數(shù)組運算與線性代數(shù)運算4.1 最關鍵的區(qū)分A*B與A.*B矩陣運算中最容易混淆的就是乘號。A * B是線性代數(shù)中的矩陣乘法要求A的列數(shù)等于B的行數(shù)結果第(i,j)個元素是A第i行與B第j列的點積。A .* B是對應元素相乘要求A和B維度完全一致結果矩陣的每個元素是A(i,j)*B(i,j)。我見過太多新手用*去做對應元素相乘結果要么維度不匹配報錯要么得到完全錯誤的結果。記住一個口訣帶點的運算是在說“按位置算”不帶點的乘除法是線性代數(shù)里的矩陣運算。這個區(qū)分同樣適用于除法和冪次A./B是逐元素相除A/B是右除A.^2是每個元素求平方A^2是矩陣自乘。從物理意義上想.*更像是對照片每個像素做亮度調整*則更像坐標系變換兩者服務于完全不同的場景。4.2 轉置與共軛轉置一個撇關系的兩個世界在MATLAB里單引號是共軛轉置點號加單引號.才是普通的非共軛轉置。當矩陣是實數(shù)時兩者結果一樣一旦涉及復數(shù)結果就完全不同了。例如Z [12i, 3-4i]; Z % 得到 [1-2i; 34i] Z. % 得到 [12i; 3-4i]如果你在做傅里葉變換、濾波器設計這類復數(shù)運算忘記加地方實會用錯。建議在沒有復數(shù)參與的時候也統(tǒng)一使用.這樣代碼更安全讀代碼的人也不會誤解你的意圖。這是一個非常小但很專業(yè)的習慣。4.3 左除和右除解線性方程組的最佳姿勢很多人初學MATLAB時會找逆矩陣inv(A)來解線性方程Axb書上有時也這么寫。但實際做數(shù)值計算時更好用的是左除運算符\x A \ b;它的意思是求解線性方程組A*x b。MATLAB內部會根據(jù)A的特點自動選擇高斯消去、LU分解、Cholesky分解等方法比直接用inv(A)*b更快、更穩(wěn)定。右除/則對應x*A b形如的方程比如x b / A。你可能會問為什么不推薦inv因為求逆本身計算量更大而且當矩陣接近奇異時數(shù)值誤差會被放大。用\相當于讓MATLAB悄悄選一條更安全的路徑。我平時幾乎不寫inv除非是要在公式里顯式寫出逆矩陣。4.4 常用函數(shù)速查從行列式到特征值矩陣基礎必須掌握的函數(shù)大概有這么幾類線性代數(shù)類det行列式、rank秩、trace跡、eig特征值和特征向量、svd奇異值分解、norm范數(shù)、cond條件數(shù)統(tǒng)計類sum求和、mean均值、median中位數(shù)、std標準差、max最大值、min最小值、sort排序數(shù)組操作類diff差分、cumsum累積和、flip翻轉、unique去重。舉個例子A [1 2 3; 4 5 6; 7 8 10]; sum(A) % 默認對每一列求和結果[12 15 19] sum(A, 2) % 對每一行求和結果[6; 15; 25] max(A, [], 1) % 對每一列取最大值 max(A, [], 2) % 對每一行取最大值注意max和min的寫法很特別想按維度求值時要傳一個空數(shù)組[]作為第二個參數(shù)否則max(A)會把整個矩陣當成一個向量來比較。這個小細節(jié)踩一次就記住了。5. 矩陣形變與批量處理技巧5.1 reshape的列優(yōu)先陷阱與permute的正確用法reshape可以把矩陣改成任意維度只要元素總數(shù)不變。例如A 1:12; B reshape(A, 3, 4);結果B是3行4列但順序是列優(yōu)先填充的第一列是1、2、3第二列是4、5、6如此類推。如果你想按行優(yōu)先的順序重新排布必須先對原數(shù)據(jù)做一些處理。一個常見需求是你有一個按行采集的信號序列想把它重排成多行多列但reshape默認按列排可能會打亂信號順序。這時可以用轉置或permute來解決。比如把行向量1:12變成3行4列且按行順序排列可以先reshape再轉置B reshape(A, 4, 3);permute是高維版的“轉置”適合處理三維矩陣的維度交換。比如permute(M, [2 1 3])會把前兩維互換。對于初學者我建議先用permute做實驗觀察維度變化理解后再用于更復雜的數(shù)據(jù)處理。5.2 用find與邏輯矩陣做批量篩選在數(shù)據(jù)處理中經常需要知道哪些元素滿足條件以及它們的位置。find是解決這個問題的利器。A [10 20 30; 40 50 60]; [row, col] find(A 25); linearIndex find(A 25);不帶賦值時find(A 25)返回滿足條件的線性索引向量。如果你寫成[row, col] find(...)則返回行列下標。這個行為與邏輯索引互補邏輯索引直接給你值find給你位置兩者結合可以高效地對滿足條件的位置做修改。比如你想把矩陣中所有大于100的替換為100可以直接用A(A100)100但如果你想以某個特定順序逐個修改這些位置那就需要先拿到線性索引再循環(huán)處理。5.3 隱式擴展自動廣播維度的新特性從R2016b開始MATLAB支持隱式擴展讓不同維度的矩陣在操作時自動擴展。比如A [1 2 3; 4 5 6]; b [10; 20]; C A b;這里A是2×3b是2×1。按照老版本規(guī)則維度不一致會報錯新版本里b會沿著列方向自動復制成2×3然后再相加。同理A - [1 2 3]會把行向量自動復制到每一行。這個功能很方便但也帶來隱患如果代碼里有維度不匹配的bug不再報錯而是默默做廣播產生一個很大但很蠢的結果。我的建議是利用隱式擴展時務必先用size確認原來矩陣的形狀特別是做矩陣和向量混合運算時最好心里明確哪個維度會被擴展。5.4 利用all、any和邏輯運算簡化判斷當你想檢查矩陣中是否所有元素滿足條件、是否存在某些元素滿足條件all和any是必用的。A [1 2 3; 4 5 6]; all(A 0) % 返回 [1 1 1]因為每列所有元素都大于0 any(A 5, 2) % 對每一行判斷返回 [0; 1]配合方向參數(shù)可以處理行和列兩個維度。比如你想篩選出所有行中至少含有一個大于5的行可以用rows any(A 5, 2); selected A(rows, :);這種寫法比循環(huán)簡潔得多也是“向量化思維”的典型例子。矩陣基礎學得好不好一個很重要的標志就是你寫代碼時還會不會動不動就寫for循環(huán)。6. MATLAB矩陣編程常見問題與避坑6.1 索引越界和維度不一致我見過最多的報錯就是“Index exceeds the number of array elements”和“Matrix dimensions must agree”。索引越界通常是因為你用了固定大小去訪問動態(tài)矩陣比如在循環(huán)里只判斷了行數(shù)沒判斷列數(shù)。解決這類問題的關鍵是每次操作前想清楚矩陣當前的實際大小可以用size、numel確認。維度不一致則多發(fā)生在點運算時兩個矩陣需要對位相乘結果一個2×3一個3×2直接報錯。遇到這種情況我通常會立刻檢查兩個矩陣的size輸出重點看是不是轉置忘了寫。這類錯誤在腳本運行中很常見但定位并不難難的是你是否養(yǎng)成了先查看維度的習慣。6.2 循環(huán)中動態(tài)增長矩陣導致性能驟降一個反例x 1:10000; y []; for i 1:length(x) y(i) x(i)^2; end這段代碼本身沒有錯但每次循環(huán)都要給y重新分配內存運行會特別慢。正確的寫法是y zeros(size(x)); for i 1:length(x) y(i) x(i)^2; end更好的是直接向量化y x.^2;用tic和toc測一下你會發(fā)現(xiàn)向量化比循環(huán)快幾個數(shù)量級。這不是鼓勵初學者一開始就拒絕循環(huán)而是提醒你當數(shù)據(jù)量大了矩陣化思維才是MATLAB的看家本領。6.3 浮點數(shù)比較別用等號由于浮點數(shù)的二進制表示有誤差A中某個數(shù)值即使看起來等于0.1實際存儲可能是0.10000000000000001。因此判斷矩陣中哪些元素等于某個特定小數(shù)時不建議用A 0.1而應該用tol 1e-8; A(abs(A - 0.1) tol) 0;這個技巧在做數(shù)值計算時非常重要尤其是在處理迭代法收斂判斷、濾波器系數(shù)比較時。直接等號比較往往找不齊目標元素或者誤判數(shù)值導致結果偏離預期。6.4 變量名與內置函數(shù)沖突矩陣基礎知識雖然簡單但有一個特別隱蔽的坑變量名會覆蓋函數(shù)名。比如你定義一個變量叫zeros之后調用zeros(3,3)就會出錯。因為MATLAB解析變量名的優(yōu)先級高于函數(shù)。我的建議是永遠不用內置函數(shù)名作為變量名比如不要用sum、max、mean、size。寫代碼時可以先運行which size檢查某個名字是否被占用。還有一個容易忽略的點MATLAB對大小寫敏感。A和a是兩個不同的變量。在長腳本里大小寫混亂很容易導致“Undefined function or variable”錯誤。統(tǒng)一駝峰或小寫命名能減少這類問題。6.5 分號、清屏與腳本調試的細節(jié)命令窗口直接執(zhí)行語句時如果想看到結果就不加分號如果想隱藏結果就加分號。在腳本編輯器中建議默認加避免中間過程輸出刷屏。調試矩陣相關代碼時我通常會這樣操作在關鍵行前面去掉分號運行一次或者用disp打印size(A)也可以用whos查看所有變量的大小和類型這個命令比who更詳細。如果發(fā)現(xiàn)某個變量維度不對可以用clear A單獨清除它而不是動不動就clear all否則腳本里其他變量全沒了前功盡棄。最后再分享一個小技巧在命令窗口輸入doc size或者help reshape能直接打開官方文檔里面有很多示例比網(wǎng)上搜到的二手答案靠譜得多。我當年學矩陣基礎時就是靠“遇到不懂就doc”這個習慣把repmat、reshape、logical indexing這些概念真正用熟的。矩陣基礎這塊內容看似瑣碎但它是后續(xù)所有MATLAB應用的地基哪怕你后面去做圖像處理、數(shù)值分析、深度學習理解矩陣怎么存、怎么取、怎么算永遠繞不開。