NumPyで最小二乗アルゴリズムを効率的に計算する

Question

最小自乗法で大きな線形システムを解く必要があります。私はnumpy.linalg.lstsq(a, b),np.dot(np.linalg.pinv(a), b)と数学的な実装の計算効率の違いを理解するのに苦労しています。

私は次の行列を使用：

h=np.random.random((50000,100)) 
a=h[:,:-1].copy() 
b=-h[:,-1].copy() 

およびアルゴリズムの結果は次のとおり

---

# mathematical implementation 
%%timeit 
np.dot(np.dot(np.linalg.inv(np.dot(a.T,a)),a.T),b) 

10のループ、3の最もよい：36.3 MSループ

あたり

---

# numpy.linalg.lstsq implementation 
%%timeit 
np.linalg.lstsq(a, b)[0] 

10ループ、3の最もよい：103ミリループ当たり

---

%%timeit 
np.dot(np.linalg.pinv(a), b) 

1ループ、3の最もよい：216ミリループ当たり

---

ある理由違いは？

Answer 1

ルーチンlstsqは、過度に決定された、不十分な、または十分に決定されたシステムを処理します。その出力はpinv（a）* bから得られるものですが、疑似逆を計算するよりも速いです。理由は次のとおりです。

一般的なアドバイス：実際に必要でない限り逆行列を計算しないでください。特定の右辺のシステムを解くことは、その行列を反転するよりも速いです。

しかし、^T解決を経由して、あなたのアプローチは= ^T Bは、あなたが行列を反転しているにもかかわらず、高速です。何がありますか？問題は、^Tを反転すると、aが完全な列ランクを持つ場合にのみ有効です。したがって、この特定の状況に問題を限定し、より一般性の低いトレードオフとして得られ、以下に示すように、安全性が低いためです。

しかし、行列を反転することは依然として非効率的です。あなたがフルランクであることがわかっている場合は、次はあなたの3回のいずれよりも高速である：悪い条件の行列を扱うとき

言っ

np.linalg.solve(np.dot(a.T, a), np.dot(a.T, b)) 

---

、lstsqは、上記にまだ望ましいです。製品を形成すると、^Tが基本的に条件番号を二乗するので、無意味な結果が得られる可能性が高くなります。ここで（numpyののと本質的に同等であるが、複数の方法を有する）scipyのダウンロードのlinalgモジュールを使用して、注意の例である：

import numpy as np 
import scipy.linalg as sl 
a = sl.hilbert(10) # a poorly conditioned square matrix of size 10 
b = np.arange(10)  # right hand side 
sol1 = sl.solve(a, b) 
sol2 = sl.lstsq(a, b)[0] 
sol3 = sl.solve(np.dot(a.T, a), np.dot(a.T, b)) 

ここlstsqほぼsolveと同じ出力（本システムの一意の解）を得ます。しかし、sol3は数値的な問題のために全く間違っています（これは警告さえありません）。

SOL1：

[ -9.89821788e+02, 9.70047434e+04, -2.30439738e+06, 
    2.30601241e+07, -1.19805858e+08, 3.55637424e+08, 
    -6.25523002e+08, 6.44058066e+08, -3.58346765e+08, 
    8.31333426e+07] 

SOL2：

[ -9.89864366e+02, 9.70082635e+04, -2.30446978e+06, 
    2.30607638e+07, -1.19808838e+08, 3.55645452e+08, 
    -6.25535946e+08, 6.44070387e+08, -3.58353147e+08, 
    8.31347297e+07] 

SOL3：

[ 1.06913852e+03, -4.61691763e+04, 4.83968833e+05, 
    -2.08929571e+06, 4.55280530e+06, -5.88433943e+06, 
    5.92025910e+06, -5.56507455e+06, 3.62262620e+06, 
    -9.94523917e+05]