約 11,911 件
https://w.atwiki.jp/nopu/pages/91.html
次も参照 Matrix Reference Manual The Matrix Cookbook 行列の見方 0. 環Rの元を矩形状に並べたものを行列という。 同じサイズの行列Mat(m,n)は,和に関してアーベル群である。 また同じサイズの正方行列の全体M(n)は,さらに「行列の積」を入れて非可換環(結合的多元環,代数)になる。 この積はCayleyが初めて使ったとされる。 M(n)は零因子を持ちうる。従って整域でない。 特に,M(n)の部分集合で,逆元を持つ元だけを集めた体GL(n)を正則行列という。 正則行列の重要な部分体として,直交行列(ユニタリ行列)や対称行列(エルミート行列)がある。 前者は内積を変えない線形写像であり,行列式の絶対値が1という特徴付けをもつ。正規直交基底を並べた行列と見ることもできる。 後者は二次形式の議論において重要であり,必ず実固有値を持つという強い性質を持つ。 ここから正定値性やシルベスターの慣性法則,ミニマックス定理およびその系としての単調定理などの性質が導かれる。 直交行列や対称行列はまた,正規行列の重要な例でもある。 正規行列はユニタリ相似変換によって対角化可能である。 1. 線形写像の表現行列 線形空間 線形写像 このとき,の基底を固定すると, に対応する表現行列が得られる。 「表現」というのは,抽象的な線形写像から行列への同型φを考えていることを意味する。 このときまた,共役写像 の表現行列は転置行列である。 固有値,特異値はこの視点からの議論。 A が対角化可能の場合,この手の話では固有ベクトルを基底にとる(固有ベクトル展開)と見通しが良くなる。 実際, に対して, であるから,事態を固有座標系で見れば,は単に各成分を固有値倍するスケーリング作用素(=対角行列)に他ならないことが分かる。 つまり対角化という作業は, という一般の線形写像を, 固有座標系への座標変換と, 成分毎のスカラー倍と, 元の座標系への座標変換の合成写像に分解する作業なのである。 一般の対角化できない行列に対しては,特異値分解を用いて同様の幾何学的意味付けができる。 Rank, Ker, Im などもこの視点。 1.5. 線形準同型の表現行列あるいは座標変換として 線形写像はまた,線形空間の準同型である。 特に,正則つまり全単射ならば同型になる。 特に,抽象的な n-dim 線形空間 V と, とは,次のように結び付けられる。 Vの基底ベクトルを並べた行列 として, これはまた,以下のように複数の視点で捉えることができる。 (i) 座標ベクトルを抽象的な線形空間の元に対応付ける写像 (ii) 座標ベクトルの空間をユークリッド空間と同一視したうえで,そこからVへの座標変換 ←本当? 2. ベクトルを並べたもの 線形写像Aが作用するベクトル空間の元vを並べたもの。 直交行列は,しばしば基底ベクトルを並べたものと見ると幾何学的な理解がしやすい。 また,一次独立の判定法として,行列式やグラム行列を用いる方法はこの観点に立っている。 行列式の幾何学的意味もまた、各ベクトルが張る超立体の体積と考えることができる。 3. 行列自体をベクトル空間の元として見る。 さらに内積やノルムを入れることができる。 i) 一本ベクトルとしてのpノルム p=2 のときフロベニウスノルム(ユークリッドノルム) p=∞ のとき最大値ノルム ii) 誘導ノルム(作用素ノルム) p=1 のとき列毎のl1ノルムの最大値 p=2 のとき(正方行列に限る)スペクトルノルム(=最大特異値) p=∞ のとき行毎のl1ノルムの最大値 任意の誘導ノルムは,スペクトル半径(絶対値最大の固有値)によって下から押さえられる。 4. 行列自体を環の元として見る。 行列の和を加法,行列の積を乗法として,非可換環になる。 ただし整域にはならない。つまり,零因子を持つ。 →零因子の作り方参照 特に三角行列や対角行列のなす部分環は可換環である。 さらに,完備なノルムを入れてBanach環の元と見ることができる。 5. 二次形式や内積の表現として 対称行列はこの観点の性質がよく研究されている。正定値,シルベスターの慣性法則など。 またその拡張として対称作用素の理論もよく研究されている。 6. 群の表現として GL(n)は行列の積を群演算として位相群(Lie群)(群演算が連続写像になること)になっている。 7. 変換あるいは添字付けられた何かあるいはテンソルとして 物理では添字に着目して話をしているように見える人も多い気がする。 要するに線形変換 8. 離散空間上の二変数関数として という書き方は, 有限集合から係数体Kへの写像 を外延的(つまり網羅的)に記述したものだとみなせる。 この考えによれば,行列の行ないし列を無限に飛ばすことは, 関数の定義域を無限集合に拡張することに他ならず,通常の二変数関数になる。 名前のついた行列 Gram行列 内積を成分に持つ行列 実は次のように書ける。 Aのadjoint作用素 Def. cofactors; 余因子 行列Aのi行とj列を潰して得られる行列をとする。 Aの余因子Δijとは,次で定義される行列式である。 各成分を余因子にもつ行列を,余因子行列(adjoint)といい,adj A と書く。
https://w.atwiki.jp/geogebra_kyozai/pages/27.html
行列 行列の加法,減法と実数倍 行列の乗法 行列の積1 行列の乗法の性質 逆行列 連立一次方程式 行列の対角化 点の移動と1次変換 合成変換と逆変換 回転移動と1次変換 回転移動 原点を通る直線に関する対称移動 1次変換と直線 *
https://w.atwiki.jp/center_math/pages/116.html
意かは基本的に正方行列とする。 ベクトル方程式 正方行列の乗数 対角行列 対角行列は単純に各項をn乗する。 を与えるのならば、 三角行列 大変は単純に各項をn乗する。 残りの1項は数列を解くことで得られる を与えた時、 とおくと、 これを特性方程式を経て解くと 一般行列 Δ≠0を満たす行列P, 対角行列D,を用いて、 を満たすような値を与える。 Pの候補としては、などが考えられる。 次に、これを用いて、 A-Et 上のP,Dを用いて D+tEは対角行列なので簡単にその答えを求められる。 ハーリー・ケミルトンの定理 のとき、 のとき、 が成り立つ。 逆行列が存在しないとき、 さらに、と表すと、 が成り立つ。 これらを使うと、特に、 行列の証明 まずはabcdを用いずに考える。 逆行列の証明 ad-bc=0の証明と同じ。 与えられた条件に対して、 が存在する場合を考えて、 矛盾を導く 例) が存在するとすると、 二番目の式に左からを掛けたものに矛盾 また、AB=BA=Eが成り立つ時、 例) のとき、 だから、 2乗=O を満たすとき、 とおくと、k A^2-pA+qE=Oの解 ハミルトン・ケーリーの法則より、 よって、 (1)p=a+dのとき、(q-ad+bc)E=O⇔q=ad-bc よって、p=a+d⇔q=ad-bc (2)p≠a+dのとき、より、 よって、A=αE,βE 逆行列 また、 行列P、Qが逆行列を持つのであれば、行列PQも逆行列を持つ。 (証明1)より、原理的に成り立つ。 (証明2)より、なら、である。 行列における逆行列 逆行列が存在しないとき、 反証 特殊な行列の性質
https://w.atwiki.jp/linearalgebra/pages/76.html
5.行列 5-1行列の演算 5-2転置行列 5-3複素共役行列 5-4随伴行列
https://w.atwiki.jp/opengl/pages/101.html
さて、最初に解説するのは 『単位行列』 という物です。 ですが、単位行列の前にまず正方行列を解説します。 正方行列というのは 2×2行列 とか 3×3行列 とか 4×4行列 のように 行と列の並び数が同じ行列の事を言います。 行と列の数がともに nである行列を n次の正方行列といいます。 DirectX とか OpenGL では 4次の正方行列を使います。 つまり 4×4行列 の事ですね。 そして 『単位行列』 というのは正方行列の左上から右下まで全て 1 で、 それ以外の所は全て 0 の行列の事です。わかりやすくすると次の図のような物です。 2次の正方行列の単位行列 3次の正方行列の単位行列 4次の正方行列の単位行列 そもそも、『単位行列』 って何? 何故、それが必要なの? とお思いでしょう。 この 『単位行列』 というのは glLoadIdentity(); がしている事で、 簡単に言うと初期化の意味があります。 では次回、どうして初期化なのかを説明しましょう。
https://w.atwiki.jp/opengl/pages/139.html
転置行列です。 行列の行と列を入れ替えて作られる行列をの転置行列と言い、と表す。 行列:行列の合成 のページで //合成 glMultMatrixfの代わり void multiplication(GLfloat* src1,GLfloat* src2,GLfloat* dst) { for(int y=0;y 4;y++){ for(int x=0;x 4;x++){ dst[y*4+x]=src2[y*4]*src1[x]+src2[y*4+1]*src1[x+4]+src2[y*4+2]*src1[x+8]+src2[y*4+3]*src1[x+12]; } } } となっていて、 「掛け合わせる順番が逆ジャネーノ?」 と思った方もいるでしょうが、これで良いのです。 DirectXは行ベクトルのため、左から掛け合わせますがOpenGLの場合、 列ベクトルであるために左から掛け合わせるには転置行列を作成して 計算しなければなりません。 しかし、列ベクトルの場合は右から掛け合わせると転置行列を作成しなくても 解を導き出せます。 一種の最適化ですね。
https://w.atwiki.jp/cvlec/pages/14.html
行列表現 行列は、M行xN列と行の次に列として表現する。 そのため、添え字は と表記する。 3x4の行列の要素は、次のような要素となる。 行列は、乗算する順番が重要となる。右側からかけるか、左側からかけるかが重要となる。
https://w.atwiki.jp/opengl/pages/140.html
逆行列です。 逆行列とは、あるN×N行列があったときと積をとると N×N行列の単位行列になる行列の事をの逆行列と呼び、と書く。 つまり、各種アフィン変換を打ち消して単位行列(初期状態)に戻す行列の事です。 4×4行列の逆行列の公式 の時、逆行列が存在して、 で、求められます。 b は以下のようになります。
https://w.atwiki.jp/introtopython/pages/10.html
目次 目次 行列を作る 行列の行数と列数を得る 行列の様々な情報を得る 転置行列を作る 転置行列を求める 行列の掛け算を行う(積を求める) 行列の掛け算 行列式を求める 逆行列を求める QR分解を行う 一般逆行列を求める 連立一次方程式を解く 行列を作る NumPyモジュールのarray関数を使う。以下は3行2列の行列を作成した例。 import numpy as np mx = np.array([[1, 2], [3, 4], [5, 6]]) print(mx) [[1 2] [3 4] [5 6]] mx array([[1, 2], [3, 4], [5, 6]]) mx[0, 0] 1 mx[2, 1] 6 行列の要素を取り出す際に指定するインデックスは、0から始まる0および正値の整数であることに注意。 行列の行数と列数を得る shape変数を使う。以下は、2行3列の行列から行数(2)と列数(3)を得る例。 import numpy as np mx = np.array([[1, 2, 3], [4, 5, 6]]) mx.shape (2, 3) mx.shape[0] 2 mx.shape[1] 3 行列の様々な情報を得る import numpy as np mx = np.array([[1, 2], [3, 4], [5, 6], [7, 8]]) print(mx) [[1 2] [3 4] [5 6] [7 8]] print( 次元 {} .format(mx.ndim)) 次元 2 print( 全要素数 {} .format(mx.size)) 全要素数 8 nrow, ncol = mx.shape print( 行数 {} .format(nrow)) 行数 4 print( 列数 {} .format(ncol)) 列数 2 print(mx[2, 1]) 6 print(mx[3, 1]) 8 転置行列を作る NumPyモジュールのndarrayを使う場合はT属性(「アトリビュート」とも)を使う。以下は、2行3列の行列を転置させた例(つまり3行2列の行列を作る)。 import numpy as np mx1 = np.array([[1, 2, 3], [4, 5, 6]]) mx1[1, 2] 6 np.array(mx1).T array([[1, 4], [2, 5], [3, 6]]) mx2 = np.array(mx1).T mx2 array([[1, 4], [2, 5], [3, 6]]) mx2[2, 1] 6 転置行列を求める Tプロパティを使う。 import numpy as np mx = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]]) print(mx) [[1 2 3] [4 5 6] [7 8 9]] print(mx.T) [[1 4 7] [2 5 8] [3 6 9]] mx = np.array([[1, 2, 3], [4, 5, 6]]) print(mx) [[1 2 3] [4 5 6]] print(mx.T) [[1 4] [2 5] [3 6]] 行列の掛け算を行う(積を求める) NumPyモジュールのndarrayを使用した例。 import numpy as np a1 = np.array([[1, 2], [3, 4]]) a2 = np.array([[5, 6, 7], [8, 9, 10]]) print(a1) [[1 2] [3 4]] print(a2) [[ 5 6 7] [ 8 9 10]] np.dot(a1, a2) array([[21, 24, 27], [47, 54, 61]]) 行列の掛け算 dotメソッドを使う。最後の例のとおり、掛けられる行列の列数と掛ける行列の行数が同じでなければ計算することはできない。 mx1 = np.array([[-1, 0, 3], [8, 1, -5]]) mx2 = np.array([[3, 4, 0], [-2, 1, 2]]) mx3 = np.array([[4, -3, 7], [2, 0, -1], [1, 5, 0]]) print(mx1) [[-1 0 3] [ 8 1 -5]] print(mx2) [[ 3 4 0] [-2 1 2]] print(mx3) [[ 4 -3 7] [ 2 0 -1] [ 1 5 0]] np.dot(mx1, mx3) array([[ -1, 18, -7], [ 29, -49, 55]]) np.dot(mx2, mx3) array([[ 20, -9, 17], [ -4, 16, -15]]) np.dot(mx1, mx2) Traceback (most recent call last) File " stdin ", line 1, in module File " __array_function__ internals ", line 5, in dot ValueError shapes (2,3) and (2,3) not aligned 3 (dim 1) != 2 (dim 0) 行列式を求める linalg.detメソッドを使う。 mx = np.array([[3]]) mx array([[3]]) np.linalg.det(mx) 3.0000000000000004 mx = np.array([[1, 2], [3, 4]]) mx array([[1, 2], [3, 4]]) np.linalg.det(mx) -2.0000000000000004 mx = np.array([[0, 3, 2, -1], [-4, 6, 1, 5], [0, -2, 3, -1], [-1, 0, 0, 2]]) mx array([[ 0, 3, 2, -1], [-4, 6, 1, 5], [ 0, -2, 3, -1], [-1, 0, 0, 2]]) np.linalg.det(mx) 27.999999999999986 逆行列を求める linalg.invメソッドを使う。 mx = np.array([[3]]) mx array([[3]]) np.linalg.inv(mx) array([[0.33333333]]) mx = np.array([[1, 2], [3, 4]]) mx array([[1, 2], [3, 4]]) np.linalg.inv(mx) array([[-2. , 1. ], [ 1.5, -0.5]]) mx = np.array([[3, -3, 1], [3, 2, 0], [-1, -5, 1]]) mx array([[ 3, -3, 1], [ 3, 2, 0], [-1, -5, 1]]) np.linalg.inv(mx) array([[ 1. , -1. , -1. ], [-1.5, 2. , 1.5], [-6.5, 9. , 7.5]]) 正則ではない行列(行列式の値が0)の場合は、linalg.invメソッドはエラーを返す。 mx = np.array([[1, 4, 7], [2, 5, 8], [3, 6, 9]]) mx array([[1, 4, 7], [2, 5, 8], [3, 6, 9]]) np.linalg.det(mx) 0.0 np.linalg.inv(mx) Traceback (most recent call last) (表示省略) numpy.linalg.LinAlgError Singular matrix QR分解を行う 行列AをQR分解で行列Qと行列Rに分解することを考える。 A = Q R 行列QはQTQ=I(Iは単位行列)となる直交行列。行列Rは上三角行列(対角成分と対角成分より上の要素はすべて0以外で、対角成分より下の要素はすべて0)。 linalg.qrを使うと、それぞれQR分解した際の行列Qと行列Rが得られる。 import numpy as np import scipy.linalg as linalg a = np.array([[1, 2, 3], [4, 5, 6], [9, 9, 8]]) print(a) [[1 2 3] [4 5 6] [9 9 8]] q, r = linalg.qr(a) print(q) [[-0.10101525 -0.71840915 -0.6882472 ] [-0.40406102 -0.60253671 0.6882472 ] [-0.90913729 0.34761733 -0.22941573]] print(r) [[ -9.89949494 -10.40457121 -10.00051019] [ 0. -1.32094586 -2.98950905] [ 0. 0. 0.22941573]] QTQ = I となるか試す。 print(np.dot(q.T, q)) [[ 1.00000000e+00 -1.50980479e-16 -9.33624473e-17] [-1.50980479e-16 1.00000000e+00 1.17405362e-16] [-9.33624473e-17 1.17405362e-16 1.00000000e+00]] 一般逆行列を求める linalg.pinvメソッドを使う。 import numpy as np mx = np.array([[1, 2, 3, 4], [5, 6, 7, 8], [9, 10, 11, 12], [13, 14, 15, 16]]) print(mx) [[ 1 2 3 4] [ 5 6 7 8] [ 9 10 11 12] [13 14 15 16]] np.linalg.det(mx) -1.820448242817726e-31 np.linalg.inv(mx) array([[ 1.50119988e+15, -3.75299969e+14, -3.75299969e+15, 2.62709978e+15], [-1.95155984e+16, 1.95155984e+16, 1.95155984e+16, -1.95155984e+16], [ 3.45275971e+16, -3.79052969e+16, -2.77721977e+16, 3.11498974e+16], [-1.65131986e+16, 1.87649984e+16, 1.20095990e+16, -1.42613988e+16]]) np.linalg.pinv(mx) array([[-0.285 , -0.145 , -0.005 , 0.135 ], [-0.1075, -0.0525, 0.0025, 0.0575], [ 0.07 , 0.04 , 0.01 , -0.02 ], [ 0.2475, 0.1325, 0.0175, -0.0975]]) 一般逆行列の定義を満たしているか否か確認。 mxi = np.linalg.pinv(mx) np.dot(np.dot(mx, mxi), mx) array([[ 1., 2., 3., 4.], [ 5., 6., 7., 8.], [ 9., 10., 11., 12.], [13., 14., 15., 16.]]) 当然、正方行列でなくとも求めることができる。 mx = np.array([[1, 2, 3, 4], [5, 6, 7, 8]]) print(mx) [[1 2 3 4] [5 6 7 8]] np.linalg.det(mx) (表示省略) numpy.linalg.LinAlgError Last 2 dimensions of the array must be square np.linalg.inv(mx) (表示省略) numpy.linalg.LinAlgError Last 2 dimensions of the array must be square np.linalg.pinv(mx) array([[-5.50000000e-01, 2.50000000e-01], [-2.25000000e-01, 1.25000000e-01], [ 1.00000000e-01, -1.38777878e-17], [ 4.25000000e-01, -1.25000000e-01]]) mxi = np.linalg.pinv(mx) np.dot(np.dot(mx, mxi), mx) array([[1., 2., 3., 4.], [5., 6., 7., 8.]]) linalg.pinvメソッドはムーア・ペンローズ一般逆行列を数値的に求めるメソッドのため、求まった一般逆行列は近似値であることに注意。ムーア・ペンローズ一般逆行列の定義の4つの式について、比較をした例は以下のとおり。すべての要素がTrueにならないことがわかる。 np.dot(np.dot(mx, mxi), mx) == mx array([[False, False, False, True], [False, False, False, False]]) np.dot(np.dot(mxi, mx), mxi) == mxi array([[False, False], [False, True], [False, False], [False, False]]) np.dot(mx, mxi).T == np.dot(mx, mxi) array([[ True, False], [False, True]]) np.dot(mxi, mx).T == np.dot(mxi, mx) array([[ True, False, False, False], [False, True, False, True], [False, False, True, True], [False, True, True, True]]) 連立一次方程式を解く solveメソッドを使う。以下の連立一次方程式を解く(解はx=3,y=2)。 x + 2y = 7 3x - 4y = 1 import numpy as np mxaa = np.array([[1, 2], [3, -4]]) mxy = np.array([7, 1]) mxaa array([[ 1, 2], [ 3, -4]]) mxy array([7, 1]) np.linalg.solve(mxaa, mxy) array([3., 2.]) 同じく以下の連立一次方程式を解く(解はx=2,y=-3,z=-14)。 3x - 3y + z = 1 3x + 2y = 0 -1x - 5y + z = -1 mxaa = np.array([[3, -3, 1], [3, 2, 0], [-1, -5, 1]]) mxy = np.array([1, 0, -1]) mxaa array([[ 3, -3, 1], [ 3, 2, 0], [-1, -5, 1]]) mxy array([ 1, 0, -1]) np.linalg.solve(mxaa, mxy) array([ 2., -3., -14.]) 名前 コメント
https://w.atwiki.jp/r-intro/pages/18.html
目次 目次 基本操作 行列を作る 行列を作る 行列から成分を取り出す 行列から特定の行か特定の列のみ取り出す 行列をベクトルに変換する 行列から非対角成分だけを取り出す 行列に行を追加する 行列に列を追加する 線形代数 転置行列を求める 単位行列を作る 行列のトレース(固有和)を求める 行列式を求める 行列式を求める 行列の掛け算 逆行列を求める 一般逆行列を求める 連立一次方程式を解く 計算 QR分解をする コレスキー分解する 固有値と固有ベクトルを求める 固有値と固有ベクトルを求める 特異値分解を行う 余因子を求める 余因子行列を求める 基本操作 行列を作る matrix関数を使う。行列の要素を格納したベクトルを第一引数に与えて、nrowオプションで行数を、ncolオプションに列数を与える。デフォルトでは列方向に順番に要素が格納される。これを行方向に順に格納するにはbyrowオプションにTRUEを指定する。 d - 1 12 mx - matrix(d, nrow = 3, ncol = 4) print(mx) [,1] [,2] [,3] [,4] [1,] 1 4 7 10 [2,] 2 5 8 11 [3,] 3 6 9 12 mx - matrix(d, ncol = 4) print(mx) [,1] [,2] [,3] [,4] [1,] 1 4 7 10 [2,] 2 5 8 11 [3,] 3 6 9 12 mx - matrix(d, ncol = 5) 警告メッセージ matrix(d, ncol = 5) で データ長 [12] が列数 [5] を整数で割った、もしくは掛けた値ではありません mx - matrix(d, ncol = 4, byrow = TRUE) print(mx) [,1] [,2] [,3] [,4] [1,] 1 2 3 4 [2,] 5 6 7 8 [3,] 9 10 11 12 mx - matrix(1 4, nrow = 2, ncol = 4, byrow = TRUE) print(mx) [,1] [,2] [,3] [,4] [1,] 1 2 3 4 [2,] 1 2 3 4 mx - matrix(0., 2, 5) print(mx) [,1] [,2] [,3] [,4] [,5] [1,] 0 0 0 0 0 [2,] 0 0 0 0 0 上記の例のとおりに、特に何も指定しなければ第二引数が行列、第三引数が列数になる。ベクトルの要素数で行列どちらかが決まればもう片方も決まるため、行または列の省略は可能。ただし、ベクトルの要素数に応じた正しい値でなければならない。 最後の例のとおりに、指定した行数×列数の要素数に対してベクトルの要素数が足らないときは、そのベクトルの要素数が行数×列数の約数であれば、自動的に繰り返し使われる。 行列を作る 行列はmatrix関数を使用して、ベクトルから作ることができる。以下は、2行3列の行列を新規に作成している。 d - 1 6 mx - matrix(d, 2, 3) mx [,1] [,2] [,3] [1,] 1 3 5 [2,] 2 4 6 matrix関数の第一引数には、成分(並べられている各々の数)に代入する値が含まれているベクトルを指定し、特に指定をしなければ第二引数には行数を、第三引数には列数を指定する。行数はnrowオプション、列数はncolオプションで明示的に指定することができ(順不同)、例えば、以下のようにしても同じ。 mx - matrix(d, ncol = 3, nrow = 2) mx [,1] [,2] [,3] [1,] 1 3 5 [2,] 2 4 6 特にオプションを指定しない場合、第一引数に与えたベクトルの要素について、列方向(縦方向)を優先に値を埋めるように代入する。これは見た目の直感となかなか合わない。列方向ではなく、行方向(横方向)を優先に値を埋めるように代入させるには、byrowオプションにTRUEを与える。 mx - matrix(d, 2, 3, byrow = TRUE) mx [,1] [,2] [,3] [1,] 1 2 3 [2,] 4 5 6 行列から成分を取り出す 行ごと列ごとに角括弧(ブラケット)をつかってインデックス(1~)を指定する。インデックスは1から始まることに注意。 mx - matrix(1 6, 3, 2) print(mx) [,1] [,2] [1,] 1 4 [2,] 2 5 [3,] 3 6 mx[1, 2] [1] 4 mx[3, 2] [1] 6 mx[4, 2] mx[4, 2] でエラー 添え字が許される範囲外です 最後の例のように、行列の各長さより大きい値を指定すると、エラーが発生する。 行列から特定の行か特定の列のみ取り出す インデックスを指定する際に、数値を与えずに空欄のままにする。以下の例ではそれぞれ2行目のみ、3列目のみを取り出している。 mx - matrix(1 9, 3, 3) print(mx) [,1] [,2] [,3] [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9 mx[2, ] [1] 2 5 8 mx[, 3] [1] 7 8 9 is.vector(mx[, 3]) [1] TRUE 最後のとおりに、1行(列)だけ取り出したものはベクトルになる。これを行列のままにする場合はdropオプションにFALSEを与える。 mx[2, , drop = FALSE] [,1] [,2] [,3] [1,] 2 5 8 複数行(列)取り出したい場合は、インデックスをベクトルで与える。以下は2行目と3行目のみを取り出した例。 mx[c(2, 3), ] [,1] [,2] [,3] [1,] 2 5 8 [2,] 3 6 9 行列をベクトルに変換する as.vector関数を使用する。 mx1 - matrix(6, nrow = 1, ncol = 1) mx2 - matrix(1 9, nrow = 3, ncol = 3) mx1 [,1] [1,] 6 as.vector(mx1) [1] 6 mx2 [,1] [,2] [,3] [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9 as.vector(mx2) [1] 1 2 3 4 5 6 7 8 9 as.vector(t(mx2)) [1] 1 4 7 2 5 8 3 6 9 1列目の要素すべて、2列目の要素すべて、・・・、という順にベクトルの要素になる。これを1行目の要素すべて、2行目の要素すべて、・・・、としたい場合は、最後の例のとおりにt関数を使って元の行列の転置行列を指定する。 行列から非対角成分だけを取り出す row関数とcol関数を使うと、ベクトル形式で簡単に取り出すことができる。 mx - matrix(1 16, 4, 4) mx [,1] [,2] [,3] [,4] [1,] 1 5 9 13 [2,] 2 6 10 14 [3,] 3 7 11 15 [4,] 4 8 12 16 mx[row(mx) != col(mx)] [1] 2 3 4 5 7 8 9 10 12 13 14 15 応用として、上の例で2、7、12(対角成分のそれぞれ1行下の成分)を取り除くには、以下のようにする。 mx[row(mx) != col(mx) + 1] [1] 1 3 4 5 6 8 9 10 11 13 14 15 16 行列に行を追加する rbind関数を使う。最後の例のとおり、列数が異なる場合はエラーが発生する。 mx1 - matrix(1 9, ncol = 3) mx2 - matrix(10 15, ncol = 3) print(mx1) [,1] [,2] [,3] [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9 print(mx2) [,1] [,2] [,3] [1,] 10 12 14 [2,] 11 13 15 mx3 - rbind(mx1, mx2) print(mx3) [,1] [,2] [,3] [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9 [4,] 10 12 14 [5,] 11 13 15 mx4 - matrix(16 19, ncol = 4) print(mx4) [,1] [,2] [,3] [,4] [1,] 16 17 18 19 rbind(mx1, mx4) rbind(mx1, mx4) でエラー 行列の列数は一致していなければなりません (2 番目の引数を参照) 行列に列を追加する cbind関数を使う。最後の例のとおり、行数が異なる場合はエラーが発生する。 mx1 - matrix(1 9, nrow = 3) mx2 - matrix(10 15, nrow = 3) print(mx1) [,1] [,2] [,3] [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9 print(mx2) [,1] [,2] [1,] 10 13 [2,] 11 14 [3,] 12 15 mx3 - cbind(mx1, mx2) print(mx3) [,1] [,2] [,3] [,4] [,5] [1,] 1 4 7 10 13 [2,] 2 5 8 11 14 [3,] 3 6 9 12 15 mx4 - matrix(16 19, nrow = 4) print(mx4) [,1] [1,] 16 [2,] 17 [3,] 18 [4,] 19 cbind(mx1, mx4) cbind(mx1, mx4) でエラー 行列の行数は一致していなければなりません (2 番目の引数を参照) 線形代数 転置行列を求める t関数を使う。 mx - matrix(1 9, 3, 3, byrow = TRUE) print(mx) [,1] [,2] [,3] [1,] 1 2 3 [2,] 4 5 6 [3,] 7 8 9 t(mx) [,1] [,2] [,3] [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9 mx - matrix(1 6, 2, 3, byrow = TRUE) print(mx) [,1] [,2] [,3] [1,] 1 2 3 [2,] 4 5 6 t(mx) [,1] [,2] [1,] 1 4 [2,] 2 5 [3,] 3 6 単位行列を作る diag関数を使う。引数を一つだけ指定した場合は、n次の単位行列を作成する。 diag(3) [,1] [,2] [,3] [1,] 1 0 0 [2,] 0 1 0 [3,] 0 0 1 diag(4) [,1] [,2] [,3] [,4] [1,] 1 0 0 0 [2,] 0 1 0 0 [3,] 0 0 1 0 [4,] 0 0 0 1 行列のトレース(固有和)を求める Rには行列Aのトレース(固有和)Tr(A)を求めるための関数が標準では組み込まれていない。diag関数とsum関数を組み合わせることで、計算を実現することができる。以下の例では、Tr(A) = 1+5+9=15である。 mx - matrix(1 9, nrow = 3, byrow = TRUE) mx [,1] [,2] [,3] [1,] 1 2 3 [2,] 4 5 6 [3,] 7 8 9 tr - sum(diag(mx)) tr [1] 15 行列式を求める det関数を使う。 d - c(1, 2, 3, 4) mx - matrix(d, 2, 2, byrow = TRUE) mx [,1] [,2] [1,] 1 2 [2,] 3 4 det(mx) [1] -2 行列式を求める det関数を使う。 mx - matrix(3, 1, 1) print(mx) [,1] [1,] 3 det(mx) [1] 3 mx - matrix(1 4, 2, 2, byrow = TRUE) print(mx) [,1] [,2] [1,] 1 2 [2,] 3 4 det(mx) [1] -2 mx - matrix(c(0, -4, 0, -1, 3, 6, -2, 0, 2, 1, 3, 0, -1, 5, -1, 2), 4, 4) print(mx) [,1] [,2] [,3] [,4] [1,] 0 3 2 -1 [2,] -4 6 1 5 [3,] 0 -2 3 -1 [4,] -1 0 0 2 det(mx) [1] 28 行列の掛け算 %*% 演算子を使う。最後の例のとおり、掛けられる行列の列数と掛ける行列の行数が同じでなければ計算することはできない。 mx1 - matrix(c(-1, 8, 0, 1, 3, -5), 2, 3) mx2 - matrix(c(3, -2, 4, 1, 0, 2), 2, 3) mx3 - matrix(c(4, 2, 1, -3, 0, 5, 7, -1, 0), 3, 3) print(mx1) [,1] [,2] [,3] [1,] -1 0 3 [2,] 8 1 -5 print(mx2) [,1] [,2] [,3] [1,] 3 4 0 [2,] -2 1 2 print(mx3) [,1] [,2] [,3] [1,] 4 -3 7 [2,] 2 0 -1 [3,] 1 5 0 mx1 %*% mx3 [,1] [,2] [,3] [1,] -1 18 -7 [2,] 29 -49 55 mx2 %*% mx3 [,1] [,2] [,3] [1,] 20 -9 17 [2,] -4 16 -15 mx1 %*% mx2 mx1 %*% mx2 でエラー 適切な引数ではありません 逆行列を求める solve関数を使う。 mx - matrix(3, 1, 1) print(mx) [,1] [1,] 3 solve(mx) [,1] [1,] 0.3333333 mx - matrix(1 4, 2, 2, byrow = TRUE) print(mx) [,1] [,2] [1,] 1 2 [2,] 3 4 solve(mx) [,1] [,2] [1,] -2.0 1.0 [2,] 1.5 -0.5 mx - matrix(c(3, 3, -1, -3, 2, -5, 1, 0, 1), 3, 3) print(mx) [,1] [,2] [,3] [1,] 3 -3 1 [2,] 3 2 0 [3,] -1 -5 1 solve(mx) [,1] [,2] [,3] [1,] 1.0 -1 -1.0 [2,] -1.5 2 1.5 [3,] -6.5 9 7.5 正則ではない行列(行列式の値が0)の場合は、solve関数はエラーを返す。 mx - matrix(1 9, 3, 3) print(mx) [,1] [,2] [,3] [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9 det(mx) [1] 0 solve(mx) solve.default(mx) でエラー Lapack routine dgesv システムは正確に特異です U[3,3] = 0 一般逆行列を求める MASSパッケージのginv関数を使う。 library(MASS) mx - matrix(1 16, 4, 4, byrow = TRUE) print(mx) [,1] [,2] [,3] [,4] [1,] 1 2 3 4 [2,] 5 6 7 8 [3,] 9 10 11 12 [4,] 13 14 15 16 det(mx) [1] 4.733165e-30 solve(mx) solve.default(mx) でエラー システムは数値的に特異です 条件数の逆数 = 2.95756e-18 ginv(mx) [,1] [,2] [,3] [,4] [1,] -0.2850 -0.1450 -0.0050 0.1350 [2,] -0.1075 -0.0525 0.0025 0.0575 [3,] 0.0700 0.0400 0.0100 -0.0200 [4,] 0.2475 0.1325 0.0175 -0.0975 一般逆行列の定義を満たしているか否か確認。 mxi - ginv(mx) mx %*% mxi %*% mx [,1] [,2] [,3] [,4] [1,] 1 2 3 4 [2,] 5 6 7 8 [3,] 9 10 11 12 [4,] 13 14 15 16 当然、正方行列でなくとも求めることができる。 mx - matrix(1 8, 2, 4, byrow = TRUE) print(mx) [,1] [,2] [,3] [,4] [1,] 1 2 3 4 [2,] 5 6 7 8 det(mx) determinant.matrix(x, logarithm = TRUE, ...) でエラー x は正方行列でなければなりません solve(mx) solve.default(mx) でエラー a (2 x 4) は正方行列である必要があります ginv(mx) [,1] [,2] [1,] -0.550 2.500000e-01 [2,] -0.225 1.250000e-01 [3,] 0.100 -2.081668e-17 [4,] 0.425 -1.250000e-01 mxi - ginv(mx) mx %*% mxi %*% mx [,1] [,2] [,3] [,4] [1,] 1 2 3 4 [2,] 5 6 7 8 ginv関数はムーア・ペンローズ一般逆行列を数値的に求める関数のため、求まった一般逆行列は近似値であることに注意。ムーア・ペンローズ一般逆行列の定義の4つの式について、比較をした例は以下のとおり。すべての要素がTRUEにならないことがわかる。 mx %*% mxi %*% mx == mx [,1] [,2] [,3] [,4] [1,] FALSE FALSE TRUE TRUE [2,] FALSE FALSE FALSE TRUE mxi %*% mx %*% mxi == mxi [,1] [,2] [1,] FALSE FALSE [2,] FALSE FALSE [3,] FALSE FALSE [4,] FALSE FALSE t(mx %*% mxi) == mx %*% mxi [,1] [,2] [1,] TRUE FALSE [2,] FALSE TRUE t(mxi %*% mx) == mxi %*% mx [,1] [,2] [,3] [,4] [1,] TRUE FALSE FALSE FALSE [2,] FALSE TRUE FALSE FALSE [3,] FALSE FALSE TRUE FALSE [4,] FALSE FALSE FALSE TRUE 連立一次方程式を解く solve関数を使う。以下の連立一次方程式を解く(解はx=3,y=2)。 x + 2y = 7 3x - 4y = 1 mxaa - matrix(c(1, 2, 3, -4), 2, 2, byrow = TRUE) mxy - matrix(c(7, 1), 2, 1) print(mxaa) [,1] [,2] [1,] 1 2 [2,] 3 -4 print(mxy) [,1] [1,] 7 [2,] 1 solve(mxaa, mxy) [,1] [1,] 3 [2,] 2 同じく以下の連立一次方程式を解く(解はx=2,y=-3,z=-14)。 3x - 3y + z = 1 3x + 2y = 0 -1x - 5y + z = -1 mxaa - matrix(c(3, 3, -1, -3, 2, -5, 1, 0, 1), 3, 3) mxy - matrix(c(1, 0, -1), 3, 1) print(mxaa) [,1] [,2] [,3] [1,] 3 -3 1 [2,] 3 2 0 [3,] -1 -5 1 print(mxy) [,1] [1,] 1 [2,] 0 [3,] -1 solve(mxaa, mxy) [,1] [1,] 2 [2,] -3 [3,] -14 計算 QR分解をする 行列XをQR分解で行列Qと行列Rに分解することを考える。 X = Q R 行列QはQ^T Q=I(Iは単位行列)となる直交行列。行列Rは上三角行列(対角成分はすべて0以外で、対角成分より下の要素はすべて0)。 qr.Q関数、qr.R関数を使うと、それぞれQR分解した際の行列Qと行列Rが得られる。引数にはqrオブジェクトを与える必要があり、qr関数を使ってあらかじめ分解したい行列Xからqrオブジェクトを作成しておく必要がある。 d - c(1, 2, 3, 4, 5, 6, 9, 9, 8) x - matrix(d, 3, 3, byrow = TRUE) x [,1] [,2] [,3] [1,] 1 2 3 [2,] 4 5 6 [3,] 9 9 8 objqr - qr(x) q - qr.Q(objqr) q [,1] [,2] [,3] [1,] -0.1010153 -0.7184092 -0.6882472 [2,] -0.4040610 -0.6025367 0.6882472 [3,] -0.9091373 0.3476173 -0.2294157 r - qr.R(objqr) r [,1] [,2] [,3] [1,] -9.899495 -10.404571 -10.0005102 [2,] 0.000000 -1.320946 -2.9895090 [3,] 0.000000 0.000000 0.2294157 Q^T Q = I となるか試す。 t(q) %*% q [,1] [,2] [,3] [1,] 1.000000e+00 -1.110223e-16 -1.665335e-16 [2,] -1.110223e-16 1.000000e+00 -9.714451e-17 [3,] -1.665335e-16 -9.714451e-17 1.000000e+00 qr.X関数を使うことで、qrオブジェクトの元の行列を取り出すことができる。 qr.X(objqr) [,1] [,2] [,3] [1,] 1 2 3 [2,] 4 5 6 [3,] 9 9 8 コレスキー分解する chol関数を使う。 mx0 - matrix(c(1, 0, 1, 0, 2, 0, 1, 0, 3), 3, byrow = TRUE) mx0 [,1] [,2] [,3] [1,] 1 0 1 [2,] 0 2 0 [3,] 1 0 3 mx - chol(mx0) mx [,1] [,2] [,3] [1,] 1 0.000000 1.000000 [2,] 0 1.414214 0.000000 [3,] 0 0.000000 1.414214 得られた結果が正しいか、確認する。 t(mx) %*% mx [,1] [,2] [,3] [1,] 1 0 1 [2,] 0 2 0 [3,] 1 0 3 上記と同じ計算(A T A)はcrossprod関数を使うと同じ結果を得ることができる。 crossprod(mx) [,1] [,2] [,3] [1,] 1 0 1 [2,] 0 2 0 [3,] 1 0 3 固有値と固有ベクトルを求める eigen関数を使う。 # 例1 d - c(1, 3, -2, -4) mx - matrix(d, 2, 2) print(mx) [,1] [,2] [1,] 1 -2 [2,] 3 -4 eig - eigen(mx) print(eig) eigen() decomposition $values [1] -2 -1 $vectors [,1] [,2] [1,] 0.5547002 0.7071068 [2,] 0.8320503 0.7071068 # 固有値(合計2つ) print(eig$value) [1] -2 -1 # 固有値 -2 のときの固有ベクトル print(eig$vector[, 1]) [1] 0.5547002 0.8320503 # 固有値 -1 のときの固有ベクトル print(eig$vector[, 2]) [1] 0.7071068 0.7071068 # 例2 d - c(1, 1, 1, 0, 1, 0, 2, 1, 0) mx - matrix(d, 3, 3) print(mx) [,1] [,2] [,3] [1,] 1 0 2 [2,] 1 1 1 [3,] 1 0 0 eig - eigen(mx) print(eig) eigen() decomposition $values [1] 2 1 -1 $vectors [,1] [,2] [,3] [1,] 0.5345225 0 -0.7071068 [2,] 0.8017837 1 0.0000000 [3,] 0.2672612 0 0.7071068 # 固有値(合計3つ) print(eig$value) [1] 2 1 -1 # 固有値 2 のときの固有ベクトル print(eig$vector[, 1]) [1] 0.5345225 0.8017837 0.2672612 # 固有値 1 のときの固有ベクトル print(eig$vector[, 2]) [1] 0 1 0 # 固有値 -1 のときの固有ベクトル print(eig$vector[, 3]) [1] -0.7071068 0.0000000 0.7071068 固有値は数値的に求めているため、基本的に戻り値は近似値である。そのため、解析的には整数で得られる場合でも戻り値が複素数になる場合がある。この場合、固有ベクトルも複素数になる。以下の例では、解析的には固有値は3と-1(二重根)と求められるが、戻り値は複素数になっている。虚部が無視できるくらい小さいため、虚部を切り捨てて実数ととして扱ってかまわないが、ケースバイケースであることに注意。 d - c(1, 2, -1, 2, 1, -1, -4, 4, -1) mx - matrix(d, 3, 3) eig - eigen(mx) print(eig) eigen() decomposition $values [1] 3+0i -1+0i -1-0i $vectors [,1] [,2] [,3] [1,] -0.9045340+0i -7.071068e-01+0.000000e+00i -7.071068e-01+0.000000e+00i [2,] -0.3015113+0i 7.071068e-01-0.000000e+00i 7.071068e-01+0.000000e+00i [3,] 0.3015113+0i 0.000000e+00+2.637113e-09i 0.000000e+00-2.637113e-09i 固有値と固有ベクトルを求める eigen関数を使う。戻り値はリストで、順番にvaluesに固有値が、vectorsにその固有値のときの固有ベクトルが含まれる。 mxaa - matrix(c(3, 1, 1, 1, 2, 0, 1, 0, 2), nrow = 3) print(mxaa) [,1] [,2] [,3] [1,] 3 1 1 [2,] 1 2 0 [3,] 1 0 2 lam - eigen(mxaa) print(lam) eigen() decomposition $values [1] 4 2 1 $vectors [,1] [,2] [,3] [1,] 0.8164966 0.0000000 0.5773503 [2,] 0.4082483 -0.7071068 -0.5773503 [3,] 0.4082483 0.7071068 -0.5773503 この例では固有値と固有ベクトルは3つあり、例えば3番目の固有値と固有ベクトルを取り出すには、以下のようにする。 lam$values[3] [1] 1 lam$vectors[, 3] [1] 0.5773503 -0.5773503 -0.5773503 特異値分解を行う svd関数を使う。特異値分解の等式は以下のとおり。 A = U D V^T UとVは直交行列。DはAの特異値を対角成分にもつ対角行列。 a - matrix(1 9, 3, 3) a [,1] [,2] [,3] [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9 r - svd(a) r $d [1] 1.684810e+01 1.068370e+00 5.543107e-16 $u [,1] [,2] [,3] [1,] -0.4796712 0.77669099 0.4082483 [2,] -0.5723678 0.07568647 -0.8164966 [3,] -0.6650644 -0.62531805 0.4082483 $v [,1] [,2] [,3] [1,] -0.2148372 -0.8872307 0.4082483 [2,] -0.5205874 -0.2496440 -0.8164966 [3,] -0.8263375 0.3879428 0.4082483 svd関数の戻り値はリスト。要素dに特異値が代入されている。 r$d [1] 1.684810e+01 1.068370e+00 5.543107e-16 上記の恒等式における、行列U(戻り値では要素u)と行列V(戻り値では要素v)が直交行列であることを確認する。 t(r$u) %*% r$u [,1] [,2] [,3] [1,] 1.000000e+00 1.110223e-16 5.551115e-17 [2,] 1.110223e-16 1.000000e+00 -1.110223e-16 [3,] 5.551115e-17 -1.110223e-16 1.000000e+00 t(r$v) %*% r$v [,1] [,2] [,3] [1,] 1.000000e+00 0.000000e+00 -1.110223e-16 [2,] 0.000000e+00 1.000000e+00 1.110223e-16 [3,] -1.110223e-16 1.110223e-16 1.000000e+00 上記の恒等式が成り立っているか確認する。 r$u %*% diag(r$d) %*% t(r$v) [,1] [,2] [,3] [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9 行列Aの特異値は、A^T Aの正の固有値の正の平方根でも求められる。 sqrt(eigen(t(a) %*% a)$values) [1] 1.684810e+01 1.068370e+00 7.942757e-08 3つ目は実質0ということだが、高精度に計算されていないことがわかる。 余因子を求める polyMatrixパッケージのcofactor関数を使う。以下、使用例。計算に用いた行列はここを参考にしており、先に余因子行列をadjoint関数を使用して求めている。 d - c(8, 3, 4, 1, 5, 9, 6, 7, 2) mx - matrix(d, 3) print(mx) [,1] [,2] [,3] [1,] 8 1 6 [2,] 3 5 7 [3,] 4 9 2 adjoint(mx) [,1] [,2] [,3] [1,] -53 52 -23 [2,] 22 -8 -38 [3,] 7 -68 37 cofactor(mx, 1, 2) [1] 22 cofactor(mx, 2, 3) [1] -68 余因子行列を求める polyMatrixパッケージのadjoint関数を使う。以下、使用例。計算に用いた行列はここhttps //jp.mathworks.com/help/symbolic/adjoint.htmlを参考にしている。 d - c(8, 3, 4, 1, 5, 9, 6, 7, 2) mx - matrix(d, 3) mx [,1] [,2] [,3] [1,] 8 1 6 [2,] 3 5 7 [3,] 4 9 2 adjoint(mx) [,1] [,2] [,3] [1,] -53 52 -23 [2,] 22 -8 -38 [3,] 7 -68 37 余因子行列の性質 A A~ = A~ A ~ = |A|・En を確認する。 print(mx %*% adjoint(mx)) [,1] [,2] [,3] [1,] -3.600000e+02 2.842171e-13 -5.684342e-14 [2,] -7.815970e-14 -3.600000e+02 -5.684342e-14 [3,] -8.704149e-14 1.421085e-13 -3.600000e+02 print(adjoint(mx) %*% mx) [,1] [,2] [,3] [1,] -3.600000e+02 2.842171e-14 -1.136868e-13 [2,] 5.684342e-14 -3.600000e+02 5.684342e-14 [3,] 8.526513e-14 1.705303e-13 -3.600000e+02 print(det(mx) * diag(nrow(mx))) [,1] [,2] [,3] [1,] -360 0 0 [2,] 0 -360 0 [3,] 0 0 -360 上2つの計算結果の非対角成分は非常に小さいので0と見なしてよく、上記の性質を満たしていることが確認できた。 名前 コメント