第
1
回
R”
アール
”
勉強会
片所 強
2007
年
8
月
19
日
1
Windows
版
R
のインストール
1.1
R
のダウンロード
The R projectというwebサイト(http://www.r-project.org/)より最新版のRを無 料でダウンロードして使用することができる。このwebサイトはすべて英語であるが、単 にダウンロードするだけであれば気にすることはない。左側のメニューより、Download の項目にある [CRAN]を選択すると、CRAN Mirrorsというページが表示される。そ こで Japanの項目にある大学のミラーページを選択する。そうしたら Download and Install Rという項目から[Windows(95 and later)]という部分をクリックし、さらに [base]という部分のリンクをクリックする。現在では、最新版は2.5.1なので、そのペー ジのタイトルもR - 2.5.1 for Windowsとなっているが、これは最新版のRが提供される たびに変更される。[R-2.5.1-win32.exe]の部分をクリックすればダウンロードが開始さ れる。
1.2
R
のインストール
Windowsにおいてソフトウェアのプログラムは難しいものではなく、Rにおいても、 原則として[次へ]か[OK]のボタンをどんどんクリックしていけばよい。途中で特に何か を指定する必要はないので、迷うこともないと思われる。万が一、インストール方法が分1.3
R
に関する有用な
web
サイトの紹介
1. RjpWiki(http://www.okada.jp.org/RWiki/) 2. R-Tips(http://cse.naro.affrc.go.jp/takezawa/r-tips/r.html) 3. Rによる統計処理(http://aoki2.si.gunma-u.ac.jp/R/) 4. R、R言語、R環境・・・・・・(http://www1.doshisha.ac.jp/ mjin/R/index.html) 原則として Rに関する情報のすべては参考URL1より得ることができる。前述したダ ウンロードからインストールをして使い始めるまでに生じる問題についても、この web サイトの情報を参考に解決できるだろう。参考URL2はRの基本的な処理が一通りまと められているので、1度は通して読んでおくとよい。参考URL3は群馬大学の青木繁信先 生のwebサイトの1部で、推定や検定、多変量解析、グラフィックスなどほぼすべての統 計処理はここにある関数によって行うことができる。参考URL4にはRに関するPDF 資料が豊富に用意されており、それらの内容も実に分かりやすいのでぜひ参照してもらい たい。2
とりあえず
R
を動かしてみよう
2.1
実に簡単な計算を
R
で実行する
RはMicrosoft Excelのように、ワークシートにデータを入力してクリックをしながら データ解析を行っていくようなものとは異なり、コマンド画面にユーザーがさまざまな命 令をコマンドしていくことによってデータを解析していくことになる。最初はキーボード 主体の動作に若干の違和感や不安を感じるかもしれないが、使用していくうちにR の便 利さを実感するだろうから、その点については安心してよい。 さっそく Rを使用して簡単な計算をしてみよう。例えば、普通の電卓で計算できる四 則演算や「ルートの計算」も行うことができる。画面の左側には常に”>”という記号があ るのが確認できるだろうが、これはR側が「命令を打ち込んでください」とユーザーがコ マンドするのを待っている状態なのである。この状態で次のようにコマンドしてみよう。> 12 + 3 - 67 #足し算と引き算 命令を打ち込んだらリターンキーを押す [1] -52 > 3 * 90 / 45 #掛け算と割り算 [1] 6 > sqrt(2) #ルートの計算 関数sqrt()を使用 [1] 1.414214 µ ´ 少し補足をしておくと、”#”はコメントであることをコンピュータに教えるための記号 で、このマーク以降の命令は無視されて、実行されることはない。また、今回は2の平方 根を求めるためにsqrt() という関数を用いた。関数とは引数にある値を指定することに よって、その指定した値に対して一定の処理を施した結果を返すためのものである。これ は5章で詳しく説明するので、今はコマンド方法のみを覚えておけば問題ない。 Rには最初からさまざまな関数が実装されているので、通常の電卓ではなく、関数電卓 と同じ機能を発揮して、より複雑な計算も簡単に行える。そのうちのいくつかを簡単に紹 介しておくが、もしかしたら見慣れないものもあるかもしれない。その場合はコメントに 記してあるキーワードで検索するなりしてほしい。もっとも、ここではこれらの数値が表 す意味を理解する必要はない。 ¶ ³ > 12 ^ 2 #累乗の計算 [1] 144 > sin(1) #sin"サイン"の計算 [1] 0.841471 > cos(5) #cos"コサイン"の計算 [1] 0.2836622 > tan(3) #tan"タンジェント"の計算 [1] -0.1425465 > log10(12) #常用対数(底が10の対数)の計算 [1] 1.079181 > exp(1) #ネイピア数
算に用いる+, -, *, /の記号の前後には半角スペースを入れたほうがよい。
2.2
代入、ベクトルと行列、データフレームの作成
代入とは「ある入れ物のなかに数値や文字を入れること」である。この「ある入れ物」 のことを変数という。通常、プログラミング言語には変数の型を宣言しなければならな い。例えば、整数を代入できる変数や文字を代入できる変数というものをあらかじめ宣言 しておかないと、その変数の中に特定の数値や文字列を代入することができない仕組みに なっている。しかし、Rにおいてはそのようなことを気にする必要はなく、たとえ整数と 小数だろうと、文字列だろうとお構いなしに代入することができる。 ¶ ³ > a <- 244 #整数値を変数aに代入 > a #変数の中身の表示 [1] 244 > b <- 0.12345 #小数値を変数bに代入 > b [1] 0.12345 > c <- "Katasho" #文字列を変数cに代入 > c [1] "Katasho" µ ´ また他のプログラミング言語においてはスカラー変数と配列変数とが区別されており、 前者は既にやったようにある変数に1つだけの数値や文字列が代入できる変数のことであ る。一方で配列変数とはある変数に1つ以上の複数の数値や文字列を代入することができ る変数のことである。これもあらかじめ変数宣言をしておかなければならないが、Rでは 特に宣言をする必要はない。> a <- c(1, 2, 3, 4, 5) #ベクトルとして変数に数値を代入 > a [1] 1 2 3 4 5 > b <- c("a", "b", "c", "d", "e") #ベクトルとして変数に文字を代入 > b [1] "a" "b" "c" "d" "e" > c <- c(1, "b", 3, "d", 5) #ベクトルとして数値と文字を代入? > c [1] "1" "b" "3" "d" "5" µ ´ 要するに複数の数値を1つの変数に代入したい場合は関数c()を用いて上のようにコマ ンドすればよい。ただし、気をつけなければいけないのは1つの変数に数値と文字を混同 させることはできないということである。上のcという変数には、命令としては数値と文 字を代入しているようであるが、実際に中身を表示させてみると、すべて文字列として扱 われているのである。もっとも、通常はこのような使い方をすることはまず在り得ない。 また、1回目に244という値を変数aに代入したが、今回新たに変数aに[1, 2, 3, 4, 5] という5つの数値を代入したので、以前に代入された244という数値は消えていることに も着目してもらいたい。同じ変数に新たな数値や文字列が代入されてしまったら、それ以 前に代入されていた値は消えて、元には戻らないので注意が必要である。それから、既に 使用されている関数名などのは変数として扱わない方が無難である。 なお蛇足ではあるが、数学的には配列とはいわずにベクトルというのが正しい。スカ ラーとは要素とか成分と訳されるが、ベクトルの中身のことを指す。例えば、今はaとい う変数に[1, 2, 3, 4, 5]という数値が格納されているが、この場合、「aは長さ5のベク トルである」という言い方をする。そして、「1番目の要素は1である」という言い方を する。 さて、ベクトルの次に学ぶべきことは行列である。行列といっても、ラーメン屋にでき る行列のことではなく、ここでいう行列とは「ベクトルが組み合わさったもの」とでも理 解しておけばよいだろう。下手に論ずるより以下の例を見たほうが速いだろう。
¶ ³ > dat <- matrix(c(
+ 1, 2, 3, + 4, 5, 6,
+ 7, 8, 9), nrow=3, ncol=3, byrow=T) > dat [,1] [,2] [,3] [1,] 1 2 3 [2,] 4 5 6 [3,] 7 8 9 > dat2 <- matrix(c( + 1, 2, 3, + 4, 5, 6, + 7, 8, 9), ncol=3) > dat2 [,1] [,2] [,3] [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9 µ ´ 行列を作成するための方法はいくつかあるが、もっとも単純なのは関数 matrix()を利 用することである。関数の機能は理解できなくてもよいから、どのようにしたらデータ行 列が作成できるのかを上の2つを見比べて理解してほしい。それでも少し説明しておく と、matrix()の引数であるnrow=3というのは3行という意味であり、ncol=3というの は3列という意味である。つまり、「3行3列のデータ行列を作成せよ」と命令している わけであるが、実質的にはncolのみを指定するだけで用は足りてしまう。byrowは転置 にするかどうかということなのだが、簡単にいえば、画面に入力したとおりにデータ行列 を作成したい場合はbyrow=T とすればよい。そうではなく、dat2のようにしたい場合 はこの引数を省略すればよい。 行列を作成するためのもう1つの方法を提示しておく。
> dat4 <- cbind(x1, x2, x3) #関数cbind()を利用 > dat4 x1 x2 x3 [1,] 1 4 7 [2,] 2 5 8 [3,] 3 6 9
> dat5 <- rbind(x1, x2, x3) #関数rbind()を利用 > dat5 [,1] [,2] [,3] x1 1 2 3 x2 4 5 6 x3 7 8 9 µ ´ さらにデータフレームの作成方法を提示しておく。データフレームを作成して「パスを 通しておく」ことによって、解析がより効率的になる。
¶ ³ > my.data <- data.frame(X1 = x1, X2 = x2, Y = x3) #デ ー タ フ レ ー ム の 作成 > my.data X1 X2 Y 1 1 4 7 2 2 5 8 3 3 6 9 > attach(my.data) #my.dataにパスを通す > X1 [1] 1 2 3 > X2 [1] 4 5 6 > Y [1] 7 8 9 > detach(my.data) #パスを解除 > X1 エラー: オブジェクト "X1" は存在しません µ ´ データフレームを作成したあとに、関数 attach()を利用すれば、データフレームにお ける変数名に直接的にアクセスが可能となる。通したパスを解除するためには detach() を用いるが、パスを解除したあとはデータフレームにおける変数名を表示させようとして も上例のようになる。基本的にそのデータフレームに対する解析が終了したら、すぐに detach()を用いてパスを解除するクセをつけたほうがよい。
2.3
基本統計量の算出
Rには平均値、中央値、分散、標準偏差、共分散に相関係数など様々な基本統計量を算 出するための関数を実装している。> dat <- rnorm(100, mean=50, sd=10) > mean(dat) #平均値 [1] 49.10156 > median(dat) #中央値 [1] 50.10624 > var(dat) #不偏分散 [1] 104.4191 > sd(dat) #不偏標準偏差 [1] 10.21857
> dat2 <- rnorm(100, mean=100, sd=20) > var(dat, dat2) #共分散 [1] -34.2106 > cor(dat, dat2) #相関係数 [1] -0.1832966 µ ´
3
心理学分野における古典的ないくつかの検定
3.1
t
検定
t検定というのはt分布を利用する検定手法の総称であって、正確には「独立2群の平 均値の差の検定」と「対応のある2条件間の差の検定」とに区別されている。しかし、こ こではこれら2種類の検定法についての説明はせずに、Rの操作方法のみを紹介する。 R でt 検定を行うためには t.test() という関数を用いる。対応ありかなしかは引数 pairedで指定することができる。省略した場合はpaired=FALSE、すなわち対応なしの t検定であるとされる。また、対応なしのt検定を行う場合は常にvar.equal=FALSEと することを奨める。¶ ³ > group1 <- c(56, 49, 59, 42, 58, 61, 53) #対照群のデータ
> group2 <- c(73, 86, 54, 60, 65, 70, 53) #実験群のデータ > t.test(group1, group2, var.equal=FALSE) #対応のないt検定
Welch Two Sample t-test
data: group1 and group2
t = -2.3395, df = 9.518, p-value = 0.04259
alternative hypothesis: true difference in means is not equal to 0 95 percent confidence interval:
-23.2278060 -0.4864797 sample estimates:
mean of x mean of y 54.00000 65.85714
> t.test(group1, group2, paired=TRUE) #対応のあるt検定
Paired t-test
data: group1 and group2
t = -2.2632, df = 6, p-value = 0.06427
alternative hypothesis: true difference in means is not equal to 0 95 percent confidence interval:
-24.6769481 0.9626624 sample estimates:
mean of the differences -11.85714
き学年によって走り幅跳びの記録の平均値に差が認められるかどうかを検定してみる。 ¶ ³ > #1年生のデータ > class.one <- c(3.9, 3.8, 3.6, 4.4, 4.7, 3.5, 3.7, 4.6, 3.9, 3.6) > #2年生のデータ > class.two <- c(5.8, 4.4, 4.5, 6.4, 4.9, 7.2, 4.8, 5.7, 5.6, 3.4) > #3年生のデータ > class.thr <- c(6.0, 5.8, 6.1, 5.6, 6.2, 6.0, 6.1, 5.7, 6.4, 5.9) > group <- as.factor(rep(c(1, 2, 3), c(10, 10, 10))) #グループを表す カテゴリカルデータ > #データフレームの作成 > my.data <- data.frame(
+ Y = c(class.one, class.two, class.thr), X = group)
> result <- summary(aov(Y ~ X, my.data)) #分散分析の結果を変数result へ代入
> result #結果の表示
Df Sum Sq Mean Sq F value Pr(>F) X 2 20.7807 10.3903 21.55 2.548e-06 *** Residuals 27 13.0180 0.4821 ---Signif. codes: 0 ’***’ 0.001 ’**’ 0.01 ’*’ 0.05 ’.’ 0.1 ’ ’ 1 µ ´ 分散分析には得られたデータ形式によってさまざまなタイプがあり、それによって解析 の仕方も異なる。また、今回と同じ形式のデータを分析する場合でもいくつかの方法があ る。実をいうと、分散分析の前にやったt検定も含めて、すべて同一の方法で(同じ線形 モデルとして)解析することができるが、ここではそれに触れることはしない。 Rにはこれ以外の検定を行うための関数が多く用意されており、最初から実装されてい る関数よりも場合によっては効率的な関数をweb上からダウンロードすることもできる
4
グラフィックス入門
RではExcelで描けるような棒グラフや折れ線グラフ、散布図などを描けることもでき る。そこでまずは散布図を描いてみよう。 ¶ ³ > x <- rnorm(100) > y <- rnorm(100) > plot(x, y) #もっとも単純な散布図のプロット> plot(x, y, xlim=c(-3,3), ylim=c(-3,3)) #x軸とy軸の最小・最大値を指 定
> plot(x, y, xlab="横軸", ylab="縦軸") #x 軸のラベルとy 軸のラベルを指 定
> plot(x, y, col="red", pch=18) #デ ー タ ポ イ ン ト の 色 (col) と 形 を 変 更 (pch) µ ´ R にはグラフを描くための関数がいくつか用意されているが、基本的にはplot()だけ で事足りるだろう。関数plot()の引数は上で紹介したもの以外にもあるが、複雑な図を 描くのでなければ、これだけで十分な場合も多い。 続いて今度は折れ線グラフを描いてみよう。これも散布図のときと同じく plot()を使 用して描くことができる。 ¶ ³ > x <- 1:6 > y <- c(1, 5, 6, 3, 9, 10) > plot(x, y, type="b") #もっとも単純な折れ線グラフ > plot(x, y, type="b", lty=2) #線のスタイルを指定
µ ´
今回はplot()の引数に、新たにtypeを指定したが、この他にもデータポイントを描か ずに線のみを描く(その場合はtype=”l”とする)ことや、イベントレコードのように階 段状に線を結ぶようなグラフを描く(その場合はtype=”s”とする)こともできる。また、 グラフの枠のみを描く(type=”n”)ことなどもでき、いったん、枠のみを描いてから点
#カテゴリが1つの棒グラフ
> x <- c(12, 34, 72, 56, 20)
> names(x) <- c("A1", "A2", "A3", "A4" , "A5") > barplot(x, col="blue") > abline(h=0) #カテゴリが2つ以上の棒グラフ > y <- c(17, 34, 69, 32, 12) > mat <- rbind(x, y) > mat A1 A2 A3 A4 A5 x 12 34 72 56 20 y 17 34 69 32 12
> barplot(mat, beside=TRUE, col=c("blue", "red")) > legend(2, 60, + paste(c("努力家","天才肌")), + fill=c("blue", "red")) #凡例を追加 µ ´ グラフィックスの補足で、心理学でよく用いられる群間比較(k群の平均値の差を検討 する場合)のプロット方法を提示しておく。例えば、今回は3条件間の平均値の違いを表 すグラフを描いてみる。 ¶ ³ > g1 <- rnorm(10, 10, 1) > g2 <- rnorm(10, 20, 1) > g3 <- rnorm(10, 5, 1) > #x軸を描かずにプロット
> plot(x, c(g1,g2,g3), xlab="GROUP", ylab="SCORE", xaxt="n") > axis(1, at=1:3, label=c("g1", "g2", "g3")) #x軸を新たに追加 > M <- c(mean(g1), mean(g2), mean(g3)) #各群の平均値
5
関数の定義
5.1
関数とは何か
既に述べたことだが、関数とは「引数にある値を指定すると、特定の処理に従ってある 値を返す機能をもったもの」のことである。関数は(1)関数名、(2)引数、(3)関数の中身 の3つで構成されていると考えると分かりやすい。例えば、平均を求める関数はmean() であるが、この場合、関数名は”mean”であり、”()”のなかに指定するものが引数である。 関数の中身とは要するに平均値を算出するためのプログラムコードのことである。なお、 引数は「いんすう」と発音するのではなく、「ひきすう」と発音するのが正しい。 ¶ ³ 関数名 ( 引数 ) { 関数の中身=プログラムコード } µ ´5.2
簡単な関数を定義してみる
まずは実に簡単な関数を定義してみよう。定義する関数は”pow”という名前の関数で、 引数にある1つの値を指定するとその値を自乗したものを返すというものである。だから 例えば、pow(2)としたら、22 = 4という値が返ってくるわけである。 ¶ ³ #関数pow()のプログラムコード pow <- function(x){ res <- x ^ 2 res } #実行例 > pow(2)#2の自乗 [1] 4 > pow(12) #12の自乗 [1] 144#関数pow2のプログラムコード pow2 <- function(x, y){
res <- x ^ y res } #実行例 > pow2(2, 4) #2の4乗 [1] 16 > pow2(2, -4) #2の-4乗 [1] 0.0625 µ ´ また、関数を定義する際にはRのコンソールウィンドウに対して直接プログラムを書 くのではなく、Windowsならメモ帳(Rでも[ファイル][新しいスクリプト]でメモ張のよ うなものが現れる)などのテキストエディタに書いてから、それをコピー・アンド・ペー ストしたほうが効率的である。これは複雑なグラフを作成するときなのにも同じことがい える。 そして Rに最初から実装されている関数を利用してプログラミングができるというの も、Rの優れているところである。例えば、sd()は不偏分散を求める関数であるが、不偏 分散ではない通常の分散を求めたい場合はこのsd()を利用して次のようにすればよい。 ¶ ³ sd2 <- function(x){ var2 <- function(x){ #引数xがベクトルならベクトルの長さを取得する if(is.vector(x)) n <- length(x) #引数xが行列なら行数を取得する if(is.matrix(x)) n <- nrow(x) var2 <- var(x) * (n-1) / n #通常の分散を求める }
これは1つの関数の中に、さらに新たな関数を定義している少し特殊といえば特殊な ケースであるが、別にプログラムコードの意味を理解する必要はない。もちろん、理解し ようと努めれば好ましいのだが、ここではRの便利さを主張するためにプログラムコー ドを提示した。もしこれをRに実装されている関数を一切、使わずに定義しようとする とけっこうな手間になるのである。
5.3
関数のプロット
Rが優れているもう1つの点として関数を描くことができるということが挙げられる。 例えば、中学・高校で習ったような1次関数や2次関数、3次関数のグラフが簡単に描け る。Excelではいちいちセルに式を書き込んで、それをオートフィルさせたりなどする必 要があるが、Rでは上のように簡単な関数を定義しておいて、それをplot()の引数に指定 するだけでよい。 ¶ ³ #1次関数を定義 one.func <- function(x){ 2*x + 1 } #上で定義した関数をプロットする > plot(one.func) µ ´ もう1つの例を挙げておこう。以下の3次式はかなり複雑でパッと見ただけでは、どの ような曲線が描かれるのか想像もつかないが、このような複雑な式でも簡単にプロットす ることができるのである。ちなみに、関数の定義の際に、関数の中身が1行で終わる場合 は””でくくる必要はない。 ¶ ³ > thr.func <- function(x) 3*x^3 - 0.5*x^2 + sqrt(x) + (3/4) > plot(thr.func) µ ´ 直線や曲線を描く方法はplot()を用いる以外にもある。それはpoints()を利用すると> plot(thr.func, xlim=c(-2,2), ylim=c(-5,5)) #thr.funcをplot()で描 く
> abline(h=0) #原点を通るx軸を描く > abline(v=0) #原点を通るy軸を描く
> x.dat <- seq(-2, 2, 0.1) #one.funcで定義した直線のx軸の値
> points(x.dat, 2*x.dat+1, type="l", col="pink") #points() で 直 線 を 描く
> legend(-1.5, 4,
+ c("one.func", "thr.func"), lwd=2, + col=c("black", "pink")) #凡例の追加