Sponsored Links
対応のないt検定
データの読み込みと準備
サンプルデータは、統計教育推進委員会が提供してくださっている運動習慣と脈拍[1]を使います。
R Console
[1] 統計教育推進委員会 – 統計活用提案のための教材サイト
> x <- read.table(pipe("pbpaste"), header=TRUE) # クリップボードからデータを読み込み、xに格納。header=TRUEで先頭行を変数名と設定。ちなみにこちらはMac用のコードでWindowsは、read.table("clipboard", header=TRUE)。
> str(x) # xの構造を確認、intは整数を示すが、exerciseの0,1は本当は順序なし因子（0：運動習慣あり、1：運動習慣なし）。
'data.frame': 102 obs. of 2 variables:
$ HR : int 59 59 81 70 68 66 65 71 75 61 ...
$ exercise: int 0 0 0 0 1 1 1 1 0 0 ...
> x$exercise <- as.factor(x$exercise) # exerciseのデータ型を順序なし因子に変換し、上書き。
> str(x) # 再びxの構造を確認。exerciseがfactorに変換されている。
'data.frame': 102 obs. of 2 variables:
$ HR : int 59 59 81 70 68 66 65 71 75 61 ...
$ exercise: Factor w/ 2 levels "0","1": 1 1 1 1 2 2 2 2 1 1 ...
> summary(x) # xの要約を算出
HR exercise
Min. : 40.00 0:55
1st Qu.: 61.25 1:47
Median : 70.00
Mean : 70.21
3rd Qu.: 76.75
Max. :122.00
基本統計量の算出
運動習慣あり（0）、運動習慣なし（1）別に、平均値、標準偏差、標準誤差を求める。
R Console
> l0 <- nrow(subset(x, exercise=="0")) # 運動習慣ありの対象者数をl0に代入。
> l0
[1] 55
> l1 <- nrow(subset(x, exercise=="1")) # 運動習慣なしの対象者をl1に代入。
> l1
[1] 47
> m0 <- mean(subset(x, exercise=="0")$HR) # 運動習慣ありの心拍数の平均値。
> m0
[1] 69.63636
> m1 <- mean(subset(x, exercise=="1")$HR) # 運動習慣なしの心拍数の平均値。
> m1
[1] 70.87234
> d0 <- sd(subset(x, exercise=="0")$HR) # 運動習慣ありの心拍数の標準偏差。
> d0
[1] 9.972859
> d1 <- sd(subset(x, exercise=="1")$HR) # 運動習慣なしの心拍数の標準偏差。
> d1
[1] 13.94651
> e0 <- d0/sqrt(l0) # 運動習慣ありの心拍数の標準誤差。
> e0
[1] 1.34474
> e1 <- d1/sqrt(l1) # 運動習慣なしの心拍数の標準誤差。
> e1
[1] 2.034307
以下の方法でも同様に算出できます。こちらの方が簡便ですね。
R Console
> l <- table(x$exercise) # 2群の対象者数
> l
0 1
55 47
> m <- tapply(x$HR, x$exercise, mean) # 2群の平均値
> m
0 1
69.63636 70.87234
> d <- tapply(x$HR, x$exercise, sd) # 2群の標準偏差
> d
0 1
9.972859 13.946508
> e <- d/sqrt(l) # 2群の標準誤差
> e
0 1
1.344740 2.034307
等分散性の検定
等分散性の検定としてF検定の一種であるlevene検定を行います。levene検定は等分散であることを帰無仮説としますので、p値が0.05以上であれば帰無仮説は棄却されず、等分散であると言えます。levene検定にはcarパッケージを使います。
R Console
> library(car) # carを読み込む。
> leveneTest(x$HR, x$exercise) #levene検定の実施。
Levene's Test for Homogeneity of Variance (center = median)
Df F value Pr(>F)
group 1 1.0037 0.3188
100 p値が0.05より大きいので、運動習慣あり群と運動習慣なし群は等分散といえます。
等分散の場合はStudentのt検定、等分散でない場合はWelchのt検定
を行います。
MEMO
ちなみにStudentのt検定もWelchのt検定も正規分布を仮定した検定です。正規分布を仮定できないノンパラメトリックなデータにはMann-WhitneyのU検定を行うことになります。しかし、Mann-WhitneyのU検定も実は不等分散では上手く検定できないらしく、そういう時には、Brunner-Munzel検定という手法が有効だとのことです[1],[2]。 MEMO
ちなみにが多くて申し訳ありませんが、このF検定からのStudentもしくはWelch検定の手順については、二段階の多重検定であり、問題視されることも増えているようです[3]。
対応のないt検定の実行
t.test(変数~群分け水準)の式になります。 var.equal=を記述しなければデフォルトではFALSEになっており、Welchのt検定が行われます。
t.test(A群の変数, B群の変数, paired=FALSE)でも同様の結果になります。
R Console
> t.test(x$HR~x$exercise, var.equal=TRUE) # var.equal=TRUEで等分散を仮定したStudentのt検定を行うよう指示。
Two Sample t-test
data: x$HR by x$exercise
t = -0.52, df = 100, p-value = 0.6042
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-5.951673 3.479720
sample estimates:
mean in group 0 mean in group 1
69.63636 70.87234 今回は残念ながら運動経験あり群となし群の間で心拍数に有意差は認められませんでした。
効果量の計算
パッケージcompute.esを用いて計算します。
mes(平均1, 平均2, 標準偏差1, 標準偏差2, サンプルサイズ1, サンプルサイズ2)
で計算できます。
MEMO
効果量には数十種類もあるそうで、t検定では群間差についての効果量であるd族が用いられる。有名な効果量としてCohenのdやHedgesのgがあり、Cohenのdは平均値の差を標本標準偏差で割った値、Hedgesのgは平均値の差を不偏標準偏差で割った値である。従って、Cohenのdは記述統計的な意味合い、Hedgesのgは推測統計的な意味合いで使われる。推測統計であるt検定ではHedgesのgが推奨されている（Cohen’s d、Hedges’ g）[1],[2]。
R Console
[1] 統計Web – Cohen’s d
[2] 統計Web – Hedges’ g
> library(compute.es) #compute.esパッケージを読み込む。
> mes(m0, m1, d0, d1, l0, l1) # 上で定めた変数を使う。
Mean Differences ES:
d [ 95 %CI] = -0.1 [ -0.5 , 0.29 ]
var(d) = 0.04
p-value(d) = 0.6
U3(d) = 45.89 %
CLES(d) = 47.09 %
Cliff's Delta = -0.06
g [ 95 %CI] = -0.1 [ -0.49 , 0.29 ]
var(g) = 0.04
p-value(g) = 0.6
U3(g) = 45.92 %
CLES(g) = 47.11 %
Correlation ES:
r [ 95 %CI] = -0.05 [ -0.25 , 0.15 ]
var(r) = 0.01
p-value(r) = 0.61
z [ 95 %CI] = -0.05 [ -0.25 , 0.15 ]
var(z) = 0.01
p-value(z) = 0.61
Odds Ratio ES:
OR [ 95 %CI] = 0.83 [ 0.41 , 1.7 ]
p-value(OR) = 0.6
Log OR [ 95 %CI] = -0.19 [ -0.9 , 0.53 ]
var(lOR) = 0.13
p-value(Log OR) = 0.6
Other:
NNT = -36.17
Total N = 102
Sponsored Links
対応のあるt検定
データの読み込みと準備
サンプルデータは、統計教育推進委員会が提供してくださっているいろいろな樹木の葉の長さと幅のデータ[1]を使います。ダウンロードしたエクセルファイルにはいろいろな樹木のデータがシート別にまとめられています。今回はハナミヅキの葉の長さと幅に差があるかどうか比較することにします。葉の長さをleaf_l、葉の幅をleaf_wと名前を付けてデータを読み込みます。
R Console
[1] 統計教育推進委員会 – 統計活用提案のための教材サイト
> x <- read.table(pipe("pbpaste"), header=TRUE) # クリップボードからデータを読み込み、xに格納。header=TRUEで先頭行を変数名と設定。ちなみにこちらはMac用のコードでWindowsは、read.table("clipboard", header=TRUE)。
> str(x) # xの構造を確認。numは実数のこと。
'data.frame': 50 obs. of 2 variables:
$ leaf_l: num 8 7 7.5 7.5 7 7.5 6.5 7.5 6.5 5.5 ...
$ leaf_w: num 6.5 5 6.5 7 4 5.5 3.5 5 4.5 3 ...
> summary(x) # xの要約を算出。
leaf_l leaf_w
Min. :4.300 Min. :2.200
1st Qu.:5.625 1st Qu.:3.425
Median :6.500 Median :4.000
Mean :6.378 Mean :4.164
3rd Qu.:7.000 3rd Qu.:5.000
Max. :8.000 Max. :7.000
基本統計量の算出
leaf_lとleaf_wのサンプルサイズ、平均値、標準偏差、標準語さを計算します。
R Console
> attach(x) # 「x$」の部分を記述しなくて良いようにするおまじない。
> l <- c(length(leaf_l), length(leaf_w)) # length関数でleaf_lとleaf_wのサンプルサイズ（ベクトルの長さ）を求める。c()で1つのベクトルにくくる。このデータは１つの葉に対して長さと幅を測定しているので、当然ながらサンプルサイズは一致する。
> l
[1] 50 50
> m <- c(mean(leaf_l), mean(leaf_w)) # 2変数の平均値
> m
[1] 6.378 4.164
> d <- c(sd(leaf_l), sd(leaf_w)) # 2変数の標準偏差
> d
[1] 0.8572119 1.1385203
> e <- d/sqrt(l) # 2群の標準誤差
> e
[1] 0.1212281 0.1610111
対応のあるt検定の実行
R Console
> t.test(leaf_l, leaf_w, paired = TRUE) # paired=TRUEで対応のあるt検定になる。
Paired t-test
data: leaf_l and leaf_w
t = 16.091, df = 49, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
1.937505 2.490495
sample estimates:
mean of the differences
2.214 p値がとても小さいので有意差が認められました。平均値からハナミヅキ の葉は縦長だと分かります。
効果量の計算
t検定の効果量gの計算方法は、
g = t検定のt値 * sqrt( 1/(サンプルサイズ))
です。
R Console
> t <- t.test(leaf_l, leaf_w, paired = TRUE) # 上述の対応あるt検定をt定める。
> t$statistic # tに含まれるt値を確認。
t
16.09142
> g <- t$statistic * sqrt( 1/ length(leaf_l)) # gを計算。サンプルサイズはleaf_lとleaf_wは同数なので、ここではleaf_lを指定。
> g
t
2.27567
効果量は非常に大きいことがわかります。
[1] 図の学のWebPage 総合案内所 目次 梶山喜一郎 – 対応のあるt検定の効果量計算 [2] 水本篤 – 効果量と検定力分析入門 ―統計的検定を正しく使うために―（PDF）
Sponsored Links
