numpyの精度：数値の比較中に問題が発生する

Question

最初に背景のビット。具体的には、固有ベクトルを見つけたら、$ argsort $を使って、固有値の1つをソートし、その値に並べ替えを適用するパーミュテーションを探します。元の行列。

ここでは、numpyパッケージを使用してPythonでコードを実装しました。コード自体は再帰的であり、固有ベクトル内の等価な値の集合を見つけると、等しい値を持つ添え字に対応する対称部分行列を抽出し、この行列にアルゴリズムを再び適用します。

これはすべて非常にうまくいっていますが、大部分は大変な仕事ですが、固有ベクトルの等しいエントリに対応する必要のあるインデックスの束が等しい値を持つとは認識されませんでした。問題は、値があるアルゴリズム（おそらくLanczos、しかし私は完全にnumpyに慣れていない）によって機械精度に計算されたということでした。 ）ルーチンseriation（

>>> T=spectral.seriation(A,index) 

    columns [ 0 1 2 3 4 5 6 7 8 9 10 11] 

    [ 3.30289130e-01 -2.75240941e-01 -2.75240941e-01 3.30289130e-01 
    -2.75240941e-01 3.30289130e-01 -2.75240941e-01 3.30289130e-01 
    3.30289130e-01 -2.75240941e-01 -1.69794463e-16 -2.75240941e-01] 

    [ 4 6 9 1 2 11 10 0 5 7 8 3] 

    difference -5.55111512313e-17 

再帰関数である：これは私が明示的に固有ベクトルで2つのエントリの違いをチェックしたサンプル出力、です。浮動小数点数の配列は考慮中の固有ベクトルであり、その下の配列は列のソート順を示します。列[4,6,9,1,2,11]は同じ値を持つことに注意してください。ただし、固有ベクトルと固有値の計算は常に近似値であり、実際には9桁目と2桁のエントリの差を出力するときはゼロではありません。アルゴリズムが[4,6,9,1,2,11]をグループ化する必要がある場合、それは[4,6,9]のみをグループ化し、他のグループに残りを入れてレンチを投射します。

これは質問です。numpyで任意の精度計算を実行する方法はありますか？これに失敗した場合、この問題に対する「良い」回避策は何でしょうか？

また、これらのエントリが等しくなければならないことが数学的に証明できるとも言えるでしょう。それは行列の性質ですが、うまくいけば問題に密接に関連していません。

Answer 1

ダブルは正確には実数ではありません[合理的ではありません]。すべての範囲に無限の数の有理数が存在します（正確には、少なくとも2つの要素を持つすべての範囲）が、それらを表す有限数のビットのみです。
したがって、「正確な」計算には丸め誤差がいくつかあるはずです。より多くのinforamtionのために

、あなたは同等の大きさの2つの浮動小数点数の減算を行う場合には、精度、すなわち場合、問題になることはありませんwhat every computer scientist should know about floating-point arithmetic

Answer 2

を読みたいかもしれません[2]、[9]実際にあります同じ場合、差はゼロになります。

実際には、デフォルトで出力には小数点以下8桁まで表示されますが、数値は異なります。通常、2桁の精度は小数点以下16桁です（実行するにはnumpy.finfo(numpy.float).eps可能な限り最小の数を与える機械イプシロン）

"%.16f\n%.16f" % myarray[[2, 9]]の出力フォーマットを使用して検査してください。

違いはありますが、類似性が7d.pに満足すれば、numpy.around(differences, 7)などの結果を切り捨てることができます。

また、データを前処理したい場合は、次のようなものを使用することもできます（これを行うより良い方法がありますが）。

sigcnd, expn = numpy.frexp(myarray) 
sigcnd = numpy.around(sigcnd, 7) 
truncated_myarray = numpy.ldexp(sigcnd, expn)

Answer 3

あなたが何か行うことができます与えられた許容値とほぼ同等の要素のインデックスをしたい場合：（numpy.allclose()の用途と同一のおおよその比較を使用して）

def almost_matches(x, array, rtol=1e-05, atol=1e-08): 
    answer = [] 
    for y in xrange(len(array)): 
     if abs(x-array[y]) <= (atol + rtol * abs(array[y])): 
      answer.append(y) 
    return answer 

を

>>> a = [3.30289130e-01, -2.75240941e-01, -2.75240941e-01, 3.30289130e-01, -2.75240941e-01, 3.30289130e-01, -2.75240941e-01, 3.30289130e-01, 3.30289130e-01, -2.75240941e-01, -1.69794463e-16, -2.75240941e-01] 
>>> almost_matches(min(a), a) 
[1, 2, 4, 6, 9, 11] 

Answer 4

numpy.allcloseとnumpy.iscloseに、公差内の等価性をテストする関数をチェックします。