ベクトル・行列の演算と関数

(Press ? for help, n and p for next and previous slide)

村田昇
2020.04.24

ベクトルの計算

ベクトルの表記

ベクトルは太字で，要素は下付き添字で表す

\begin{equation} \boldsymbol{a} =(a_{1},a_{2},\dotsc,a_{k}) \end{equation}
```
a <- c(a1,a2,...,ak) # k次元ベクトルの作成
```
別の書き方:

\begin{equation} (\boldsymbol{a})_{i} =\text{(ベクトル $\boldsymbol{a}$ の第 $i$ 成分)} \end{equation}

ベクトルの加法

同じ長さのベクトル の和および差:
数値の和と差のように扱うことができる

\begin{equation} \boldsymbol{a}\pm\boldsymbol{b} =(a_{1}\pm b_{1},a_{2}\pm b_{2},\dotsc,a_{k}\pm b_{k}) \end{equation}

\begin{equation} (\boldsymbol{a}\pm\boldsymbol{b})_{i} =a_{i}\pm b_{i} \end{equation}
```
a + b # 同じ長さのベクトル a,b の和
a - b # ベクトルの差
```

ベクトルの乗法

同じ長さの2つのベクトル の乗法:
(2種類ある)
- 成分ごとの積 (Hadamard積; 要素積)
- ベクトルの内積
データ解析ではどちらも良く用いられる

Hadamard積

同じ長さのベクトル の成分ごとの積

\begin{equation} \boldsymbol{a}\circ\boldsymbol{b} =(a_{1}b_{1},a_{2}b_{2},\dotsc,a_{k}b_{k}) \end{equation}

\begin{equation} (\boldsymbol{a}\circ\boldsymbol{b})_{i} =a_{i}b_{i} \end{equation}
```
a * b # ベクトルの成分ごとの積
a / b # 成分ごとの商も計算可
```

内積

同じ長さのベクトル の内積

\begin{align} \boldsymbol{a}\cdot\boldsymbol{b} &=a_{1}b_{1}+a_{2}b_{2}+\dotsb+a_{k}b_{k}\\ &=\sum_{i=1}^{k}a_{i}b_{i} \end{align}
```
a %*% b # ベクトルの内積
```

行列の計算

行列の表記

行列は大文字で，要素は下付き添字で表す

\begin{equation} A = \begin{pmatrix} a_{11}&a_{12}&\dots&a_{1n}\\ a_{21}&a_{22}&\dots&a_{2n}\\ \vdots&&\ddots&\vdots\\ a_{m1}&a_{m2}&\dots&a_{mn} \end{pmatrix} \end{equation}
```
A <- matrix(c(a11,a21,...,amn),m,n) # m x n 行列の作成
```
別の書き方:

\begin{equation} (A)_{ij} =\text{(行列 $A$ の $ij$ 成分)} \end{equation}

行列の加法

同じ大きさの行列 の和および差:
ベクトルと同じように記述することができる

\begin{equation} (A\pm B)_{ij}=a_{ij}\pm b_{ij} \end{equation}
```
A + B # 同じサイズの行列の和
A - B # 行列の差
```

行列の乗法

2つの行列の乗法:
(2種類ある)
- 同じ大きさの行列 の成分ごとの積 (Hadamard積; 要素積)
- $n\times m$ 型行列と $m\times l$ 型行列 の積
データ解析ではどちらも良く用いられる

Hadamard積

同じ大きさの行列 の成分ごとの積

\begin{equation} (A\circ B)_{ij}=a_{ij}b_{ij} \end{equation}
```
A * B # 行列の成分ごとの積
A / B # 成分ごとの商も計算可
```

行列の積

$n\times m$ 型行列 $A$ と $m\times l$ 型行列 $B$ の積

\begin{equation} (AB)_{ij}=\sum_{k=1}^{m}a_{ik}b_{kj} \quad\text{($AB$は$n\times l$行列)} \end{equation}
```
A %*% B # 行列の積
```

行列式とトレース

$n$ 次正方行列 $A$ の行列式 $\det(A)$
```
det(A) # 行列式
```
$n$ 次正方行列 $A$ のトレース(対角成分の総和)

\begin{equation} \mathrm{trace}(A)=\sum_{i=1}^{n}a_{ii} \end{equation}

関数は用意されていないので以下を利用:
- 関数 diag(): 行列の対角成分を取り出す(ベクトル)
- 関数 sum(): ベクトルの総和を計算する
```
sum(diag(A)) # 行列のトレース
```

演習

例題

適当な2次正方行列 $A$ で Hamilton-Cayleyの定理

\begin{equation} A^2-\mathrm{trace}(A)A+\det(A)E_{2}=O_{2} \end{equation}

の成立を確認せよ．
ただし $E_{2}$ は2次単位行列，$O_{2}$ は2次正方零行列

解答例

## 行列を作成 (好きに設定してよい)
(A <- matrix(1:4,2,2)-diag(rep(3,2)))


     [,1] [,2]
[1,]   -2    3
[2,]    2    1

## 左辺を計算
A%*%A - sum(diag(A)) * A + det(A) * diag(rep(1,2))


             [,1]         [,2]
[1,] 1.776357e-15 0.000000e+00
[2,] 0.000000e+00 1.776357e-15

練習問題

1から10の2乗値からなるベクトルを作成せよ．
1から10までの和を計算せよ．
行列を用いて九九の表を作成せよ．
30度の回転行列を2回乗ずると60度の回転行列となることを確認せよ．

\begin{equation} \text{(回転行列)} = \begin{pmatrix} \cos(\theta)&-\sin(\theta)\\ \sin(\theta)& \cos(\theta) \end{pmatrix} \end{equation}

ベクトルと行列の計算

ベクトルと行列の乗法

列(縦)ベクトル・行(横)ベクトルという 区別はない
行列とベクトルの順序で適切に判断される

計算結果は行列で表現される

A <- matrix(1:4,2,2); b <- c(5,6) # 行列とベクトルを作成

A %*% b # 行列 x ベクトル = 列ベクトル

     [,1]
[1,]   23
[2,]   34

b %*% A # ベクトル x 行列 = 行ベクトル

     [,1] [,2]
[1,]   17   39

連立1次方程式の解法

連立1次方程式
- $A$ : $n$ 次正則行列
- $\boldsymbol{b},\boldsymbol{x}$ : $n$ 次元列ベクトル
\begin{align} A\boldsymbol{x}&=\boldsymbol{b} &&\text{(連立1次方程式)}\\ \boldsymbol{x}&=A^{-1}\boldsymbol{b} &&\text{(解;$A$が正則な場合)} \end{align}
解を求めるには関数 solve() を利用する
```
x <- solve(A, b)    
```
ベクトル $\boldsymbol{b}$ の代わりに行列も扱える

逆行列

正則な $n$ 次正方行列 $A$ の逆行列 $A^{-1}$

\begin{equation} AA^{-1}=A^{-1}A=E_{n} \quad\text{($E_{n}$は$n$次単位行列)} \end{equation}
関数 solve() を利用して求めることができる

\begin{align} AX&=E_{n}\\ X&=A^{-1}E_{n}=A^{-1} \end{align}
```
solve(A,B) # AX=B の解Xを求める
solve(A) # 逆行列 (Bが単位行列の場合省略できる)
```
(他にもいくつか方法は用意されている)

関数の適用

ベクトルや行列に関数( $\sin,\exp,\dots$ など)を適用すると成分ごとに計算した結果が返される
ベクトル $\boldsymbol{a}$ ，行列 $A$ に関数 $\sin$ を適用

\begin{equation} (\sin(\boldsymbol{a}))_{i}=\sin(a_{i}) \end{equation}
```
sin(a) # 成分ごとに計算される．sin(a)[i]=sin(a[i])
```
\begin{equation} (\sin(A))_{ij}=\sin(a_{ij}) \end{equation}
```
sin(A) # 成分ごとに計算される．sin(A)[i,j]=sin(A[i,j])
```

演習

例題

適当な3次正方行列 $A$ と3次元ベクトル $\boldsymbol{b}$ を作成して $\boldsymbol{x}$ に関する連立1次方程式

\begin{equation} A\boldsymbol{x}=\boldsymbol{b} \end{equation}

を解け

解答例

## 行列とベクトルを作成 (好きに設定してよい)
## rnorm(9) は正規乱数を9つ作成する(第5回で詳しく説明)
(A <- matrix(rnorm(9),3,3)+diag(rep(1,3)))
(b <- 1:3)


           [,1]      [,2]      [,3]
[1,]  0.9259176 0.4112264 1.0078438
[2,] -0.7882785 2.0553422 0.5155167
[3,] -0.8777801 0.7914346 0.8761933

[1] 1 2 3

## 解を計算
(x <- solve(A,b))
A%*%x # 結果の確認(b になるはず)


[1] -1.28469800 -0.07189792  2.20182011

     [,1]
[1,]    1
[2,]    2
[3,]    3

練習問題

1から10の2乗値からなるベクトルを作成せよ．
2次元ベクトルを回転行列で変換しても長さが変わらないことを確かめよ．
例題の $A$ と $\boldsymbol{b}$ を用いて
```
A %*% b + b %*% A   
```
を計算するとエラーになるが，何故そうなるか理由を考えよ．

関数の定義

自作関数

他の言語と同様にRでも関数を定義できる
関数の定義には関数 function() を利用する

半径 r から球の体積と表面積を求める関数

myfunc <- function(r){
    V <- (4/3) * pi * r^3 # 球の体積
    S <- 4 * pi * r^2     # 球の表面積
    out <- c(V,S) # 返り値のベクトルを作る
    names(out) <- c("volume", "area") # 返り値の要素に名前を付ける
    return(out) # 値を返す
}
myfunc(1) # 実行


 volume     area 
4.18879 12.56637

制御構造

制御文

最適化や数値計算などを行うためには，条件分岐や繰り返しを行うための仕組みが必要となる
R言語を含む多くの計算機言語では
- if (条件分岐)
- for (繰り返し・回数指定)
- while (繰り返し・条件指定)
これを 制御文 と言う

`if` 文

条件Aが真のときプログラムXを実行する
```
if(条件A) プログラムX
```
上記の if 文に条件Aが偽のときプログラムYを実行することを追加する
```
if(条件A) プログラムX else プログラムY    
```

`if` 文の例

20200724が19で割り切れるか?

if(20200724 %% 19 == 0) {# %% は余りを計算
    print("割り切れます")
    print(20200724 %/% 19) # 商を表示
} else { # {}で囲まれたブロックが1つのプログラム
    print("割り切れません")
    print(20200724 %% 19) # 余りを表示
}


[1] "割り切れます"
[1] 1063196

`for` 文

ベクトル V の要素を順に変数 i に代入してプログラムXを繰り返し実行する
```
for(i in V) プログラムX
```
プログラムXは変数 i によって実行内容が変わる

`for` 文の例

アルファベットの20,15,11,25,15番目を表示

print(LETTERS) # LETTERS ベクトルの内容を表示
for(i in c(20,15,11,25,15)) {
    print(LETTERS[i]) # 順番に表示
}

 [1] "A" "B" "C" "D" "E" "F" "G" "H" "I" "J" "K" "L" "M"
[14] "N" "O" "P" "Q" "R" "S" "T" "U" "V" "W" "X" "Y" "Z"

[1] "T"
[1] "O"
[1] "K"
[1] "Y"
[1] "O"

`while` 文

条件Aが真である限りプログラムXを繰り返す
```
while(条件A) プログラムX
```
プログラムXは実行するたびに実行内容が変わり，いつか条件Aが満たされなくなるように書く
(repeat 文というものもある)

`while` 文の例

20200809を素因数分解する

n <- 20200809 # 分解の対象
p <- 2 # 最初に調べる数
while(n != 1){ # 商が1になるまで計算する
    for(i in p:n){ # pからnまで順に調べる
        if(n%%i == 0) { # 余りが0か確認
            print(i) # 割り切った数を表示
            n <- n/i # 商を計算して分解の対象を更新
            p <- i # 最初に調べる数を更新
            break # for文を途中で終了
        }  
    }
}


[1] 3
[1] 31
[1] 281
[1] 773

演習

例題

三角形の3辺の長さ $x,y,z$ を与えると面積 $S$ を計算する関数を作成せよ．
参考: ヘロンの公式 より

\begin{equation} S=\sqrt{s(s-x)(s-y)(s-z)},\qquad s=\frac{x+y+z}{2} \end{equation}

が成り立つ．

解答例

area <- function(x,y,z){
    s <- (x+y+z)/2
    S <- (s*(s-x)*(s-y)*(s-z))^(1/2)
    ## S <- sqrt(s*(s-x)*(s-y)*(s-z)) # 平方根を求める関数を用いても良い
  return(S)
}
area(3,4,5) # 直角三角形で検算
area(12,13,5)


[1] 6

[1] 30

演習

非負の整数 $n$ の階乗 $n!$ を求める関数を作成せよ．
整数 $n$ の Fibonacci数を求める関数を作成せよ．

\begin{align} F_{0}&=0\\ F_{1}&=1\\ F_{n}&=F_{n-1}+F_{n-2} \end{align}
行列 $X$ が与えられたとき，各列の平均を計算する関数を作成せよ．
前問で $X$ がベクトルの場合にはその平均を計算するように修正せよ．
(関数 is.vector() が利用できる)

ベクトル・行列の演算と関数

ベクトルの計算

ベクトルの表記

ベクトルの加法

ベクトルの乗法

Hadamard積

内積

行列の計算

行列の表記

行列の加法

行列の乗法

Hadamard積

行列の積

行列式とトレース

演習

例題

解答例

練習問題

ベクトルと行列の計算

ベクトルと行列の乗法

連立1次方程式の解法

逆行列

関数の適用

演習

例題

解答例

練習問題

関数の定義

自作関数

制御構造

制御文

if 文

if 文の例

for 文

for 文の例

while 文

while 文の例

演習

例題

解答例

演習

`if` 文

`if` 文の例

`for` 文

`for` 文の例

`while` 文

`while` 文の例