Rのベクトルと行列、演算について(2)

2023年8月23日 2023年8月20日

FPshima

前回の続きです。行列の演算を見ていきましょう

行列を準備する

行列は特殊なものが多くあって、色々と面倒です。個別の話はしませんが、のちのことを考えて、対角成分が連番になるような行列を作ります。　ちなみにdiag()はn行n列の単位行列を、作ります。

diag(4)　　　　　　　　　　#　4行4列の単位行列　(大きさが1)
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    0    1    0    0
[3,]    0    0    1    0
[4,]    0    0    0    1

a <- diag(1:4)　　　　　　#　対角成分が1:4の対角行列
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    0    2    0    0
[3,]    0    0    3    0
[4,]    0    0    0    4

逆行列を求めます。逆行列とは元の行列との積を取ると単位行列になるような行列のことです。

  I <- solve(a)
> I
     [,1] [,2]      [,3] [,4]
[1,]    1  0.0 0.0000000 0.00
[2,]    0  0.5 0.0000000 0.00
[3,]    0  0.0 0.3333333 0.00
[4,]    0  0.0 0.0000000 0.25
> I%*%a                 # もとの行列aとの積が4行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
> a%*%I　　　　　　　　　　#　どっちから積を求めても単位行列になる
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    0    1    0    0
[3,]    0    0    1    0
[4,]    0    0    0    1

対角成分を除く上三角成分をすべて置換する。または下三角成分を置換する関数

> a[upper.tri(a)] <- 1　　#　上三角成分を1に置換
> a
     [,1] [,2] [,3] [,4]
[1,]    1    1    1    1
[2,]    0    2    1    1
[3,]    0    0    3    1
[4,]    0    0    0    4

> a<-diag(1:4)　　　　　#　aを作り直した
> a[lower.tri(a)] <- 1　#　下三角成分を1に置換
> a
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    1    2    0    0
[3,]    1    1    3    0
[4,]    1    1    1    4

> lower.tri(a)              # つまりこういう事。で、aの要素を[]内で指定して格納してます
      [,1]  [,2]  [,3]  [,4]
[1,] FALSE FALSE FALSE FALSE
[2,]  TRUE FALSE FALSE FALSE
[3,]  TRUE  TRUE FALSE FALSE
[4,]  TRUE  TRUE  TRUE FALSE

行列の(i,j)を(j,i)に入れ替える、いわゆる転地行列はｔ()でできる。

 
 > a
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    1    2    0    0
[3,]    1    1    3    0
[4,]    1    1    1    4
> t(a)
     [,1] [,2] [,3] [,4]
[1,]    1    1    1    1
[2,]    0    2    1    1
[3,]    0    0    3    1
[4,]    0    0    0    4

次に、行列式を求めてみましょう。　行列式とは、n次元の単位空間が、どの程度変化しているかを示します。（何を言っているかわからないけど、n次元は人間には捉えられないので気にしないようにしましょう。。。）

行列式はdet()で求められます。　(英語でdeterminantから）

他にcrossprod(x,y)という関数も用意されています。これはt(x)%*%yとおなじです。

 a
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    1    2    0    0
[3,]    1    1    3    0
[4,]    1    1    1    4
det(a)
[1] 24

> crossprod(a,I)
     [,1] [,2]      [,3] [,4]
[1,]    1  0.5 0.3333333 0.25
[2,]    0  1.0 0.3333333 0.25
[3,]    0  0.0 1.0000000 0.25
[4,]    0  0.0 0.0000000 1.00
> t(a)%*%I
     [,1] [,2]      [,3] [,4]
[1,]    1  0.5 0.3333333 0.25
[2,]    0  1.0 0.3333333 0.25
[3,]    0  0.0 1.0000000 0.25
[4,]    0  0.0 0.0000000 1.00

最後に行列の固有値と、固有ベクトルを求めましょう。ただし、固有値も固有ベクトルも、あったりなかったりします。

固有ベクトルとは、ある行列で線形変換しても、方向が変わらないようなベクトルのこと…

固有値とは、固有ベクトルが何倍になったかなーって値です。　ある行列をA　固有ベクトルをa　固有値をλとしたとき、

Ａa =λa

導出方法とかは数学系のサイトに譲るとして、Rでは以下のように書きます。

> a
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    1    2    0    0
[3,]    1    1    3    0
[4,]    1    1    1    4
> eigen(a)
eigen() decomposition
$values
[1] 4 3 2 1

$vectors
     [,1]       [,2]       [,3]       [,4]
[1,]    0  0.0000000  0.0000000  0.7071068
[2,]    0  0.0000000  0.7071068 -0.7071068
[3,]    0  0.7071068 -0.7071068  0.0000000
[4,]    1 -0.7071068  0.0000000  0.0000000

> det(a)
[1] 24
> e <- eigen(a)
> e$values[1]*e$values[2]*e$values[3]*e$values[4]
[1] 24　　#　固有値の積はdet(a)行列式と一致します

> b <- c(3,0,1,2)
> c <- matrix(b,nrow = 2,ncol = 2)
> c
     [,1] [,2]
[1,]    3    1
[2,]    0    2
> eigen(c)
eigen() decomposition
$values
[1] 3 2

$vectors
     [,1]       [,2]
[1,]    1 -0.7071068
[2,]    0  0.7071068

#　Ａa =λa　(定義の式)を確認してみましょう
> c%*%e$vectors[,1]
     [,1]
[1,]    3
[2,]    0
> c%*%e$vectors[,2]
          [,1]
[1,] -1.414214
[2,]  1.414214
> e$values[1]%*%e$vectors[,1]
     [,1] [,2]
[1,]    3    0
> e$values[2]%*%e$vectors[,2]
          [,1]     [,2]
[1,] -1.414214 1.414214
#　どちらも一致しました。なお、固有値の並び(順番)と、固有ベクトルの[,並び(順番)]が対応しています。

まぁ　直接値を求めることは、あまりないかもしれませんが、一応こんな感じの結果になります。

行列を準備する

おすすめFXアイテム