最小二乘法在很多地方都會用到。比如形狀公差中，如果設計者採用下面的標註：
圖1 圖紙標註
圖1中的標註，平面度的理想要素就必須採用最小二乘法擬合出來（如果沒有那個G的修飾符號，默認是最大值最小法擬合平面）。還有，很多基準的擬合，也會近似的採用最小二乘法擬合。（GZHl:智慧汽車供應鏈）
之前在德輝學堂介紹過最小二乘法，但是有很多好學的小夥伴總是追問，最小二乘法的數學公式究竟是怎麼樣的？
本期的這一篇文章，我們將介紹一個簡潔的最小二乘法數學公式，慢慢剖析它，爭取讓好學的小夥伴們能認識它，然後再結合Excel利用它來做一些計算。
本期文章將分3個部分來講解：
1. 3個方程組
2. 矩陣的初步概念
3. 最小二乘法的公式和案例
本期的文章可能會讓很多小夥伴感到頭疼不適，因為涉及的數學概念較多（主要是線性代數），如果你對線性代數比較了解，可以直接跳到第三節，如果不太了解，建議耐心把它看完。我們儘量從實用的角度出發，避重就輕，盡最大限度幫助小夥伴們快速理解這個公式的含義並學會使用它。
1. 3個方程組
1.1 不定方程組
方程組，我們初中就學過解方程。
我們知道，平面上的直線方程為：
圖2 直線方程
比如有一天，老師出的題目是這樣的：已知平面中的一個點P（x0, y0)，求過這個點P的直線方程？（GZHl:智慧汽車供應鏈）
注意，直線方程有兩個重要的參數，k和b，只要我們知道這兩個參數，這個直線方程就確定了。所以我們把這個點P（x0,y0)帶入直線方程（1）可得：
上式其實就是一個方程組（只不過方程組裡只有一個方程），這裡的未知數是k和b。如果根據方程組（2）能夠求出k和b, 那麼直線方程就求出來了。
問題是，兩個未知數，一個方程組，怎麼解？顯然在方程組(2)中， k和b的解有很多個，不唯一。想像一下就可以理解，平面上，過一個點的直線有無數條，所以k和b的解有很多個。
圖3 過一個點有無數條直線
可以看出來，方程組（2）的特點是，方程的個數（也就是已知點的個數）小於未知數的個數（2個）。這樣的方程組，我們叫不定方程組，它沒有唯一解，有無數個解。
不定方程組對我們的意義在於，如果你用三坐標只採一個點，在沒有其它約束的條件下，三坐標是沒有辦法擬合出一條直線的，因為有無數條直線的可能，三坐標會直接懵逼出錯。
當然，不定方程組，有一個唯一的最小範數解，和我們的主題無關，不去管它。
1.2 正定方程組
如果在平面上，有兩個已知點：P1(x1, y1), P2(x2,y2), 求過這兩個點的直線方程？相信初中班的小傢伙們就可以解出來的，因為直接把兩個已知點帶入直線方程（1），就可以得有兩個方程的方程組：（GZHl:智慧汽車供應鏈）
根據方程組（3），我們可以直接解出k和b, 再次強調， k和b是未知數。只要P1和P2不在同一個點上，方程組（3）一定有唯一解k和b。從幾何上講，過P1和P2兩點一定有一條唯一直線。（GZHl:智慧汽車供應鏈）
來看看圖吧：
圖4 過兩個點能確定一條唯一的直線
方程組（3）是我們上學的時候經常做的習題，在做作業的時候，只要未知數的個數和方程的個數相等，我們就可以把它給解出來。（GZHl:智慧汽車供應鏈）
這種方程的個數和未知數的個數相等的方程組，我們把它叫住正定方程組。它的特點是有唯一確定解（當然這是有前提條件的，比如圖4中的兩個點不能在同一個位置上）。
在測量的時候，我們是經常有意識或者無意識的在用正定方程組，比如建基準系的時候，所謂的「3-2-1」原則中的「3」，指的是第一基準建立的時候，要在一個平面上採3個點（3個點不在同一條直線上），就能確定一個唯一的平面。因為平面方程是
顯然平面方程有3個未知的參數a,b,c，三坐標採了3個點，軟體就能列出有3個方程的方程組，也就是正定方程組。然後解出出a ,b,c這三個參數，得到平面方程，然後就可以想幹嘛就幹嘛。
「2」用在第二基準的建立，要求在產品表面上採2個點。它也是一個正定方程組，因為除了過兩個點（有了兩個方程），還要加一個和第一基準的垂直約束（第三個方程），也能列出一個正定方程組。
「1」也是一樣，過一個點，加上和第一基準，第二基準，兩個垂直的約束，也是一個正定方程組。
站在設計工程師的角度來說，我個人是懷疑3-2-1原則的合理性的，這個原則可以使軟體的計算變得簡單，但是和實際裝配的符合性不好。實際上，我們希望測量工程師在第一，第二，第三基準上的採的點越多越好，才和實際裝配接近。
如此一來，我們就不得不面臨一個苦逼的問題，超定方程組。
1.3 超定方程組
還是用直線舉例，假如測量工程師在一個面上採了3個點P1(x1,y1)，P2(x2,y2), P1(x3,y3)要我們求出過這3個點的直線方程。根據我們前面套路，我們可以把這3個點分別帶入直線方程（1），得到一下方程組：
這時，你憑直覺，是不是有一種不好的預感呢？未知數k和b可能沒有解, 因為方程組（4）中很有可能兩個方程之間的關係是矛盾的。大家可以想像一下，如果P1, P2, P3這3個點不在同一條直線上，如何去求過這3個點的直線呢？
圖4 過三個點能確定一條唯一的直線？（GZHl:智慧汽車供應鏈）
顯然，去尋找一條同時過這3個點的直線肯定是扯，因為它根本就不存在。也就是說，不存在一個k和b, 使得三個點都在同一條直線y=kx+b上。
所以方程組（4）的特點是，它的方程的個數大於未知數的個數，這樣的方程組，我們把它叫住超定方程組。
超定方程組，聽起來是不是很牛逼的樣子？
事實上是苦逼，因為它沒有解。（GZHl:智慧汽車供應鏈）
更苦逼的是，在我們現實的測量中（甚至公差分析中），面臨的都是超定方程組。比如，在一個實際平面上採8個點（放心，這8個點一定不在同一個平面上），非要讓你解出一個平面，你面臨的就是一個解超定方程組的問題。
比如，前面講的直線，變態的是，測量工程師可能採了12個點，需要軟體算出一條直線，軟體怎麼辦？
所以在測量的時候，需要解超定方程組簡直就是家常便飯。（GZHl:智慧汽車供應鏈）
明明知道沒有解，非要去解，那怎麼辦呢？沒有辦法，日子總得過下去，我們只能找出一個近似解來「將就」一下。
上帝並沒有讓我們完全絕望，超定方程組有個特點，那就是它有唯一的最小二乘解。
最小二乘解就是近似解的一種。最小二乘法的特點在本公眾微信號的《來認識傳說中的最小二乘法》中已經介紹過，有興趣的小夥伴，可以點擊文章最後的連結。
如何求得這個最小二乘解呢？它有固定的公式，這正是本篇文章的目的。在介紹最小二乘法的公式之前，我們還有個難題必須要克服，那就是要認識矩陣和矩陣相關的概念，我們馬上進入第2節。（GZHl:智慧汽車供應鏈）
2. 矩陣的初步概念
說到矩陣，可能有一半的小夥伴被嚇跑了。
矩陣本身是個非常好用的數學工具，大約是因為上學那會兒為了應付考試，沒有看到它的真正用處，可能很多小夥伴和我一樣開始對它有心理陰影了。
本篇文章儘量只簡單介紹和最小二乘法數學公式相關的矩陣的3個概念，1. 矩陣的概念和矩陣的乘法，2. 矩陣的轉置，3. 矩陣的逆矩陣。
2.1 矩陣概念和矩陣的乘法
矩陣，其實就是表格裡的一堆數。我們上學那會兒，學的123，xyz等等都是單個數單個數來認識的，然後進行處理和運算。而矩陣，則是多個數，放在一起，很簡單。
舉個例子，共享單車的興起，據說救活了很多瀕臨倒閉的自行車廠，因為他們收到了很多訂單。比如下面是某自行車的兩個分廠2018年供給各共享單車品牌的銷售清單，單位萬輛。
以上的這個自行車銷售清單的表格，我們就可以用矩陣來表達, 給矩陣取個名字叫A， 就有：
上圖中用中括號表達的這個東西就是高大上的矩陣A, 它就包含了表1中的所有信息。它是不是就是一張表格？所以它就像表格一樣，可以有很多行，也可以有很多列。
我們再來，因為每個品牌單車的樣子和型號不一樣，所以自行車單價也就不一樣，下面是每個品牌單車的單價清單（單位RMB）。
上面這個表格，也可以寫成矩陣的形式, 給矩陣取個名字叫B. 則有：
所以，矩陣其實就是把一堆數字，按特定順序放在特定位置的一張表格裡。
好了，我們對矩陣有了一個初步的認識，我們再來粗粗的認識一下矩陣的乘法。
繼續我們的故事，到2018年年底，各品牌共享單車一分錢都沒有給自行車廠家，流言傳來，說共享單車不行了，這時自行車廠家慌了，趕緊核算一下各分廠這一年的損失。
這個算法小學生都會做的，我們把表格放在一起，很容算出來杭州分廠和上海分廠的損失（數量x單價=金額）。
圖5 自行車廠家的損失
圖3中第三張表格，顯示的就是自行車廠家損失的計算方法，損失的結果也可以用矩陣C來表達，則有：
好了，從數學上講，我們還可以用矩陣的乘法來表達上面的結果，很簡單：
AB=C
即為：
公式（5）就是矩陣的乘法，這個就像我們平時用的加減乘除一樣，就是一種算法，通過矩陣的乘法，我們就可以直接計算出自行車廠家各分廠的損失。
大家最細要仔細觀察等式（5）和圖5中第三張表中的11960和10720是怎麼來的。矩陣乘法，對初學者來說，稍微麻煩一點，要記住結果的規則，我列一個公式：
等式（5）的結果就是根據公式（6）計算出來的，公式（6）就是矩陣乘法的遊戲規則。有興趣的小夥伴最好記住，因為很多時候，我們需要根據方程組，也就是等式（6）的右邊的樣子，反推回等式左邊的樣子的。
大家注意到沒有，矩陣的好處是可以「批處理」同類型的計算（當然，矩陣的好處遠遠不止這些）。
如果你實在記不住，也不願意記，怎麼辦？沒有關係，用Excel裡的一個函數MMULT（）它可以直接幫你計算，你記住這個函數就好了。見圖6：
圖6 利用Excel的函數計算矩陣乘法
2.2 矩陣的轉置
矩陣裡有行，也有列，如果同行變成同列就是轉置（相當於翻轉）, 比如把第一行變成第一列，第二行變成第二列...。矩陣A的轉置後叫轉置矩陣，轉置矩陣就記著：
舉個例子吧，假設有：
把A矩陣的第一行變成第一列，第二行變成第二列，就可以得到A的轉置矩陣。那麼矩陣A的轉置矩陣為：
所謂轉置是不是很easy?
還是怕弄錯？好吧，沒有關係，Excel裡邊有個函數Transpose()可以幫你完成矩陣的轉置。
2.3 矩陣的逆
矩陣的逆，我們不深入了解，暫時把它想像成矩陣的倒數（其實矩陣沒有這個倒數這個概念）。
任何一個非零數乘以它的倒數，結果是1，而這個1其實就是一個單位。所以有時候我們用乘以倒數的方法來代替除法。
矩陣裡邊也有類似1一樣的單位矩陣，比如：
上面的矩陣E就是一個2階的單位矩陣。所謂逆矩陣的特點，就是任何一個矩陣乘以它的逆矩陣等於單位矩陣。比如：
根據前面講的矩陣的乘法規則，我們會發現：
所以，B就是A的逆矩陣，A也是B的逆矩陣。記為：
這裡要強調一下，只有方的矩陣（方陣），即行數等於列數的矩陣才有「逆」這一說，而且就算是方的矩陣， 也不一定有逆矩陣（它的充要條件是行列式不等於0），類似不是任何數都有倒數的原理一樣（比如0就沒有倒數）。
問題是，如果我已經知道了一個矩陣，如何求它的逆矩陣呢？其實手工求比較麻煩，還要引入伴隨矩陣之類的概念，不去深入了。我們還是偷懶，直接用Excel吧。
Excel裡邊有一個函數叫Minverse(), 它直接可以求出逆矩陣。如果逆矩陣不存在，他會提示出錯。
3. 最小二乘法的公式和案例
好了，我們做了很多鋪墊，講了超定方程組，講了矩陣和矩陣的乘法，轉置矩陣，逆矩陣等概念，千呼萬喚，現在終於可以把最小二乘法的公式給牽出來了。
定義對於超定方程（矩陣方程）
如果滿足
可逆（這個條件大家不用擔心，只要用三坐標採的點不再同一個位置上，肯定是可逆的），那麼該超定方程（7）有唯一最小二乘解，且解的公式為：
這個公式（8）就是我們今天的主角。
需要說明的是，公式（7）中，A, X, b他們都是矩陣（所以它們被加粗），A和b是已知的係數或者值，x是不知道的未知數矩陣（如果你在這裡開始懵圈不用急，後邊我們還會舉例來認識這個公式）。
公式（8）就是最小二乘法的解法，看起來很複雜的樣子，但是對我們來講，只要知道超定方程（7）中的係數矩陣A和b, 利用公式（8）就可以輕鬆解出最小二乘解x.
這個公式怎麼來的？為什麼成立？這解釋起來太痛苦，必須要具備比較全面線性代數的知識，從基到張量空間再到坐標去解釋，為了不讓最後幾個讀者跑掉，不再講了。
反正，有興趣的小夥伴可以去查閱本文後邊的文獻。不過我可以把手按在手機上發誓，這公式不是我杜撰的，而且確定是正確的，算出來的這個x一定能滿足平方和最小。
還要補充一下，公式（7）是一個矩陣方程，而我們前面講的是線性的超定方程組，所以我們要學會自己把方程組轉化成矩陣方程。舉個例子，前面的超定方程組為：
稍微做一個變化：
這樣就可以轉化成矩陣方程了：
上面這個矩陣可能有點突然，如果大家把方程（10）中等號左邊的矩陣乘法按照乘法的遊戲規則（6）打開，就是方程組（9）了。
這一點還比較重要，有興趣的小夥伴最好能熟練變換。
矩陣方程（10）中的x1,y1,x2,y2,x3,y3都是已知點的坐標數據，k和b是未知數。我們再給每個矩陣取個名字。令：
則矩陣方程（10）就可以寫成：
看到沒有，上面這個就是把超定方程組（9）轉化成標準的矩陣方程了，它的好處在於我們可以直接套公式。所以，我們能再利用公式
把k,b給解出來。
我就把它的長公式寫出來吧。
公式（11）看起來很複雜，其實右邊那一堆就像加減乘除一樣，最終得出一個2行1列的矩陣。如果你耐心把它輸入Excel裡邊，也就是一步計算的事啊。
到了這裡，k和b的最小二乘解就解出來了。
為了能夠Excel來幫我們計算，我們再回憶一下Excel中的三個處理矩陣的函數：
如果你已經有了矩陣方程Ax=b, 也就是說你有已知矩陣A和b, 那麼用Excel計算最小二乘解x就是分分鐘的事情。你只要在Excel裡邊輸入一個長長的公式就可以解決。當然你也可以一步一步輸入公式計算。
數學公式為：
利用Excel的函數公式為：
x=MMULT(MMULT(MINVERSE(MMULT(TRANSPOSE(A),A)),TRANSPOSE(A)),b)
這個公式耐心看，吃透，它和數學公式是一模一樣的。如果你還是不明白，那就拷貝黏貼吧。
好了，到這裡，我們還是舉一個實際的例子吧。我在平面上採了5個點，坐標如下圖所示，請把這5個點擬合成最小二乘直線，並寫出直線方程：
圖7 5個點的分布
已知直線的方程為y=kx+b, 將A,B,C,D,E五個點帶入直線方程，可以得到一個超定方程組（有些文獻也把它叫做觀測方程）：
將超定方程組（12）整理成矩陣方程，如下：
顯然有：
然後可以解得x的最小二乘解為：
所以可以得到直線方程的參數
k=0.731, b=0.551。
即最小二乘直線的方程為：
y=0.731x+0.551
Excel中的表格計算和公式輸入，見圖8：
圖8 Excel表格計算
我們再把最小二乘直線方程 y=0.731x+0.551和原來的5個點放在一起，做一下比較。見圖9.
圖9 最小二乘直線
顯然，根據最小二乘法的規則，只有y這條直線滿足殘差的平方和e最小：
其中，也就是說，e值在所有可能的直線中，只有y=0.731x+0.551這條直線滿足e是最小的。也就是說,儘管有點「將就」（也不得不將就），但是y=0.731x+0.551的直線已經是折衷後的優選方案了。
好了，最小二乘直線是擬合出來了，幹嘛用呢？
該幹嘛就幹嘛，可以做基準去評價其它被測要素，可以做理想要素去夾被測要素（比如直線度），隨你的便。
為了加深印象，我們再舉一個例子，我們利用最小二乘法再來計算一個面輪廓度。
已知圖紙標註和實際零件如圖10所示：
圖10 圖紙和零件上採的點
如圖10，我們首先在基準要素（下表面）上採了六個點，點坐標如下：
被測要素（上表面）上採的點也有6個，數據如下：
因為平面方程可以寫成下面的形式：
z=ax+by+c （14）
公式（14）中，a,b,c是平面方程的參數，只要知道a,b,c，我們就知道最小二乘法擬合出來的基準平面了。同樣的方法，把D1，D2...D6的x,y,z坐標值分別代入平面方程（14），可以得到下面的超定方程組:
接下來的思路是如何把它轉化成矩陣方程。我們令
則方程組（14）就可以寫成矩陣方程：
Ax=b
顯然A是已知的係數矩陣（代入坐標值就已知），x包含3個未知數a,b,c， b也是一個數據已知的矩陣（所有的已知z）。我們就可以套公式啦。
我們將原始的數據整理成A和b, 然後利用Excel的函數：
x=MMULT(MMULT(MINVERSE(MMULT(TRANSPOSE(A),A)),TRANSPOSE(A)),b)
可以直接求出a,b,c。
Excel的具體數據如下：
根據表5的計算，可以得到基準A的方程是：
為了方便後邊直接套用公式，需要將上面這個基準平面的方程直接轉化標準平面方程：Ax+By+Cz+D=0, 轉化後為：
顯然，標準平面方程，我們可以得出: A=0.002, B=-0.005,C=-1, D=0.026, 這四個標準平面方程的參數在算距離的時候，馬上要用到的。
然後再求被測要素上每一個點到該基準面的距離，就可以算出輪廓度（根據ISO標準，輪廓度為和理論尺寸24差值最大的距離的2倍，）。這裡需要利用點到面的距離公式（A,B,C,D四個參數剛好可以在這裡用上）：
將被測要素每點的坐標代入公式（16），用點的實際坐標代替公式（16）中的x0,y0,z0。計算出每點到基準面的距離，最後可以計算出輪廓度（按照ISO標準）。最後的計算結果參考下面的表6。
表6 輪廓度的計算
根據上面的計算，可以得到輪廓度的實測值為0.344。
好了，我們的文章到這裡快要結束了。
稍微回顧一下思路，我們求最小二乘解的時候，思路如下：
1. 寫出目標函數，比如直線方程，平面方程；
2. 將已知點帶入目標函數，構建超定方程組；
3. 再將超定方程組轉化成矩陣方程Ax=b的形式（這一步有點難度，需要有興趣的小夥伴們反覆練習掌握）；
4. 利用公式
求出最小二乘解；
本期文章介紹的最小二乘解的公式，對線性問題適用（不僅僅是對測量GD&T, 其它地方也可以用），對非線性問題不能直接使用，如果要擬合最小二乘圓或圓柱之類的特徵就沒有辦法直接採用這個公式。對於非線性的公式，如有機會，我們以後和大家討論。
最後申明一下，本文的目的純粹是為GD&T的發燒友們準備的技普知識，涉及到相關數學概念的解釋一定有諸多不嚴謹之處，更加詳細的，嚴謹的概念解釋最好參考專業的文獻，比如本文最後的參考文獻。
本期小結
本期文章介紹了通過矩陣來計算最小二乘法的一個公式（注意，這並非最小二乘解的唯一方法），它非常簡潔，但要求讀者對線性代數有一定基礎。為了照顧到把線性代數已經還給老師的小夥伴們，我們前期做了鋪墊。第1節，講了三種方程組，說明我們在擬合的時候，面臨最多的是超定方程組。第2節，講了矩陣的一些相關知識，為了幫助小夥伴們理解公式的含義。第3節才正式介紹最小二乘法的公式和相關例子，並用輪廓度計算作為例子，來解釋了這個公式。
【後記】
寫本期文章的過程非常痛苦，一方面在抖音，快手等超短，超爽視頻流行的年代，公眾號長文註定不受歡迎，而且還是針對專業化的小眾文章；另外一方面，技普文章一旦遇到數學問題，馬上就變得誠惶誠恐，因為數學既嚴謹，水又很深，一不小心就發現自己在膚淺的胡說八道。不管怎樣，終於寫完，不為別的，只為真正有好奇心的同行，哪怕寥寥無幾。
這裡要感謝OGP的伊揚威老師，他幫我用OGP的軟體驗證最小二乘法的計算結果的正確性，還熱心的通過視頻教我如何使用這款軟體，讓我有信心把這篇文章發出來。
如果你有任何問題，或者發現本文的錯誤之處，歡迎留言給我們。我也給你留一個思考題，在本文最後一個案例中，如果要求被測面相對於基準面A面的平行度，你可以通過今天的講的公式計算出來嗎？你可以留言告訴我們哦。
另外文中案例已經做成Excel表格，歡迎加微信索要。
參考文獻
精密測量的數學方法 熊有倫 中國計量出版社 線性代數的幾何意義 任光千 謝聰 胡翠芳 西安電子科技大學出版社 線性代數與空間解析幾何 科學出版社