2016-10-24 13 views
2

私は、PythonでA = Bxという形式の線形代数問題を扱い、これをMATLABとMathematicaの同僚のコードと比較しています。 Bが特異行列であるときに、Pythonと他のものとの違いに気付きました。 numpy.linalg.solve()を使用すると、私は特異行列エラーを投げるので、代わりに.pinv()(Moore Penrose擬似逆行列)を実装しました。特異行列の方程式を扱うときにMathematicaとPythonの答えが異なるのはなぜですか?

逆行列の格納は計算効率が悪く、Pythonで特異行列を扱うより良い方法があるかどうか不思議です。しかし、私の疑問は、Pythonが無限の解空間から答えを選ぶ方法、そしてなぜMATLABとMathematica以外のものを選ぶのかということにあります。

B = np.array([[2,4,6],[1,0,3],[0,7,0]]) 
A = np.array([[12],[4],[7]]) 
BI = linalg.pinv(B) 
x = BI.dot(A) 

私にはPythonの出力がある答え:ここで

は私のおもちゃの問題である

[[ 0.4] 
[ 1. ] 
[ 1.2]] 

これが正しい答えは確かですが、それは私が意図したものではありません:(1,1,1)。なぜPythonはこの特定のソリューションを生成しますか?可能な解決策ではなく、ソリューションのスペースを返す方法はありますか?私の同僚のコードが返されました(1、1、1) - PythonがMathematicaとMATLABと異なるのはなぜですか?

+0

より便利な比較のために、代わりに 'BI'を印刷しないでください。 – Eric

+0

_ "私の同僚のコードが返されました(1,1,1)" _ - これはmatlab、mathematica、またはその両方でしたか?そのコードを表示できますか? – Eric

+0

'pinv(B)* A'は私にmatlabと同じ結果を与えます。だからあなたのMATLABコードは何ですか? – Eric

答えて

2

つまり、あなたのコード(そして明らかにnp.linalg.lstsq)は、np.linalg.pinvで実装されているMoore-Penrose pseudoinverseを使用しています。 MATLABとMathematicaは、システムを解くためにガウス消去を使用する可能性が高いです。我々は、LU分解を用いて、Pythonでこの後者のアプローチを複製することができ:

B = np.array([[2,4,6],[1,0,3],[0,7,0]]) 
y = np.array([[12],[4],[7]]) 
P, L, U = scipy.linalg.lu(B) 

これはB = P L UとしてBを分解、Uは現在上位対角行列であり、そしてP Lは可逆です。特に、我々は見つける:

>>> U 
array([[ 2., 4., 6.], 
     [ 0., 7., 0.], 
     [ 0., 0., 0.]]) 

>>> np.linalg.inv(P @ L) @ y 
array([[ 12.], 
     [ 7.], 
     [ 0.]]) 

を目的は、この下に決定された、変換された問題、U x = (P L)^{-1} yを解決することを目的とします。解決策セットは元の問題と同じです。解をx = (x_1, x_2, x_3)と書いてみましょう。それから、すぐに解決策にはx_2 = 1が必要であることがわかります。それでは2 x_1 + 4 + 6 x_2 = 12が必要です。 x_1を解くと、x_1 = 4 - 3 x_2となります。それで、どんな解決策も、形式は(4 - 3 x_2, 1, x_2)です。

上記のソリューションを生成する最も簡単な方法は、単にx_2 = 1を選択することです。次にx_1 = 1と入力すると、MATLABから(1、1、1)という答えが返されます。一方

np.linalg.pinvBためpseudionverse特性を満足ユニーク行列であるムーア・ペンローズ疑似逆を計算します。ここでの重点はユニークです。したがって、あなたが言う:

私の質問は、Pythonは無限の解空間から答えを選択する方法である答えはあなたが逆行列を使用するときに選んだのすべてが実際にあなたによって行われていることである

np.linalg.pinv(B)はユニークな行列なので、np.linalg.pinv(B) @ yはユニークです。

解決策の完全なセットを生成するには、上の@ali_mのコメントを参照してください。

+0

ありがとう@ jme - あなたの答えは私がコードの下で何が起こっているのか理解するのを助けました!この問題のLU分解を実装します。 – SMal

関連する問題