統計データ解析I

第2講サンプルコード

Published

April 25, 2025

ベクトルの操作

問題

以下に示すベクトルを作成しなさい．

1から10までの自然数のベクトル
1以上30以下の奇数を昇順に並べたベクトル
すべての要素が1からなる長さ10のベクトル

作成したベクトルを操作しなさい．

ベクトルの長さを求める．
3番目の要素を取り出す．
最後の要素を取り出す．

解答例

1から10までの自然数のベクトルを作成する．

(foo <- 1:10) # () は表示．"foo <- 1:10; print(foo)" と等価

 [1]  1  2  3  4  5  6  7  8  9 10

1以上30以下の奇数を昇順に並べたベクトルを作成する．

(bar <- seq(from = 1, to = 30, by = 2)) # 変数名を含めた丁寧な書き方

 [1]  1  3  5  7  9 11 13 15 17 19 21 23 25 27 29

すべての要素が1からなる長さ10のベクトルを作成する．

(baz <- rep(1, length.out = 10)) # 関数rep_len() を用いても良い

 [1] 1 1 1 1 1 1 1 1 1 1

rep(1, len = 10)                 # 変数名は見分けられれば短くてもOK

 [1] 1 1 1 1 1 1 1 1 1 1

ベクトルの長さを求める．

length(foo)

[1] 10

3番目の要素を取り出す．

bar[3]

[1] 5

最後の要素を取り出す．

bar[length(bar)] # 長さを用いて最後の要素番号を指定すれば良い

[1] 29

tail(bar, n = 1) # 関数tail() を用いる方法もある

[1] 29

Note

関数 rep() や関数 seq() には高速化された派生型がいくつかあるので調べてみよ．

行列の操作

問題

以下に示す行列を作成しなさい．

\begin{equation} M = \begin{pmatrix} 1 & 2 & 3 \\ 4 & 5 & 6 \end{pmatrix} \end{equation}

行列を操作しなさい．

2行2列成分を取り出す．
転置行列を作成する．
行名をつける．

解答例

行列を作成する．

(M <- matrix(1:6, 2, 3, byrow = TRUE))

     [,1] [,2] [,3]
[1,]    1    2    3
[2,]    4    5    6

(M2 <- rbind(1:3, 4:6))                # cbind を使って作ることもできる

     [,1] [,2] [,3]
[1,]    1    2    3
[2,]    4    5    6

2行2列成分を取り出す．

M[2,2]

[1] 5

転置行列を作成する．

matrix(1:6, 3, 2)  # 愚直に作成し直す

     [,1] [,2]
[1,]    1    4
[2,]    2    5
[3,]    3    6

t(M)               # 関数t()を調べてみよ

     [,1] [,2]
[1,]    1    4
[2,]    2    5
[3,]    3    6

行名をつける．

rownames(M) <- c("first","second")
M                                  # print(M) と同義

       [,1] [,2] [,3]
first     1    2    3
second    4    5    6

列名をつけることもできる．

colnames(M2) <- c("alpha","beta","gamma")
M2

     alpha beta gamma
[1,]     1    2     3
[2,]     4    5     6

ベクトルと行列の計算

問題

以下の問いに答えなさい．

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}

解答例

1から10までの2乗値からなるベクトルを作成する．

1:10        # 1から10までのベクトル

 [1]  1  2  3  4  5  6  7  8  9 10

1:10 * 1:10 # Hadamard積を利用する

 [1]   1   4   9  16  25  36  49  64  81 100

1から10までの和を計算する．

1:10 %*% rep(1, 10) # (1,2,...,10)と(1,1,...,1)の内積

     [,1]
[1,]   55

sum(1:10)           # 関数sum()はベクトルの成分の合計を計算することもできる

[1] 55

九九の表を作成する．

matrix(rep(1:9, 9), 9, 9) # 行ごとに1から9を並べる

      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
 [1,]    1    1    1    1    1    1    1    1    1
 [2,]    2    2    2    2    2    2    2    2    2
 [3,]    3    3    3    3    3    3    3    3    3
 [4,]    4    4    4    4    4    4    4    4    4
 [5,]    5    5    5    5    5    5    5    5    5
 [6,]    6    6    6    6    6    6    6    6    6
 [7,]    7    7    7    7    7    7    7    7    7
 [8,]    8    8    8    8    8    8    8    8    8
 [9,]    9    9    9    9    9    9    9    9    9

matrix(rep(1:9, 9), 9, 9, byrow = TRUE) # 列ごとに1から9を並べる

      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
 [1,]    1    2    3    4    5    6    7    8    9
 [2,]    1    2    3    4    5    6    7    8    9
 [3,]    1    2    3    4    5    6    7    8    9
 [4,]    1    2    3    4    5    6    7    8    9
 [5,]    1    2    3    4    5    6    7    8    9
 [6,]    1    2    3    4    5    6    7    8    9
 [7,]    1    2    3    4    5    6    7    8    9
 [8,]    1    2    3    4    5    6    7    8    9
 [9,]    1    2    3    4    5    6    7    8    9

matrix(rep(1:9, 9), 9, 9) * matrix(rep(1:9, 9), 9, 9, byrow = TRUE)

      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
 [1,]    1    2    3    4    5    6    7    8    9
 [2,]    2    4    6    8   10   12   14   16   18
 [3,]    3    6    9   12   15   18   21   24   27
 [4,]    4    8   12   16   20   24   28   32   36
 [5,]    5   10   15   20   25   30   35   40   45
 [6,]    6   12   18   24   30   36   42   48   54
 [7,]    7   14   21   28   35   42   49   56   63
 [8,]    8   16   24   32   40   48   56   64   72
 [9,]    9   18   27   36   45   54   63   72   81

#' 行列とその転置を利用してもよい
(M <- matrix(rep(1:9, 9), 9, 9)) # 9x9型行列の積による方法

      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
 [1,]    1    1    1    1    1    1    1    1    1
 [2,]    2    2    2    2    2    2    2    2    2
 [3,]    3    3    3    3    3    3    3    3    3
 [4,]    4    4    4    4    4    4    4    4    4
 [5,]    5    5    5    5    5    5    5    5    5
 [6,]    6    6    6    6    6    6    6    6    6
 [7,]    7    7    7    7    7    7    7    7    7
 [8,]    8    8    8    8    8    8    8    8    8
 [9,]    9    9    9    9    9    9    9    9    9

M * t(M)

      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
 [1,]    1    2    3    4    5    6    7    8    9
 [2,]    2    4    6    8   10   12   14   16   18
 [3,]    3    6    9   12   15   18   21   24   27
 [4,]    4    8   12   16   20   24   28   32   36
 [5,]    5   10   15   20   25   30   35   40   45
 [6,]    6   12   18   24   30   36   42   48   54
 [7,]    7   14   21   28   35   42   49   56   63
 [8,]    8   16   24   32   40   48   56   64   72
 [9,]    9   18   27   36   45   54   63   72   81

(M9 <- matrix(1:9,9,1))          # 9x1型行列の積による別の方法

      [,1]
 [1,]    1
 [2,]    2
 [3,]    3
 [4,]    4
 [5,]    5
 [6,]    6
 [7,]    7
 [8,]    8
 [9,]    9

M9 %*% t(M9)

      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
 [1,]    1    2    3    4    5    6    7    8    9
 [2,]    2    4    6    8   10   12   14   16   18
 [3,]    3    6    9   12   15   18   21   24   27
 [4,]    4    8   12   16   20   24   28   32   36
 [5,]    5   10   15   20   25   30   35   40   45
 [6,]    6   12   18   24   30   36   42   48   54
 [7,]    7   14   21   28   35   42   49   56   63
 [8,]    8   16   24   32   40   48   56   64   72
 [9,]    9   18   27   36   45   54   63   72   81

(V9 <- 1:9); V9 %o% V9           # 講義の範囲外の方法 (%o% はベクトルの積の一つ)

[1] 1 2 3 4 5 6 7 8 9

      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
 [1,]    1    2    3    4    5    6    7    8    9
 [2,]    2    4    6    8   10   12   14   16   18
 [3,]    3    6    9   12   15   18   21   24   27
 [4,]    4    8   12   16   20   24   28   32   36
 [5,]    5   10   15   20   25   30   35   40   45
 [6,]    6   12   18   24   30   36   42   48   54
 [7,]    7   14   21   28   35   42   49   56   63
 [8,]    8   16   24   32   40   48   56   64   72
 [9,]    9   18   27   36   45   54   63   72   81

30度の回転行列の2乗は60度の回転行列であることを確かめる．

theta <- pi/6  # 30度のラジアン値
(R30 <-        # 30度の回転行列
     matrix(c(cos(theta),-sin(theta),
              sin(theta), cos(theta)), 2, 2, byrow = TRUE))

          [,1]       [,2]
[1,] 0.8660254 -0.5000000
[2,] 0.5000000  0.8660254

(R60 <-        # 60度の回転行列
     matrix(c(cos(2*theta),-sin(2*theta),
              sin(2*theta), cos(2*theta)), 2, 2, byrow = TRUE))

          [,1]       [,2]
[1,] 0.5000000 -0.8660254
[2,] 0.8660254  0.5000000

R30 %*% R30    # 30度の回転行列の2乗を表示

          [,1]       [,2]
[1,] 0.5000000 -0.8660254
[2,] 0.8660254  0.5000000

問題

以下の問いに答えなさい．

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

解答例

1から10までの2乗値からなるベクトルを作成する．

(1:10)^2 # ^2も関数として成分ごとに計算される

 [1]   1   4   9  16  25  36  49  64  81 100

例題と同様に A と \boldsymbol{b} を作成する．

(A <- matrix(rnorm(9), 3, 3) + diag(rep(1, 3)))

          [,1]      [,2]        [,3]
[1,] 1.2944115 0.8964453 -0.34020354
[2,] 0.6400273 1.1305498 -0.08463192
[3,] 1.0028176 1.4576977  0.50870264

(b <- 1:3)

[1] 1 2 3

エラーになる理由を考察する．

A %*% b           # 列ベクトル (3x1型行列)

         [,1]
[1,] 2.066691
[2,] 2.647231
[3,] 5.444321

b %*% A           # 行ベクトル (1x3型行列)

         [,1]     [,2]     [,3]
[1,] 5.582919 7.530638 1.016641

A %*% b + b %*% A # 異なる次元(大きさ)の行列は足し算できない

Error in A %*% b + b %*% A: non-conformable arrays

Tip

計算途中でエラーが生じても，出力を停めないようにするには以下の chunk option を指定すれば良い(ソースコードを参照)．

#| error: true

回転してもベクトルの長さが変わらないことを確認する．

#' 回転行列とベクトルを作成する (好きに設定してよい)
theta <- 2*pi/3 # 120度のラジアン値
(R <- matrix(c(cos(theta),-sin(theta),
               sin(theta), cos(theta)), 2, 2, byrow = TRUE))

           [,1]       [,2]
[1,] -0.5000000 -0.8660254
[2,]  0.8660254 -0.5000000

(x <- 1:2)

[1] 1 2

(y <- R %*% x) # xを回転してyを作成 (結果は行列になるので扱いに注意が必要)

           [,1]
[1,] -2.2320508
[2,] -0.1339746

#' 長さを確認する
x %*% x                       # xの長さの2乗はベクトルの内積で計算できる

     [,1]
[1,]    5

as.vector(y) %*% as.vector(y) # yの長さの2乗 (ベクトルに変換して計算)

     [,1]
[1,]    5

t(y) %*% y                    # 行列の積として計算することもできる

     [,1]
[1,]    5

sum(y*y)                      # 行列の要素積を計算して全ての成分の和を計算

[1] 5

sum(y^2)                      # 上記は各成分の2乗和として計算することも可能(xも同様)

[1] 5