library(conflicted) # 関数名の衝突を警告
conflicts_prefer( # 優先的に使う関数を指定
dplyr::filter(),
dplyr::select(),
dplyr::lag(),
)
library(tidyverse)統計データ解析I
第3講 実習
準備
以下で利用する共通パッケージを読み込む．
関数の使い方
問題
以下の問に答えなさい．
ヘルプ機能 (Help タブの検索窓, または関数
help(),?) を用いて関数sample()を調べなさい．サイコロを 1 回振る試行を模擬せよ．
サイコロを 10 回振る試行を模擬せよ．
Tip引数
replaceを調べなさい．1 が出易いサイコロを作りなさい．
Tip引数
probを調べなさい．1 から 6 をランダムに並べ替えよ．
解答例
関数 sample() を調べる．
help(sample)サイコロを1回振る試行を模擬する．
sample(x = 1:6, size = 1) # x = 見本点(標本点)の集合，size = 試行回数[1] 3サイコロを10回振る試行を模擬する．
sample(x = 1:6, size = 10) # 復元しないのでエラーになるError in sample.int(length(x), size, replace, prob): cannot take a sample larger than the population when 'replace = FALSE'sample(x = 1:6, size = 10, replace = TRUE) [1] 1 6 1 3 3 2 5 4 1 51 が他の目の 3 倍出易いサイコロを作る．
sample(1:6, 10, # 引数名(x, size)を省略
replace = TRUE, # 復元抽出を行う
prob = c(3, rep(1,5))) # =c(3,1,1,1,1,1) 各見本点の確率(重みづけ) [1] 1 4 5 4 4 1 4 1 1 41 から 6 をランダムに並べ替える．
sample(1:6, 6) # 既定値では復元抽出しないので並べ替えになる[1] 4 2 5 3 1 6自作関数の定義
問題
以下の問に答えなさい．
1 から整数 n までの和を求める関数を作成せよ．
Tip関数
sum()を調べなさい． 等差数列の和を利用してもよい．整数 n の階乗 n! を求める関数を作成せよ
Tip関数
prod()を調べなさい．
解答例
和 (summation) を計算する関数を定義する．
my_sum <- function(n){
out <- sum(1:n) # 1からnまでの整数を生成して和を求める
return(out)
}
my_sum(10) # 1から10までの和を確認する[1] 55以下のように定義することもできる．
my_sum2 <- function(n){
out <- n*(n+1)/2 # 等差数列の和を利用した場合
return(out)
}
my_sum2(10)[1] 55階乗 (factorial) を計算する関数を定義する．
my_fact <- function(n){
out <- prod(1:n) # 1からnまでの整数を生成して積を求める
return(out)
}
my_fact(5) # 5!で確認する[1] 120
Tip
関数内でオブジェクトに代入せずに直接 return() に渡すこともできる．
my_sum <- function(n){return(sum(1:n))}また return() を省略すると，関数内の最後の計算結果が返値となるので，
my_sum <- function(n){sum(1:n)}としても同じ結果となる．このような短縮した表記は無名関数の定義などで活用される．
制御構造
問題
以下の問に答えなさい．
整数 n の Fibonacci 数を求める関数を作成せよ．
- Fibonacci 数は以下の漸化式で計算される． \begin{align} F_{0}&=0\\ F_{1}&=1\\ F_{n}&=F_{n-1}+F_{n-2} \end{align}
行列 X が与えられたとき，各列の平均を計算する関数を作成せよ．
前問で X がベクトルの場合にはその平均を計算するように修正せよ．
Tip関数
is.vector()が利用できる．
解答例
Fibonacci 数を返す関数を定義する．
my_fibo <- function(n){
f0 <- 0 # 第0項の設定
f1 <- 1 # 第1項の設定
if(n<0) {
print("計算できません")
return(NA) # 欠損値を返す
}
if(n==0) { # n=0の場合
return(f0)
}
if(n==1) { # n=1の場合
return(f1)
}
for(i in 2:n) { # n>=2の場合
fn <- f1 + f0 # fn = fn-1 + fn-2 の計算
f0 <- f1 # fn-2 の値の更新 (f0が覚えておく)
f1 <- fn # fn-1 の値の更新 (f1が覚えておく)
}
return(fn) # 計算結果を返す
}Fibonacci 数を n 項目まで計算する． n が小さい場合には以下のように愚直に行えばよい．
c(my_fibo(1),my_fibo(2),my_fibo(3),my_fibo(4),my_fibo(5),my_fibo(6))[1] 1 1 2 3 5 8同じ関数に複数の値を代入する方法はいくつか用意されているので， n が大きい場合には例えば以下のように行えばよい．
sapply(1:30, my_fibo) # sapplyを調べてみよ [1] 1 1 2 3 5 8 13 21 34 55
[11] 89 144 233 377 610 987 1597 2584 4181 6765
[21] 10946 17711 28657 46368 75025 121393 196418 317811 514229 832040行列の列ごとの平均を計算する関数を定義する．
my_col_ave <- function(X) {
ave <- rep(0, length = ncol(X)) # 平均を記録するベクトルを用意
for(i in 1:ncol(X)){ # 列ごとに計算
ave[i] <- sum(X[,i])/nrow(X) # 平均の定義に従って計算
## ave[i] <- mean(X[,i]) # 平均を計算する関数を用いても良い
}
return(ave)
}
(A <- matrix(1:12,3,4,byrow=TRUE)) # 適当な行列を定義する [,1] [,2] [,3] [,4]
[1,] 1 2 3 4
[2,] 5 6 7 8
[3,] 9 10 11 12my_col_ave(A) # 正しい答えを返す[1] 5 6 7 8想定しない型を入力してみる．
(x <- 1:12) # 適当なベクトルを定義する [1] 1 2 3 4 5 6 7 8 9 10 11 12my_col_ave(x) # うまく動かないError in 1:ncol(X): argument of length 0ベクトルと行列を扱えるように修正する．
my_col_ave <- function(X){
if(is.vector(X)){ # ベクトルとそれ以外の場合分け
ave <- mean(X) # ベクトルならば平均を計算
} else { # ベクトル以外は行列と想定して計算
ave <- rep(0,length=ncol(X))
for(i in 1:ncol(X)){
ave[i] <- mean(X[,i])
}
}
return(ave)
}
my_col_ave(A)[1] 5 6 7 8my_col_ave(x) # いずれも正しい答えを返す[1] 6.5