このページは編集中です。完成予定は6月26日午後4時です。　Rを使ったプログラミング演習2：　シミュレーション †

↑

第10回授業の獲得目標： [worried] †

↑

条件分岐 [smile] †

プログラミングの基本技の最後は条件分岐。if()を使って表現する。

i=4                       #iに4を代入
if(i==3){                 #iの値が3ならば
  print("三")            #   三 と表示する
} else {                  #iの値が3では無ければ
  print("三以外")   #   三以外と表示する
}

if命令の()の中の式を条件式といい、iの値を評価している。評価に使われるのは比較演算子というもので、

==　    等しければ
!=      等しく無ければ
>, >=   左辺が右辺より大きいなら、左辺が右辺以上ならば
<, <=   左辺が右辺より小さいなら、左辺が右辺以下ならば

ということを意味している

''比較演算子''である''=='' に対して、代入の時に使った ''='' は''代入演算子''といいう。

上のif命令では、次のような条件分岐が行われている。

if(i==3){
   print("三")
}
 ・もしiの値が３ならば、画面に三を表示。
 ・もしiの値が３で無いならば、何もしない
else {                  #iの値が3では無ければ
  print("三以外")   #   三以外と表示する
}

#演習：乱数を１つ発生させ、値が0.5以上なら「0.5以上」、0.5より小さいなら「0.5より小さい」と表示させる この演習では、乱数を発生させるrunif()という関数を使う。試しに

runif(5)

と入力してみよう。0から1までの乱数が5個表示されたはずだ。乱数はシミュレーションではよく使うので、もう少し慣れ親しんでおこう。 例えば、runif()を使って乱数を10個発生させ、ヒストグラムを書いてみよう。同じことを、100個、1000個、10000個、100000個についてもやってみよう。

hist(runif(10))

バラツキがだんだんと無くなってきたのがわかるだろうか。

#穴埋め式練習問題
if (runif(1) >= _ ) {   #条件式0.5以上ならば
  print("0.5以上")
} _ { 
  print("0.5より小さい")
}
↑

繰り返しと条件分岐 †

　では、上の繰り返しと条件分岐を組み合わせて使ってみよう

練習問題：for命令を使って1から10までの数字を表示させなさい、5だけは"five"と英単語で表示させなさい
↑

関数の定義 [smile] †

「関数の定義」なんていうと難しそうだが、ようするに、先ほどまでに作ったプログラムに名前をつけて、いろいろと数値を変えて解析できるようにしようというものだ。こういうときに関数を定義する関数、'function()'を使う。 例えば、円の面積を計算する関数を作るなら、

 menseki=function(r){r*r*pai}

とすれば良い。これで自分独自の関数menseki()ができた。 では、実行してみよう。

menseki()

演習：1から入力した数値までの全てを横一列に表示させるプログラムを作りdisplayという名前の関数として定義する

display=function(a){       #関数定義の始まり
 kekka=c()                 #kekkaに空ベクトルを代入して初期化
 for(i in 1:a){              #a回（iの値を1からaまで変化させる間）繰り返し
   kekka=c(kekka,i)     #kekkaというベクトルにiを要素として代入
 }
 print(kekka)              #kekkaの内容を表示
}                          #関数定義の終わり
display(10)

↑を使って何回も実行してみると、結果がいろいろ変わるのがわかる。

演習上の関数の定義方法に従って、入力した数までの合計値を計算するsumupという名前の関数を作成する。下の囲みの中の_の部分（１文字に対応するとは限らない）を埋めて、プログラムを完成しなさい。

 sumup = ______       #関数sumpuを定義
   kotae = _          #kotaeを初期化
   _ (i in ___){      #1からaまで繰り返し
     kotae = _____    #kotaeにiの値を足したものとkotaeに代入
   }                  #繰り返し終了
  _____               #kotaeを表示
  }                   #関数定義の終了
↑

授業で使ったRの基本関数 †

授業で使うR関数一覧 †

↑

演算子 †

演算子意味使用例
+足し算> 4+3
[1] 7
-引き算> 4-3
[1] 1
*掛け算> 4*3
[1] 12
/割り算> 4/3
[1] 1.333333
^累乗> 4^3
[1] 64
%/%整数商> 7%/%3
[1] 2
%%剰余> 7%%3
[1] 1
↑

定数（Rの基本システムで決まっている値):　次の５つ †

> LETTERS
  [1] "A" "B" "C" "D" "E" "F" "G" "H" "I" "J" "K" "L" "M" "N" "O" "P" "Q" "R" "S" "T" "U"
 [22] "V" "W" "X" "Y" "Z"
> letters
 [1] "a" "b" "c" "d" "e" "f" "g" "h" "i" "j" "k" "l" "m" "n" "o" "p" "q" "r" "s" "t" "u"
 [22] "v" "w" "x" "y" "z"
> month.abb
  [1] "Jan" "Feb" "Mar" "Apr" "May" "Jun" "Jul" "Aug" "Sep" "Oct" "Nov" "Dec"
> month.name
  [1] "January"   "February"  "March"     "April"     "May"       "June"      "July"     
  [8] "August"    "September" "October"   "November"  "December" 
> pi
 [1] 3.141593
↑

関数 †

Rの関数意味使用例
c()ベクトルの作成x = c(1, 2, 3, 4)
mean()ベクトルの平均mean(1, 2, 3, 4, 5)
scan()コピー・ペーストや手入力でデータをベクトルとして読み込む。x=scan()
hist()ヒストグラムを表示x=c(1,1, 3, 1, 5, 3, 4, 4); hist(x)
summary()要約統計量の表示x=c(1,1, 3, 1, 5, 3, 4, 4); summary(x)
plot()散布図の表示x=c(1,1, 3, 1, 5, 3, 4, 4); plot(x)
print()オブジェクトを表示x=c(1,1, 3, 1, 5, 3, 4, 4); print(x)
for(ループ変数 in リスト) {式}ループ変数をリストの内容に従って、１つずつ変化させる。変化の1回ごとに、{ }の式を実行。for(i in 1:10) { print(i) } #注：1:10は 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 を表す
matrix(ベクトル, 行数, 列数)ベクトルを行列に変換して表示matrix(1:6, nrow=2, ncol=3)
if (条件式) {式1} else {式2}条件式が真の時は 式1を実行、偽のときは式2を実行if (runif(1) < 0.5) { print("0.5より小さい")} else { print("0.5以上")}
runif()0-1の乱数を指定した個数発生runif(5)
function(引数1, 引数2,....){式 }自分で関数を定義するmenseki=function(r){r*r*pi}
plot.new()図を新しく書き直すという関数
↑

Rを用いたシミュレーション †

↑

円周率(π）をモンテカルロ法で求める †

simpai=function(x){
	kaisu=0;  #回数の初期値を0にする
	for (i in 1:x) {
		point=runif(2)
		distant=sqrt(point[1]^2+point[2]^2)
		if(distant<1){
			kaisu=kaisu+1
		}
	}
	answer=4*kaisu/x
	print(answer)
}
> simpai(10)
[1] 2.8
> simpai(100)
[1] 3.16
> simpai(1000)
[1] 3.156
> simpai(10000)
[1] 3.1624
> simpai(100000)
[1] 3.15652
> simpai(1000000)
[1] 3.13994
> simpai(10000000)
[1] 3.140556