「繰り返し処理によるシミュレーション」の編集履歴(バックアップ)一覧に戻る

繰り返し処理によるシミュレーション - (2008/08/13 (水) 00:58:51) のソース

関数の使い方(参考:[[Rの基本的な使い方]])とベクトルの扱い(参考:[[ベクトルと代入]])がある程度理解できてきたら、あとはそれに繰り返し処理の知識を加えることによってシミュレーションができるようになります。

*繰り返し処理
**for文
まずはfor文です。というかそもそも「繰り返し処理」とは何ぞという疑問があるかも知れません。簡単に言えば、「どこかの数値をちょっとずつ変えながら同じような命令を何度も実行させる」ということです。例を見て、実際に実行させてみれば理解できると思います。

for文はこのような構造をしています。
 for(i in M){式}
ここでiというのが繰り返しのたびに変化させる数値で、「変数」といいます。そして、inの後ろのMですが、これはベクトルでなければなりません。このベクトルの要素が最初から順番にiの中に代入されていきます。そしてそのiを使って(あるいは使わなくてもいいのですが)中括弧{}の中の式が繰り返し実行されます。この中括弧は別に無くてもいいのですが、繰り返し処理の範囲をはっきりさせるためにもちゃんと書きましょう。特に式が複数行に渡る場合は重要です。
 for(i in M){
    式その1
    式その2
    式その3
    }
ここで括弧の中が何文字か字下げされていますが、Rエディタ上でTABキーを押すとこのように何文字かまとめてスペースを空けることができます。中括弧の中をこのように字下げすると見やすくなりますし、なにやらプログラミングをしているっぽくなりますのでぜひとも活用しましょう。

抽象的な説明ばかりではイメージもしにくいと思いますので、1~100までの総和をfor文によって計算してみましょう。
 x <- 0
 for(i in 1:100){
    x <- x + i
    }
これで、xの中には1~100までの総和の答えが入ります。

ここでは、iという変数を1~100まで順につかい、x+iという数値を計算してxに代入しています。妙な書き方に見えるかもしれませんが、これは「前のx」にiをプラスして「新しいx」を作る操作だと思ってください。つまり、繰り返し処理の間「古いx」が「新しいx」に次々と置き換えられていくわけです。

このような書き方をする場合、xというオブジェクトがあらかじめ存在していなければなりません。そこでfor文の前に
 x <- 0
という代入操作をして、数値の0を持ったxを作り出しているわけです。

今の例ではiという変数を実際に数値とみなして計算の一部に組み入れましたが、変数を添字とみなす使い方もよくやります。
 for(i in 1:100){
    x[i] <- i
    }
これで、xは1~100までの数値をひとつずつ要素として持つ長さ100のベクトルになりました。これは「xのi番目にiという数値を入れておくれ」という命令です。ところで、この命令は一回ループが回るたびにベクトルの長さがひとつ増えるわけです。このような作業は若干時間がかかります(100程度では気になりませんが、100万とか1億とかのオーダーになると顕著です)。ですから、あらかじめ目的の長さのベクトルを作成しておいてそこへ代入する、といったことがよく行われます。これは1行最初に追加するだけです。
 x <- numeric(100)
これで、xは0が100個並んだベクトルになりました。numeric()関数は引数に指定された個数だけ0の並んだベクトルを作成する関数です。空っぽの容器を作るときにしばしば使います。

ところで今上げた2つの例は何も繰り返し処理を使わなくてもベクトルを使えば簡単に計算できますし、簡単に作れます。
 sum(1:100)
 1:100
繰り返し処理というのは大変時間のかかる作業です(といってもせいぜい数秒~数分の話なんですが)。ですからベクトルでまとめて処理できるような場合はなるべくベクトルでまとめて処理するようにしましょう。

**for文を使ったシミュレーション事例:ブートストラップ法
ブートストラップと呼ばれる統計の手法があります。これは標本のデータから重複を許したサンプリング(普通は標本数と同じだけ)を行い、新たな標本を作製する、といった作業を何度も繰り返す手法です。この標本からのサンプリングをリサンプリングと呼びます。こうして作られた複数の標本から計算される統計量(例えば平均、分散)のばらつき方は、母集団からサンプリングを何度も繰り返した時のばらつき方に近いという性質があります。つまり、複雑な確率密度関数や難解な中心極限定理を使うことなく、実際のデータのばらつきを観察することで推定値や信頼区間を決定することができるのです。

なにはともあれやってみましょう。まずは標本データのベクトルを用意します。
 x <- c(1, 2, 3, 5, 3, 4, 4, 7, 8, 10, 1)
全部で10の要素からなるデータです。ここからリサンプリングを行うには、sample()関数を使います。sample()関数には3つの引数を指定します。
 sample(x, 10, replace=T)
第一引数はサンプリングの対象となるデータベクトル、第二引数はそこからいくつのデータをサンプリングするか、第三引数は重複抽出を許すか否か(T or TRUEで許す、F or FALSEで許さない)を指定します。そうしてこれを何度も繰り返すわけですから、