JMPには、Inverse()、GInverse()、Sweep()という、逆行列を計算するための関数があります。Solve()関数は、連立一次方程式の式を解くのに使用します。
Inverse()関数は、正則な正方行列の逆行列を戻します。Inverse()は、省略してInvと書くこともできます。行列Aにおいて、AとInverse(A)の積(よくA(A-1)と記載される)は、結果として単位行列を戻します。
A = [5 6, 7 8];
AInv = Inv( A );
A*A Inv;
[1 0,0 1]
A = [1 2, 3 4];
AInv = Inverse( A );
A*AInv;
[1 1.110223025e-16,0 1]
メモ: 2番目の例にあるように、浮動小数点の精度の限界が原因で、値がわずかに異なる場合があります。
行列Aの(Moore-Penrose型の)一般逆行列とは、次のような行列Gです。
AGA=A
GAG=G
(AG)´ = AG
(GA)´ = GA
GInverse()関数は、非正方行列を含めた任意の行列を引数に取り、特異値分解を使ってMoore-Penrose型の一般逆行列を求めます。この関数は、フルランクでない行列の逆行列を求めるときに便利です。次の方程式系を見てください。
この方程式系の解を求めるには、次のスクリプトを使用します。
A = [1 2 2, 2 4 4, 1 1 1];
B = [6, 12, 1];
Show( GInverse( A ) * B );
GInverse(A) * B = [-4, 2.5, 2.5];
Solve()関数は、連立一次方程式の解を求めます。Solve()は、x = A-1bとなるベクトルxを求めます(ここで、Aは正則な正方行列、bはベクトル)。この行列Aとベクトルbの行数は同じである必要があります。Solve(A,b)はInverse(A)*bと同じです。
A = [1 -4 2, 3 3 2, 0 4 -1];
b = [1, 2, 1];
x = Solve( A, b );
[-16.9999999999999, 4.99999999999998, 18.9999999999999]
A*x;
[1, 2, 0.999999999999997]
メモ: 例にあるように、浮動小数点の精度の限界が原因で、値がわずかに異なる場合があります。
Sweep()関数は、正方行列の部分(ピボット)の逆行列を求めます。これをすべてのピボットに対して行うと、最終的には逆行列が求められます。通常、この行列は対角軸が0にならないよう、正定値行列(または負定値行列)である必要があります。Sweep()は、与えられた行列が正定値行列であるかどうかはチェックしません。正定値行列でない場合も、掃き出し法において、対角軸に0以外の値が現れる限り演算は続けられます。0(または0に近い値)が対角軸に現れた場合、その行と列を0にします。結果としてg2型の一般逆行列が求められます。
行列Eが、次に示すように、さらに小さな行列A、B、C、Dに分割されているとします。
Sweep()の構文は次のようになります。
Sweep( E, [...] );
[...]は行列Aの部分を示します。
これにより、次のような結果の行列が作成されます。
メモ:
• Aの位置の部分行列は逆行列になります。
• 部分行列Bの位置はAx = Bの解になります。
• 部分行列Cの位置はxA = Cの解になります。
Sweep()は逐次的で可逆的です。
• A = Sweep( A, {i,j} )はA = Sweep( Sweep( A, i ), j )と同じで、逐次的です。
• A = Sweep( Sweep( A, i ), i )はAを元の値に復元します。これは可逆的です。
次のような交差積の部分から成る行列があるとします。
X′Xを引数としたSweep演算の結果は、次のようになります。
統計的な観点からは、次のように見ることができます。
• 右上部はモデルY = Xb + eの最小2乗推定値
• 右下部は誤差の平方和
• 左上部は推定値の共分散に比例する行列
Sweep関数は、ステップワイズ回帰で必要な部分的な解を求める際に有効です。
2つ目の引数は、掃き出しを行う行(列)を列挙したベクトルです。たとえば、Eが4x4の行列で、4行すべてに掃き出し法を適用して、E-1を得るには、以下のように指定します。
E = [ 5 4 1 1, 4 5 1 1, 1 1 4 2, 1 1 2 4];
Sweep( E, [1, 2, 3, 4] );
[0.56 -0.44 -0.02 -0.02,
-0.44 0.56 -0.02 -0.02,
-0.02 -0.02 0.34 -0.16,
-0.02 -0.02 -0.16 0.34]
Inverse( E ); // これらの結果は同じ
[0.56 -0.44 -0.02 -0.02,
-0.44 0.56 -0.02 -0.02,
-0.02 -0.02 0.34 -0.16,
-0.02 -0.02 -0.16 0.34]
メモ: Sweep()の詳細と、逆行列を作成するGauss-Jordan法との関連については、Goodnight, J.H. (1979) “A tutorial on the SWEEP operator.” The American Statistician, August 1979, Vol. 33, No. 3. pp. 149-58を参照してください。
Sweep()については、分散分析(ANOVA)の例でさらに説明しています。
Det()関数は、正方行列の行列式を戻します。2 x 2行列の行列式は、下に示すように、対角線上の要素の積の差です。n x n行列の行列式は、(n - 1) × (n - 1)の行列式の重み和として再帰的に定義されます。逆行列を作成するには、行列の行列式が0以外でなければなりません。
F = [1 2, 3 5];
Det( F );
-1