Rで統計分析する

2018-10-16 (火) 19:43:20 (9h) | Topic path: Top / バイオ・データ・マイニング / Rで統計分析する

はじめに

ここでは、Rを使って平均、分散、標準偏差を求めます。

『Rによるバイオインフォマティクスデータ解析』の3.3節「基本統計関数」に少しだけ出てきます。

準備

Rのインストールについては、次のページを見てください。

ここでは、標準で使用できるirisデータセットを使います。

> data(iris)

このデータセットは、アヤメの種類(Species)を花びらの長さ(Sepal.Length)、幅(Sepal.Width)、がくの長さ(Petal.Length)、幅(Petal.Width)によって分類する問題です。 長さと幅は連続値、種類はsetosa, versicolor, virginicaのいずれかをとる離散値です。

このデータセットには、setosa, versicolor, virginicaという3種類のアヤメについて、それぞれ50個ずつ、合計150個のデータが含まれています。 ランダムに10個のデータを選択して、見てみましょう。

> iris[sort(sample(1:150,10)),]
    Sepal.Length Sepal.Width Petal.Length Petal.Width    Species
4            4.6         3.1          1.5         0.2     setosa
22           5.1         3.7          1.5         0.4     setosa
65           5.6         2.9          3.6         1.3 versicolor
97           5.7         2.9          4.2         1.3 versicolor
100          5.7         2.8          4.1         1.3 versicolor
108          7.3         2.9          6.3         1.8  virginica
116          6.4         3.2          5.3         2.3  virginica
122          5.6         2.8          4.9         2.0  virginica
136          7.7         3.0          6.1         2.3  virginica
146          6.7         3.0          5.2         2.3  virginica

記述統計

標本を要約し、標本の情報をわかりやすく記述することを記述統計といいます。

ここでは、散布図、ヒストグラム、ボックスプロット(箱ひげ図)を作成します。

散布図

横軸と縦軸にそれぞれ別の量をとり、測定値を点として表したグラフを散布図といいます。 2つの値を同時に測定し、[math]n[/math] 個の測定値の組を [math](x_1, x_2), (x_2, y_2),\dots, (x_n, y_n)[/math] としたとき、それぞれを [math]x[/math] 座標、[math]y[/math] 座標とする点として表したものです。

散布図を作成するには、plot関数を用います。

> x = iris[,3]
> y = iris[,4]
> plot(x, y, main="Petal", xlab="Length [cm]", ylab="Width [cm]")
scatter.png

列ごとの散布図を作ることもできます。

> plot(iris[,1:4])
multiscatter.png

ヒストグラム

測定値が存在する範囲をいくつかの区間に分け、各区間とその区間に属する測定値の個数との関係性を度数分布と言います。 度数分布を表すグラフをヒストグラムといい、底辺の長さが各区間の幅に比例し、その面積がその区間の度数に比例する長方形を近接して並べたものです。

ヒストグラムを表示するにはhist関数を使います。

> x = iris[51:100,3]
> hist(x, main="Histogram of Petal.Length (Versicolor)", xlab="Petal.Length [cm]")
histogram.png

ボックスプロット(箱ひげ図)

ヒストグラムの他に、測定値の分布やばらつき具合を表すグラフとしてボックスプロット箱ひげ図)があります。 長方形の箱とその両端から伸びるひげで標本の統計量を表します。

箱の両端は、第一四分位数(最小値から全体の1/4のところにある測定値、25 パーセンタイル)と第三四分位数(最小値から全体の3/4のところにある測定値、75 パーセンタイル)を表し、箱の中の線は中央値(標本の大きさが奇数のときは全体の中央にある測定値、偶数のときは中央の2つの測定値の平均値、メディアン、50% パーセンタイル)を表します。

ひげの表し方には2種類あり、一つはひげの両端が最小値と最大値を表します。 もう一つは、ひげの両端が箱の両端から第三四分位数と第一四分位数の差の1.5倍の範囲内での最小値と最大値を表します。 後者の場合、ひげの両端よりも外側にある測定値を特異値(または外れ値)として丸印で表すこともあります。

> x = iris[1:50,3]
> y = iris[51:150,3]
> z = iris[101:150,3]
> boxplot(x, y, z, main="Petal", names=c("Setosa", "Versicolor", "Virginica"), ylab="Length [cm]")
boxplot.png

推定統計

標本から母集団の情報を推定することを推定統計といいます。

平均

母平均は、母集団の値の合計を母集団の大きさ [math]N[/math] で割ったものです。 \[ \mu(X) = \frac{1}{N} \sum_{i=1}^N x_i \]

標本平均は、標本の値の合計を標本の大きさ [math]n[/math] で割ったものです。 \[\overline{X} = \frac{1}{n} \sum_{i=1}^{n} x_i \]

標本平均を求めるには、mean関数を使います。

> mean(x)
[1] 1.462
> mean(y)
[1] 4.906
> mean(z)
[1] 5.552

標本平均の期待値は母平均に等しいため、標本平均を求めることで母平均を推定できます。

分散

測定値の標本平均からのズレを偏差といいます。

母分散は、母集団の値の偏差の二乗の総和を母集団の大きさ [math]N[/math] で割ったものです。 \[\sigma(X)^2 = \frac{1}{N} \sum_{i=1}^{n} (x_i - \mu(X))^2 \]

標本分散不偏分散実験分散)は、標本に含まれる値の偏差の二乗の総和を標本数から1を引いた数で割ったものです。 \[ s(X)^2 = \frac{1}{n-1} \sum_{i=1}^{n} (x_i - \mu(X))^2 \]

標本分散を求めるには、var関数を使います。

> var(x)
[1] 0.03015918
> var(y)
[1] 0.6815798
> var(z)
[1] 0.3045878

標本分散の期待値は母分散に等しいため、標本分散を求めることで母分散を推定できます。

標準偏差

母標準偏差は、母分散の平方根です。 \[\sigma(X) = \sqrt{\sigma(X)^2}\]

標本標準偏差不偏標準偏差実験標準偏差)は、不偏分散の平方根です。 \[ s = \sqrt{s^2} \]

標本標準偏差を求めるには、sd関数を使います。

> sd(x)
[1] 0.173664
> sd(y)
[1] 0.8255785
> sd(z)
[1] 0.5518947

標本標準偏差の期待値は母標準偏差に等しいため、標本標準偏差を求めることで母標準偏差を推定できます。

検定

平均値の差の検定(t検定)

setosaのSepal.Lengthの平均値は5.006でした。 versicolorのSepal.Lengthの平均値は5.936です。

この二つの平均値にの差には、意味があるのでしょうか。

これを、t検定という方法によって調べることができます。

t検定では、「二つの平均値の差には意味がある」という仮説に対して、「二つの平均値の差には意味がない」という反対の仮説(帰無仮説)を立て、帰無仮説が成り立たないことを示す(棄却する)ことによって、最初の仮説が成り立つことを示します。

t検定を行うには、t.test関数を使います。

> t.test(y, z)

	Welch Two Sample t-test

data:  y and z
t = -5.686, df = 135.83, p-value = 7.631e-08
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -0.870676 -0.421324
sample estimates:
mean of x mean of y 
    4.906     5.552 

これは、二つの群の分散が同じであることを仮定しない、Welchのt検定といいます。

二つの群の分散が同じであることを仮定する検定はStudentのt検定といい、オプションとして var.equal=T を付けます。
#geshi(rsplus){{
> t.test(y, z, var.equal=T)

	Two Sample t-test

data:  y and z
t = -4.9985, df = 148, p-value = 1.612e-06
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -0.9013938 -0.3906062
sample estimates:
mean of x mean of y 
    4.906     5.552 
 

まず、p値(p-value)を見ます。

p値が、0.1よりも小さいときは「10%有意水準で統計的有意な差がある」、0.05よりも小さいときは「5%有意水準で統計的有意な差がある」、0.01よりも小さいときは「1%有意水準で統計的有意な差がある」といいます。

つぎに、95%信頼区間(95 percent confidence interval)を見ます。 平均値の差が95%の確率でこの区間に入ることを意味します。

この区間がゼロを含むときは、5%有意水準で有意な差があるとはいえません。 「5%有意水準で有意な差がない」ではないことに注意しましょう。

演習

Sepal.LengthとSepal.Widthの関係を表す散布図、VersicolorのSepal.Lengthの分布を表すヒストグラム、SetosaとVersicolorとVirginicaのSepal.Lengthの分布を表すボックスプロットを作成しよう。

Setosa, Versicolor, VirginicaのSepal.Lengthについて、それぞれの標本平均 [math]\overline{X}[/math] と標本標準偏差 [math]s[/math] を求めてみよう。

また、VersicolorのSepal.Lengthの平均値とVirginicaのSepal.Lengthの平均値の間に統計的有意な差があるかどうかを調べてみよう。

参考文献

トップ   編集 凍結解除 差分 バックアップ 添付 複製 名前変更 リロード   新規 一覧 単語検索 最終更新   ヘルプ   最終更新のRSS