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 でエラー: 適切な引数ではありません
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と見なしてよく、上記の性質を満たしていることが確認できた。
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
eigen関数を使う。固有値、固有ベクトルを求めたい行列を引数に与えると、リストで、固有値と固有ベクトルが一緒に返ってくる。
以下に求めた例を示す。見てのとおり、対称行列でなくても求めることができる。
> mxaa <- matrix(c(3, 1, -2, 0), 2)
> print(mxaa)
[,1] [,2]
[1,] 3 -2
[2,] 1 0
> eigen(mxaa)
eigen() decomposition
$values
[1] 2 1
$vectors
[,1] [,2]
[1,] 0.8944272 0.7071068
[2,] 0.4472136 0.7071068
valuesが固有値、vectorsが固有ベクトル。例えば1つ目の固有ベクトルとその固有ベクトルを取り出すには、以下のようにする。
> eigen(mxaa)$values[1]
[1] 2
> eigen(mxaa)$vectors[1, ]
[1] 0.8944272 0.7071068
2つ目の固有値と固有ベクトル。
> eigen(mxaa)$values[2]
[1] 1
> eigen(mxaa)$vectors[2, ]
[1] 0.4472136 0.7071068
以下、2つの行列の固有値と固有ベクトルを求めた例。
> mxaa <- matrix(c(6, -1, -3, -3, 4, 3, 5, -5, -4), 3)
> print(mxaa)
[,1] [,2] [,3]
[1,] 6 -3 5
[2,] -1 4 -5
[3,] -3 3 -4
> eigen(mxaa)
eigen() decomposition
$values
[1] 3 2 1
$vectors
[,1] [,2] [,3]
[1,] 7.071068e-01 -0.3015113 -3.400738e-15
[2,] 7.071068e-01 -0.9045340 8.574929e-01
[3,] 2.706169e-15 -0.3015113 5.144958e-01
> mxaa <- matrix(c(3, 0, 0, -1, 2, 0, 1, 1, 3), 3)
> print(mxaa)
[,1] [,2] [,3]
[1,] 3 -1 1
[2,] 0 2 1
[3,] 0 0 3
> eigen(mxaa)
eigen() decomposition
$values
[1] 3 3 2
$vectors
[,1] [,2] [,3]
[1,] 1 0.0000000 0.7071068
[2,] 0 0.7071068 0.7071068
[3,] 0 0.7071068 0.0000000
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
上記の例のとおりに、特に何も指定しなければ第二引数が行列、第三引数が列数になる。ベクトルの要素数で行列どちらかが決まればもう片方も決まるため、行または列の省略は可能。ただし、ベクトルの要素数に応じた正しい値でなければならない。
最後の例のとおりに、指定した行数×列数の要素数に対してベクトルの要素数が足らないときは、そのベクトルの要素数が行数×列数の約数であれば、自動的に繰り返し使われる。
行ごと列ごとに角括弧（ブラケット）をつかってインデックス（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] でエラー: 添え字が許される範囲外です
最後の例のように、行列の各長さより大きい値を指定すると、エラーが発生する。
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
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
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
上記と同じ計算（AT A）はcrossprod関数を使うと同じ結果を得ることができる。
> crossprod(mx)
[,1] [,2] [,3]
[1,] 1 0 1
[2,] 0 2 0
[3,] 1 0 3
