Rによる統計分析ガイド(Windows対応版)
【概要】
Rでは、記述統計量の計算、クロス集計、各種検定、グラフ描画を行える。R標準のパッケージ(base/stats/grDevices)による基本処理から、momentsパッケージによる分布解析まで対応する。データはデータフレーム(行と列からなる2次元の表)を基盤として扱う。
【本資料の前提と準備】
本資料はWindowsパソコンにインストールしたRで動作することを前提とする。利用前に次を準備する。
- R本体のインストール:CRAN(Rの公式配布サイト)からWindows版のRをインストールする。
- パッケージの導入:歪度・尖度の計算でmomentsパッケージを使用する。初回のみ
install.packages("moments")を実行する。それ以外の処理はR標準の関数で動作する。 - 日本語の表示:グラフの日本語ラベルが文字化けしないよう、
windowsFonts()でWindows標準の日本語フォントを登録し、グラフ描画時にpar(family = "JP")で使用する(「グラフによる可視化」の節を参照)。
本資料のコード例は、各節を上から順に実行することを想定している。節単独で実行する場合は、その節で使うデータ作成コードも実行する。
統計的検定は次の流れで進める。本資料の各節はこの流れのいずれかの段階に対応する。
① データの確認・記述統計 → ② 正規性の確認 → ③ 手法の選択(正規性を棄却しない=パラメトリック/正規性を棄却する=ノンパラメトリック)→ ④ 検定の実行 → ⑤ 結果(p値)の解釈
【目次】
記述統計量の概要
記述統計量は、データの特徴を数値で要約する指標である。データ全体の特性を把握し、分析手法の選択に用いる。上記の流れの①に対応する。
主な記述統計量には、平均値、標準偏差、中央値、四分位数、最大値、最小値、分散、歪度、尖度がある。Rのsd関数とvar関数は、n − 1 で割る不偏推定に基づく標本標準偏差と標本分散を計算する。各語の定義は「用語説明」を参照する。
SPSS・Python・Rの処理対応表
次の表は、SPSS、Python、Rにおける主要な統計処理の対応を示す。
| 処理内容 | SPSS | Python (pandas/scipy) | R |
|---|---|---|---|
| 記述統計量 | FREQUENCIES |
df.describe()、stats.skew()、stats.kurtosis() |
summary()、sd()、var()、skewness()、kurtosis() |
| 頻度表 | FREQUENCIES |
value_counts() |
table() |
| クロス集計表 | CROSSTABS |
pd.crosstab() |
table() |
| 集約 | AGGREGATE |
groupby().agg() |
aggregate() |
| Welchのt検定 | T-TEST |
stats.ttest_ind(equal_var=False) |
t.test(var.equal = FALSE) |
| 一元配置分散分析 | ONEWAY |
stats.f_oneway() |
oneway.test() |
| Mann-Whitneyの検定 | NPAR TESTS /M-W= |
stats.mannwhitneyu() |
wilcox.test() |
| Kruskal-Wallis検定 | NPAR TESTS /K-W= |
stats.kruskal() |
kruskal.test() |
Rによる分析の実装例
次のデータセットを使用する。
| 科目 | 受講者 | 得点 |
|---|---|---|
| 国語 | A | 90 |
| 国語 | B | 80 |
| 算数 | A | 95 |
| 算数 | B | 90 |
| 理科 | A | 80 |
データフレームの作成
ベクトルでデータを定義し、data.frame関数でデータフレームに変換する。stringsAsFactors = FALSEを指定し、文字列を文字型として扱う。
# データの作成
科目 <- c("国語", "国語", "算数", "算数", "理科")
受講者 <- c("A", "B", "A", "B", "A")
得点 <- c(90, 80, 95, 90, 80)
# データフレームの作成
df <- data.frame(科目 = 科目, 受講者 = 受講者, 得点 = 得点,
stringsAsFactors = FALSE)
print(df)
個別の統計量の計算
R標準の関数で、基本統計量を個別に算出する。sd関数は標本標準偏差、var関数は標本分散を返す。
scores <- df$得点
print(paste("平均値:", mean(scores)))
print(paste("標準偏差:", sd(scores)))
print(paste("分散:", var(scores)))
print(paste("中央値:", median(scores)))
print(paste("最大値:", max(scores)))
print(paste("最小値:", min(scores)))
print(paste("第1四分位数:", unname(quantile(scores, 0.25))))
print(paste("第3四分位数:", unname(quantile(scores, 0.75))))
summary関数とカテゴリ別集計
summary関数は複数の統計量を一括で算出する。カテゴリ別の代表値を求める場合はaggregate関数に単一の関数(平均など)を指定すると、結果がデータフレームとして返る。複数の統計量をカテゴリ別に見る場合はtapply関数を使う。
# 全体の統計量
print(summary(df$得点))
# 科目別の平均得点(結果はデータフレーム)
print(aggregate(得点 ~ 科目, data = df, FUN = mean))
# 科目別に複数の統計量を見る場合
print(tapply(df$得点, df$科目, summary))
グラフによる可視化
R標準のグラフィック機能でデータの分布を可視化する。箱ひげ図は分布の概要を、ヒストグラムは頻度分布を示す。png関数で画像ファイルへの出力を開始し、dev.off関数で出力を終了する。
Windowsで日本語ラベルが文字化けしないようにするには、windowsFonts関数でWindows標準の日本語フォント(例:MS Gothic)を登録する。画像デバイスを開いた後、par(family = "JP")でそのフォントを使用する。
# Windows標準の日本語フォントを登録(名前 JP として登録)
windowsFonts(JP = windowsFont("MS Gothic"))
# 箱ひげ図
png("score_distribution.png", width = 800, height = 600, type = "windows")
par(family = "JP")
boxplot(得点 ~ 科目, data = df,
main = "科目別得点分布",
xlab = "科目", ylab = "得点")
dev.off()
# ヒストグラム
png("score_histogram.png", width = 800, height = 600, type = "windows")
par(family = "JP")
hist(df$得点, breaks = 5,
main = "得点分布",
xlab = "得点", ylab = "頻度")
dev.off()
歪度と尖度
momentsパッケージで分布の形状を計算する。歪度が0に近いほど左右対称である。尖度は、正規分布と比べた分布の裾の重さや尖り具合を表す指標である。パッケージ未導入の場合はinstall.packages("moments")でインストールする。
尖度の定義はパッケージにより異なる。momentsパッケージのkurtosisは正規分布で値が3になる定義である。一方、e1071パッケージなど過剰尖度(kurtosis − 3)を返す実装では正規分布で0になり、値が3だけ小さくなる。どちらの定義かを確認して比較する。
# 初回のみ: install.packages("moments")
library(moments)
# 乱数シードの設定(同じ値を指定すると同じ乱数列を再現できる)
set.seed(42)
# 正規分布に従う乱数を10000個生成
normal_data <- rnorm(10000)
print(paste("歪度:", round(skewness(normal_data), 3)))
print(paste("尖度:", round(kurtosis(normal_data), 3)))
クロス集計表
table関数で2つのカテゴリ変数の関係を集計する。組み合わせごとの度数を表形式で示す。ここでは、前節までのdfを上書きしないよう、クロス集計用のデータフレーム名をcross_dfとする。
グループ1 <- c("a", "b", "c", "a", "b")
グループ2 <- c("d", "d", "e", "e", "e")
cross_df <- data.frame(グループ1 = グループ1, グループ2 = グループ2,
stringsAsFactors = FALSE)
print(table(cross_df$グループ1, cross_df$グループ2))
検定手法の選択
統計的検定では、データの特性に応じて手法を選択する(流れの③)。まず各群の正規性を確認し、その結果に応じて次のように選ぶ。
- 2群の比較:各群で正規性を棄却しない場合はWelchのt検定(
t.test)、正規性を棄却する場合はMann-Whitneyの検定(wilcox.test)を使用する。 - 3群以上の比較:各群で正規性を棄却しない場合は一元配置分散分析(
oneway.test)、正規性を棄却する場合はKruskal-Wallis検定(kruskal.test)を使用する。
ここで扱う例は、互いに独立した群の比較である。対応のあるデータ(同じ人を前後で測る場合など)では、対応のある検定を選ぶ必要がある。
正規性の検定
shapiro.test関数でデータの正規性を検定する(流れの②)。p値が有意水準(通常0.05)以上の場合、正規性を棄却しない。この結果に基づき、続く節のパラメトリック検定(t検定・分散分析)またはノンパラメトリック検定(Mann-Whitney検定・Kruskal-Wallis検定)を選択する。
使用上の注意:shapiro.testは標本サイズが3〜5000の範囲でのみ実行できる。大標本では、実用上問題のない小さな逸脱でも正規性なしと判定されやすい。検定結果だけに頼らず、ヒストグラムやQ-Qプロット(理論分布と比べる散布図。qqnorm関数で描画)も併用する。
# 乱数シードの設定(同じ値を指定すると同じ乱数列を再現できる)
set.seed(42)
# 正規分布に従う乱数を100個生成
normal_sample <- rnorm(100)
result <- shapiro.test(normal_sample)
print(paste("検定統計量:", round(result$statistic, 3)))
print(paste("p値:", format.pval(result$p.value, digits = 3, eps = 0.001)))
# Q-Qプロット
windowsFonts(JP = windowsFont("MS Gothic"))
png("normal_qqplot.png", width = 800, height = 600, type = "windows")
par(family = "JP")
qqnorm(normal_sample, main = "Q-Qプロット")
qqline(normal_sample, col = "red")
dev.off()
# 解釈:
# p値が0.05以上なら、正規性を棄却しない。
# p値が0.05未満なら、正規性を棄却する。
Welchのt検定(2群・正規性を棄却しない場合)
t.test関数で2群間の平均値の差を検定する。Welchのt検定は等分散を仮定しないため、2群の分散が異なる場合にも適用できる。t.test関数では、var.equal = FALSEを指定するとWelchのt検定を実行する。
# 乱数シードの設定(同じ値を指定すると同じ乱数列を再現できる)
set.seed(42)
# 正規分布に従う2群のデータを生成
group1 <- rnorm(100, mean = 50, sd = 10)
group2 <- rnorm(100, mean = 55, sd = 15)
test_result <- t.test(group1, group2,
paired = FALSE,
var.equal = FALSE)
print(paste("t値:", round(test_result$statistic, 3)))
print(paste("p値:", format.pval(test_result$p.value, digits = 3, eps = 0.001)))
# 解釈:
# p値が0.05未満なら、2群の平均値に有意差があると判断する(帰無仮説を棄却)。
Mann-Whitneyの検定(2群・正規性を棄却する場合)
2群の正規性を棄却する場合は、wilcox.test関数でMann-Whitneyの検定(2標本のWilcoxon順位和検定)を行う。母集団の正規性を仮定せず、順位に基づいて2群の差を評価する。2群の分布形状が大きく異ならない場合は、中心位置の差として解釈できる。
本例では、正規分布ではない指数分布のデータを使う。また、exact = FALSEを指定し、正規近似によるp値を計算する。標本サイズが大きい場合や同順位を含むデータでも、同じ計算方式で実行できる。
# 乱数シードの設定(同じ値を指定すると同じ乱数列を再現できる)
set.seed(42)
# 正規分布ではない2群のデータを生成
group1 <- rexp(60, rate = 1.0)
group2 <- rexp(60, rate = 0.5)
test_result <- wilcox.test(group1, group2,
paired = FALSE,
exact = FALSE)
print(paste("W値:", round(test_result$statistic, 3)))
print(paste("p値:", format.pval(test_result$p.value, digits = 3, eps = 0.001)))
# 解釈:
# p値が0.05未満なら、2群の順位分布または中心位置に有意差があると判断する。
一元配置分散分析(3群以上・正規性を棄却しない場合)
oneway.test関数で3群以上の平均値の差を検定する。var.equal = FALSEを指定すると等分散を仮定しないWelch型の一元配置分散分析を実行し、var.equal = TRUEを指定すると等分散を仮定する通常の一元配置分散分析を実行する。有意差が認められた場合、事後検定(どの群間に差があるかを特定する検定。多重比較ともいう)を行う。
# 乱数シードの設定(同じ値を指定すると同じ乱数列を再現できる)
set.seed(42)
# 正規分布に従う3群のデータを生成
group_a <- rnorm(30, mean = 50, sd = 10)
group_b <- rnorm(30, mean = 55, sd = 12)
group_c <- rnorm(30, mean = 60, sd = 9)
# 値とグループ名を1つのデータフレームにまとめる
# factor() はグループ名をカテゴリ変数として扱うための関数
anova_df <- data.frame(
値 = c(group_a, group_b, group_c),
グループ = factor(c(rep("A", 30), rep("B", 30), rep("C", 30)))
)
result <- oneway.test(値 ~ グループ,
data = anova_df,
var.equal = FALSE)
print(paste("F値:", round(result$statistic, 3)))
print(paste("p値:", format.pval(result$p.value, digits = 3, eps = 0.001)))
# 解釈:
# p値が0.05未満なら、少なくとも1組の群間に平均値の有意差があると判断する。
# どの群間かは事後検定で特定する。
Kruskal-Wallis検定(3群以上・正規性を棄却する場合)
3群以上の正規性を棄却する場合は、kruskal.test関数でKruskal-Wallis検定を行う。一元配置分散分析のノンパラメトリック版であり、順位に基づいて3群以上の差を評価する。有意差が認められた場合、ノンパラメトリックの事後検定(例:pairwise.wilcox.test関数)でどの群間に差があるかを調べる。
本例では、正規分布ではない指数分布のデータを使う。
# 乱数シードの設定(同じ値を指定すると同じ乱数列を再現できる)
set.seed(42)
# 正規分布ではない3群のデータを生成
group_a <- rexp(30, rate = 1.0)
group_b <- rexp(30, rate = 0.7)
group_c <- rexp(30, rate = 0.4)
# 値とグループ名を1つのデータフレームにまとめる
# factor() はグループ名をカテゴリ変数として扱うための関数
kw_df <- data.frame(
値 = c(group_a, group_b, group_c),
グループ = factor(c(rep("A", 30), rep("B", 30), rep("C", 30)))
)
result <- kruskal.test(値 ~ グループ, data = kw_df)
print(paste("カイ二乗値:", round(result$statistic, 3)))
print(paste("p値:", format.pval(result$p.value, digits = 3, eps = 0.001)))
# 解釈:
# p値が0.05未満なら、少なくとも1組の群間に順位分布の有意差があると判断する。
用語説明
- 記述統計量:データの特性を数値で要約したもの。平均値、標準偏差、中央値などを含む。
- 平均値:データの総和をデータ数で割った値。分布の中心傾向を示す。外れ値の影響を受けやすい。
- 標準偏差:データのばらつきを示す指標。Rの
sd関数は n − 1 で割る標本標準偏差を返す。 - 中央値:順序付けたデータの中央に位置する値。外れ値の影響を受けにくい。
- 四分位数:データを4等分する境界値。第1四分位数、中央値(第2四分位数)、第3四分位数からなる。
- 最大値:データ内の最も大きい値。
- 最小値:データ内の最も小さい値。
- 分散:各データ点と平均値の差の二乗に基づくばらつきの指標。Rの
var関数は n − 1 で割る標本分散を返す。標本標準偏差は標本分散の平方根である。 - 歪度:分布の非対称性を示す指標。正規分布では0となる。
- 尖度:分布の裾の重さや尖り具合を示す指標。正規分布では3となる定義と、過剰尖度として0となる定義がある。算出式はパッケージにより異なる。
- データフレーム:行と列からなる2次元のデータ構造。
- 外れ値:他のデータから離れた値。平均値や分散などの統計量に影響する。
- ヒストグラム:データの分布を表すグラフ。各区間の頻度を棒グラフで示す。
- 箱ひげ図:データの分布を四分位数と外れ値で表すグラフ。中央値、四分位範囲、外れ値を示す。
- Q-Qプロット:データの分位点を理論分布(正規分布など)の分位点と対比して描く散布図。点が直線上に並ぶほど理論分布に近い。Rでは
qqnorm・qqline関数で描画する。 - クロス集計表:2つのカテゴリ変数の関係を表形式で示したもの。組み合わせごとの度数を示す。
- t検定:2群の平均値の差の統計的有意性を評価する検定手法。
- Welchのt検定:等分散を仮定しないt検定。分散が異なる2群の比較に適用する。Rでは
t.test関数にvar.equal = FALSEを指定して実行する。 - Mann-Whitneyの検定(Wilcoxonの順位和検定):2群の差を順位に基づいて検定するノンパラメトリック手法。正規性を仮定しない。Rの
wilcox.test関数は、2標本ではMann-Whitney検定、1標本・対応ありではWilcoxon符号順位検定として動作する。 - ノンパラメトリック検定:母集団分布の形状を正規分布などに限定しない検定手法。順位や符号を用いる。
- パラメトリック検定:母集団分布に特定の形状(正規分布など)を仮定する検定手法。
- 一元配置分散分析:3群以上の平均値の差を同時に検定する手法。Rの
oneway.test関数では、var.equal = FALSEによりWelch型、var.equal = TRUEにより等分散を仮定する方法を実行できる。 - Kruskal-Wallis検定:3群以上の差を順位に基づいて検定するノンパラメトリック手法。一元配置分散分析のノンパラメトリック版である。Rでは
kruskal.test関数で実行する。 - 事後検定:分散分析などで有意差が認められた後、どの群間に差があるかを特定する検定。多重比較ともいう。
- Shapiro-Wilk検定:データの正規性を評価する検定手法。Rでは
shapiro.test関数で実行する。標本サイズは3〜5000に限られる。 - p値:帰無仮説が真の場合に、観測データ以上に極端なデータが得られる確率。有意水準(例:0.05)と比較して判断する。非常に小さいp値を表示するときは、0と丸めず、
<0.001のように表示する。 - 帰無仮説:「差がない」「効果がない」とする仮説。p値が有意水準を下回るとき、この仮説を棄却する。
- 乱数シード:乱数生成の初期値。同じ値を指定すると同じ乱数列を再現できる。
- windowsFonts:Windows環境でグラフに使用するフォントを登録するR標準関数。日本語表示には「MS Gothic」などを指定する。