ブレークダウン Last modified: Dec 12, 2012
目的
複数のカテゴリー変数に基づいて多元分類を行い，各群のデータ数，平均値，標準偏差を求める。
さらに，一元分類の場合には一元配置分散分析またはクラスカル・ウォリス検定を行うこともできる。
使用法
breakdown(i, j, df, latex=TRUE, captions=NULL, labels=NULL,
test=c("none", "parametric", "non-parametric"),
var.equal=FALSE, digits=3, output="")
引数
i 分析対象となる変数が入っている列の番号または変数名ベクトル
c 関数を使った複数指定可能
j 群を表す変数 factor の列番号または変数名ベクトル（一元分類を繰り返すのではなく，多元分類を行う）
c 関数を使った複数指定可能
df 読み込んだデータフレームの名前
latex LaTeX ソースを出力する。Word 用なら，latex=FALSE にする
captions latex=TRUE の場合，各表の表題を文字列ベクトルとして指定することができる（デフォルトではあり合わせの表題を付ける）。
labels latex=TRUE の場合，各表の label を文字列ベクトルとして指定することができる（デフォルトでは付けない）。
test 一元分類の場合には検定法を指定することができる。"parametric", "non-parametric" の何れか
"parametric" の場合，2 群の場合には t.test 関数，3 群の場合には oneway.test 関数が呼ばれる
"non-parametric" の場合，2 群の場合には wilcox.test 関数，3 群の場合には kruskal.test 関数が呼ばれる
デフォルトは検定をしない
var.equal oneway.test(t.test) の var.equal 引数として使用する
digits 平均値，標準偏差を書き出すときに小数点以下何桁までにするかを指定
デフォルトでは 3
output 出力コネクション
デフルトはコンソールに出力
ソース
インストールは，以下の 1 行をコピーし，R コンソールにペーストする
source("http://aoki2.si.gunma-u.ac.jp/R/src/breakdown.R", encoding="euc-jp")
#####
#
# 複数のカテゴリー変数により多元分類を行い，各群の平均値，標準偏差を求め，必要なら平均値・代表値の差の検定を行う
#
#####
breakdown <- function(i, # 分析対象の変数が入っているデータフレーム上の列番号または変数名ベクトル
j, # 群を表す変数が入っているデータフレーム上の列番号または変数名ベクトル
df, # データフレーム
latex=TRUE, # LaTeX 形式で集計表を出力する（デフォルトは LaTeX 形式）
captions=NULL, # latex=TRUE のときに，各表の表題を表す文字列ベクトルを指定できる（NULL のときはデフォルトの表題）
labels=NULL, # latex=TRUE のときに，各表の label を表す文字列ベクトルを指定できる（NULL のときは付けない）
test=c("none", "parametric", "non-parametric"), # デフォルト none では検定を行わない。検定を行うときはその種類を指定する
statistics=c("mean", "median"), # （平均値，標準偏差）を求めるか（中央値，四分偏差）を求めるかを指定する
var.equal=FALSE, # t-test, oneway の場合に等分散性を仮定するかどうかの引数
digits=3, # 平均値，標準偏差を表示するときの小数点以下の桁数
output="", # ファイルに出力するときはファイル名（デフォルトはコンソールに出力）
encoding=getOption("encoding")) # ファイルに出力するときのエンコーディング（デフォルトは OS による）
{
getNum <- function(str, df) { # 変数名から列番号を得る
names <- colnames(df)
seq_along(names)[names %in% str]
}
SIQ <- function(x) return(diff(fivenum(x)[c(2,4)])) # 四分偏差を求める下請け関数
if (output != "") { # 結果をファイルに出力する場合の処理
output <- file(output, open="w", encoding=encoding)
}
test <- match.arg(test) # 引数が省略形で与えられたときに，正確な名前をとる
statistics <- match.arg(statistics) # 引数が省略形で与えられたときに，正確な名前をとる
if (statistics == "mean") {
MEAN <- mean # 位置の母数を求める関数：平均値
SD <- sd # 散布度を求める関数：標準偏差
M.str <- "平均値"
S.str <- "標準偏差"
}
else {
MEAN <- median # 位置の母数を求める関数：中央値
SD <- SIQ # 散布度を求める関数：四分偏差
M.str <- "中央値"
S.str <- "四分偏差"
}
format <- paste("%.", digits, "f", sep="") # 小数点以下 3 桁で出力する書式
if (is.character(i[1])) {
i <- getNum(i, df)
}
if (is.character(j[1])) {
j <- getNum(j, df)
}
index <- 0
for (k in i) { # ベクトルで指定されたすべての変数について分析する
ok <- complete.cases(df[,k], df[,j]) # 欠損値を持たないケースを特定
df2 <- df[ok,] # 欠損値を持つケースを除く
nl <- length(j) # 何元分類にwなるか
lst <- if (nl == 1) list(df2[,j])
else as.list(append(NULL, df2[,j]))
x <- df2[, k] # 分析対象変数
nt <- length(x) # 全体のデータ数
mt <- MEAN(x) # 全体の平均値
st <- SD(x) # 全体の標準偏差
nms <- as.matrix(cbind(aggregate(df2[,k], lst, length), # 各群のデータ数
aggregate(df2[,k], lst, MEAN)[,nl+1], # 各群の平均値
aggregate(df2[,k], lst, SD)[,nl+1])) # 各群の標準偏差を取り出し，行列形式にする
nr <- nrow(nms) # 行数（分類数）
nc <- ncol(nms) # 列数（分類詳細と，データ数，平均値，標準偏差）
str <- paste(colnames(df2)[j], collapse="，") # 多元分類に使われる変数名のリスト
if (latex) { # LaTeX 形式で集計結果を出力する
index <- index+1
cat("\n\\begin{table}[htbp]\n", file=output) # \begin{table}[htbp]
if (is.null(captions)) {
cat(sprintf("\\caption{%s別の%sの集計}\n", # \caption{○○別の□□の集計}
str, colnames(df2)[k]), file=output)
}
else {
cat(sprintf("\\caption{%s}\n", captions[index]), file=output) # \caption{○○○○}
}
if (!is.null(labels)) {
cat(sprintf("\\label{%s}\n", labels[index]), file=output) # \labels{○○○○}
}
cat("\\centering\n", file=output) # \centering
cat("\\begin{tabular}{", rep("l", nc-3), # \begin{tabular}{l…ccc} \hline
"ccc} \\hline\n", sep="", file=output)
cat(rep("&", nc-3), # 分類に使用する変数分の &
sprintf(" \\multicolumn{3}{c}{%s}\\\\ \\cline{%i-%i}\n", # 最後の3列を使って分析対象の変数名と罫線
colnames(df2)[k], nc-2, nc), file=output)
cat(colnames(df2)[j], "データ数", M.str, S.str, # 分類変数名 & データ数 & 平均値 & 標準偏差
sep=" & ", file=output)
cat("\\\\ \\hline\n", file=output) # \\ \hline
for (l in 1:nr) { # 各分類ごとに，
cat(nms[l,1:(nc-2)], # 分類基準，データ数，平均値，標準偏差
sprintf(format, as.numeric(nms[l, (nc-1):nc])), sep=" & ", file=output)
cat("\\\\", file=output) # \\
if (l == nr) cat("\\hline\n", file=output)
else cat("\n", file=output) # 最後の行なら \hline
}
cat(sprintf("\\multicolumn{%i}{l}{全体}", nc-3), # 全体について，データ数，平均値，標準偏差
nt, sprintf(format, mt), sprintf(format, st), sep=" & ", file=output)
cat("\\\\ \\hline\n", file=output) # \\ \hline
cat("\\end{tabular}\n", file=output) # \end{tabular}
}
else { # tab で区切って出力する
cat("\n表 ", str, "別の", colnames(df2)[k], "の集計", # 表 ○○別の□□の集計
sep="", file=output)
cat("\n", colnames(df2)[k], sep="\t", file=output, # 分析対象変数
fill=TRUE)
cat(colnames(df2)[j], "データ数", M.str, S.str, # 分類変数名 データ数 平均値 標準偏差
sep="\t", file=output, fill=TRUE)
for (l in 1:nr) { # 各分類ごとに，
cat(nms[l,1:(nc-2)], sprintf(format, # 分類基準，データ数，平均値，標準偏差
as.numeric(nms[l, (nc-1):nc])), sep="\t",
file=output, fill=TRUE)
}
cat("全体", nt, sprintf(format, mt), # 全体について，データ数，平均値，標準偏差
sprintf(format, st), sep="\t", file=output, fill=TRUE)
}
if (nr >= 2 && nl == 1) { # 一元分類のときのみ検定を行う
if (latex && test != "none") {
cat("\\\\ \\noindent\n", file=output)
}
if (nr == 2 && test == "parametric") { # 2 群の場合には t.test を呼ぶ
res <- t.test(x~df2[,j], var.equal=var.equal)
cat(sprintf(if (latex) "$t$値 = %.3f, 自由度 = %.3f, $P$値 = %.3f\n"
else "t 値 = %.3f, 自由度 = %.3f, P 値 = %.3f\n",
res$statistic, res$parameter, res$p.value), file=output)
}
else if (nr >= 3 && test == "parametric") { # 3 群以上の場合には oneway.test を呼ぶ
res <- oneway.test(x~df2[,j], var.equal=var.equal)
cat(sprintf(if (latex) "$F$値 = %.3f, 自由度 = (%i, %.3f), $P$値 = %.3f\n"
else "F 値 = %.3f, 自由度 = (%i, %.3f), P 値 = %.3f\n",
res$statistic, res$parameter[1], res$parameter[2], res$p.value), file=output)
}
else if (nr == 2 && test == "non-parametric") { # 2 群以上の場合には wilcox.test を呼ぶ
res <- wilcox.test(x~df2[,j]) # マン・ホイットニーの U 検定
cat(sprintf(if (latex) "$U$ = %.3f, $P$値 = %.3f\n"
else "U = %.3f, P 値 = %.3f\n",
res$statistic, res$p.value), file=output)
}
else if (nr >= 3 && test == "non-parametric") { # 3 群以上の場合には kruskal.test を呼ぶ
res <- kruskal.test(x~df2[,j])
cat(sprintf(if (latex) "$\\chi^2_{kw}$値 = %.3f, 自由度 = %i, $P$値 = %.3f\n"
else "カイ二乗値(kw) = %.3f, 自由度 = %i, P 値 = %.3f\n",
res$statistic, res$parameter, res$p.value), file=output)
}
}
if (latex) { # LaTeX 形式で集計結果を出力する場合は，
cat("\\end{table}\n", file=output) # \end{table}
}
}
if (output != "") { # 結果をファイルに出力した場合の後始末
close(output)
}
}
使用例
男女 血液型 測定値1 測定値2
1 2 101.3 135.8
1 2 101.8 165.7
1 1 122.6 178.9
2 4 114.4 155.8
1 4 89.7 143.0
2 2 83.4 123.7
1 2 93.2 151.3
2 1 99.6 108.8
1 3 127.5 101.6
2 4 104.7 156.1
2 4 84.2 173.7
2 2 90.6 122.6
1 1 90.0 102.2
1 1 97.6 176.7
2 4 97.4 125.0
1 1 92.8 172.8
2 3 103.4 146.5
2 2 101.2 159.0
2 2 87.8 148.0
2 1 95.5 140.5
1 2 89.0 171.8
2 3 94.3 121.2
1 3 69.9 114.9
1 3 79.4 160.0
1 3 105.7 98.1
1 2 95.2 133.5
2 4 113.7 86.5
2 2 106.8 148.9
2 2 112.1 206.1
2 1 91.7 163.0
のようなファイル test.dat があるとする
入力と変数の定義
df <- read.table("test.dat", header=TRUE)
df[,1] <- factor(df[,1], levels=1:2, labels=c("男", "女"))
df[,2] <- factor(df[,2], levels=1:4, labels=c("A", "B", "O", "AB"))
ファイルに出力するとき
breakdown(3:4, 1, df, latex=FALSE, output="ファイル名")
出力結果例
latex=FALSE の場合
> breakdown(3:4, 1, df, latex=FALSE)
表 男女別の測定値1の集計
測定値1
男女 データ数 平均値 標準偏差
男 14 96.836 15.060
女 16 98.800 9.969
全体 30 97.883 12.413
表 男女別の測定値2の集計
測定値2
男女 データ数 平均値 標準偏差
男 14 143.307 29.535
女 16 142.838 28.217
全体 30 143.057 28.336
> breakdown(3:4, 2, df, latex=FALSE)
表 血液型別の測定値1の集計
測定値1
血液型 データ数 平均値 標準偏差
A 7 98.543 11.127
B 11 96.582 8.795
O 6 96.700 20.483
AB 6 100.683 12.465
全体 30 97.883 12.413
表 血液型別の測定値2の集計
測定値2
血液型 データ数 平均値 標準偏差
A 7 148.986 32.410
B 11 151.491 24.184
O 6 123.717 24.759
AB 6 140.017 30.792
全体 30 143.057 28.336
> breakdown(3:4, 1:2, df, latex=FALSE)
表 男女，血液型別の測定値1の集計
測定値1
男女 血液型 データ数 平均値 標準偏差
男 A 4 100.750 14.901
女 A 3 95.600 3.951
男 B 5 96.100 5.458
女 B 6 96.983 11.421
男 O 4 95.625 26.093
女 O 2 98.850 6.435
男 AB 1 89.700 NA
女 AB 5 102.880 12.570
全体 30 97.883 12.413
表 男女，血液型別の測定値2の集計
測定値2
男女 血液型 データ数 平均値 標準偏差
男 A 4 157.650 37.053
女 A 3 137.433 27.230
男 B 5 151.620 17.206
女 B 6 151.383 30.543
男 O 4 118.650 28.501
女 O 2 133.850 17.890
男 AB 1 143.000 NA
女 AB 5 139.420 34.388
全体 30 143.057 28.336
latex=TRUE（デフォルト）の場合
> breakdown(4, 1:2, df)
\begin{table}[htbp]
\caption{男女，血液型別の測定値2の集計}
\begin{center}
\begin{tabular}{lccc} \hline
& \multicolumn{3}{c}{測定値2}\\ \cline{2-4}
男女 & 血液型 & データ数 & 平均値 & 標準偏差\\ \hline
男 & A & 4 & 157.650 & 37.053\\
女 & A & 3 & 137.433 & 27.230\\
男 & B & 5 & 151.620 & 17.206\\
女 & B & 6 & 151.383 & 30.543\\
男 & O & 4 & 118.650 & 28.501\\
女 & O & 2 & 133.850 & 17.890\\
男 & AB & 1 & 143.000 & NA\\
女 & AB & 5 & 139.420 & 34.388\\\hline
全体 & 30 & 143.057 & 28.336\\ \hline
\end{tabular}
\end{center}
\end{table}
直前のページへ戻る E-mail to Shigenobu AOKI
