マクネマー検定は対応のある二値データの解析に用いられる統計手法ですが、複数の計算方法が存在し、それぞれ異なる統計的性質を持ちます。本レポートでは、以下の4つの手法を包括的に比較します:
mcnemar.test(correct=FALSE))mcnemar.test(correct=TRUE))2×2分割表:
| 2回目陽性 | 2回目陰性 | 合計 | |
|---|---|---|---|
| 1回目陽性 | a | b | a+b |
| 1回目陰性 | c | d | c+d |
| 合計 | a+c | b+d | n |
重要: マクネマー検定では不一致ペア (b, c) のみが統計量の計算に使用されます。
原理: カイ二乗分布による近似
統計量:
χ² = (b - c)² / (b + c)
特徴:
原理: Yatesの連続性補正付きカイ二乗近似
統計量:
χ² = max(0, (|b - c| - 1)²) / (b + c)
特徴:
原理: 二項分布による正確計算
p値:
p = 2 × P(Binomial(n = b+c, p = 0.5) ≤ min(b,c))
特徴:
原理: Exact testの保守性を緩和した手法
p値:
Mid-p = 2 × [P(X ≤ min(b,c)) - 0.5 × P(X = min(b,c))] where X ~ Binomial(n = b+c, p = 0.5)
特徴:
関連論文:
マクネマー検定のサンプルサイズ計算では、不一致ペア数 (b+c) が重要な要素となります。
# マクネマー検定のサンプルサイズ計算 (複数手法)
# 1. 近似計算
mcnemar_sample_size_approx <- function(alpha = 0.05, power = 0.80, p10, p01) {
z_alpha <- qnorm(1 - alpha/2)
z_beta <- qnorm(power)
effect_size <- abs(p10 - p01)
pooled_prob <- p10 + p01
n_discord <- ((z_alpha * sqrt(pooled_prob) + z_beta * sqrt(pooled_prob - effect_size^2/pooled_prob))^2) / effect_size^2
return(ceiling(n_discord))
}
# 2. 連続性補正付き
mcnemar_sample_size_corrected <- function(alpha = 0.05, power = 0.80, p10, p01) {
z_alpha <- qnorm(1 - alpha/2)
z_beta <- qnorm(power)
effect_size <- abs(p10 - p01)
pooled_prob <- p10 + p01
n_discord <- ((z_alpha * sqrt(pooled_prob) + z_beta * sqrt(pooled_prob - effect_size^2/pooled_prob))^2 + 2*z_alpha*sqrt(pooled_prob)) / effect_size^2
return(ceiling(n_discord))
}
# 3. 正確計算 (推奨)
mcnemar_sample_size_exact <- function(alpha = 0.05, power = 0.80, p10, p01) {
delta <- p10 - p01
p_discord <- p10 + p01
if (p_discord <= 0) stop("p10 + p01 は正の値である必要があります")
p0 <- 0.5
p1 <- p10 / p_discord
power_function <- function(n) {
critical_value <- qbinom(1 - alpha/2, n, p0)
power_calc <- 1 - pbinom(critical_value, n, p1) + pbinom(n - critical_value - 1, n, p1)
return(power_calc)
}
n <- 10
while (power_function(n) < power && n < 10000) {
n <- n + 1
}
return(n)
}
# 使用例
p10 <- 0.20 # 1回目陽性→2回目陰性
p01 <- 0.10 # 1回目陰性→2回目陽性
cat("サンプルサイズ計算結果 (不一致ペア数):\n")
cat("近似計算:", mcnemar_sample_size_approx(p10 = p10, p01 = p01), "\n")
cat("連続性補正付き:", mcnemar_sample_size_corrected(p10 = p10, p01 = p01), "\n")
cat("正確計算:", mcnemar_sample_size_exact(p10 = p10, p01 = p01), "\n")
# 不一致ペア数から総サンプルサイズの推定
estimate_total_sample_size <- function(discord_pairs, p10, p01, p11) {
discord_proportion <- p10 + p01
total_n <- ceiling(discord_pairs / discord_proportion)
return(total_n)
}
# 例: 不一致ペア数50の場合
p11 <- 0.30 # 両方陽性
p00 <- 1 - p10 - p01 - p11 # 両方陰性
total_n <- estimate_total_sample_size(50, p10, p01, p11)
cat("推定総サンプルサイズ:", total_n, "\n")
各手法のType I error制御と検出力を包括的に比較
# マクネマー検定: 各方法の性能比較シミュレーション
set.seed(123) # 再現性のため
# 1. Mid-p test の実装
midp_mcnemar <- function(b, c) {
if (b + c == 0) return(1)
n <- b + c
min_bc <- min(b, c)
p_midp <- 2 * (pbinom(min_bc, n, 0.5) - 0.5 * dbinom(min_bc, n, 0.5))
return(min(p_midp, 1))
}
# 2. Exact test の実装
exact_mcnemar <- function(b, c) {
if (b + c == 0) return(1)
min_bc <- min(b, c)
p_exact <- 2 * pbinom(min_bc, b + c, 0.5)
return(min(p_exact, 1))
}
# 3. データ生成関数
generate_mcnemar_data <- function(n, p10, p01, p11) {
p00 <- 1 - p10 - p01 - p11
if (p00 < 0) stop("確率の合計が1を超えています")
outcomes <- rmultinom(1, n, c(p11, p10, p01, p00))
return(list(a = outcomes[1], b = outcomes[2], c = outcomes[3], d = outcomes[4]))
}
# 4. 単一シミュレーション関数
run_single_simulation <- function(n_pairs, p10, p01, p11, n_sim = 5000) {
results <- data.frame(
r_standard = numeric(n_sim),
r_corrected = numeric(n_sim),
exact = numeric(n_sim),
midp = numeric(n_sim)
)
for (i in 1:n_sim) {
data <- generate_mcnemar_data(n_pairs, p10, p01, p11)
table_2x2 <- matrix(c(data$a, data$c, data$b, data$d), nrow = 2)
# R標準のmcnemar.testを使用 (エラーハンドリング付き)
tryCatch({
mcnemar_false <- mcnemar.test(table_2x2, correct = FALSE)
mcnemar_true <- mcnemar.test(table_2x2, correct = TRUE)
results$r_standard[i] <- mcnemar_false$p.value
results$r_corrected[i] <- mcnemar_true$p.value
}, error = function(e) {
results$r_standard[i] <<- NA
results$r_corrected[i] <<- NA
})
results$exact[i] <- exact_mcnemar(data$b, data$c)
results$midp[i] <- midp_mcnemar(data$b, data$c)
}
return(results)
}
# 5. Type I Error シミュレーション
type1_error_simulation <- function(n_pairs, p10_p01, p11, n_sim = 5000, alpha = 0.05) {
cat(sprintf("Type I Error シミュレーション (n=%d, p10=p01=%.2f, p11=%.2f)\n",
n_pairs, p10_p01, p11))
results <- run_single_simulation(n_pairs, p10_p01, p10_p01, p11, n_sim)
# NAを除いて計算
type1_rates <- sapply(results, function(x) {
valid_values <- x[!is.na(x)]
if(length(valid_values) == 0) return(NA)
mean(valid_values < alpha)
})
# NAの割合も計算
na_rates <- sapply(results, function(x) mean(is.na(x)))
result_df <- data.frame(
Method = c("R標準(correct=FALSE)", "R標準(correct=TRUE)", "Exact", "Mid-p"),
Type1_Error_Percent = round(type1_rates * 100, 2),
NA_Rate = round(na_rates * 100, 1),
Expected = rep(5.0, 4)
)
print(result_df)
return(list(type1_rates = type1_rates, na_rates = na_rates))
}
# 6. 検出力シミュレーション
power_simulation <- function(n_pairs, p10, p01, p11, n_sim = 5000, alpha = 0.05) {
effect_size <- p10 - p01
cat(sprintf("検出力シミュレーション (n=%d, p10=%.2f, p01=%.2f, 効果量=%.2f)\n",
n_pairs, p10, p01, effect_size))
results <- run_single_simulation(n_pairs, p10, p01, p11, n_sim)
# NAを除いて計算
power_rates <- sapply(results, function(x) {
valid_values <- x[!is.na(x)]
if(length(valid_values) == 0) return(NA)
mean(valid_values < alpha)
})
result_df <- data.frame(
Method = c("R標準(correct=FALSE)", "R標準(correct=TRUE)", "Exact", "Mid-p"),
Power_Percent = round(power_rates * 100, 1),
Rank = rank(-power_rates, na.last = "keep")
)
print(result_df)
return(power_rates)
}
# 7. メインシミュレーション実行
cat("マクネマー検定 各方法の性能比較シミュレーション\n")
cat(paste(rep("=", 60), collapse = ""), "\n\n")
# パート1: Type I Error シミュレーション
cat("【パート1: Type I Error シミュレーション】\n")
cat(paste(rep("-", 40), collapse = ""), "\n")
# 中確率パターン
cat("\n◆ 中確率パターン (p10=p01=0.15, p11=0.30)\n")
cat("小サンプル (n=30):\n")
type1_mid_small <- type1_error_simulation(30, 0.15, 0.30, 5000)
cat("\n大サンプル (n=200):\n")
type1_mid_large <- type1_error_simulation(200, 0.15, 0.30, 5000)
# パート2: 検出力シミュレーション
cat("\n\n", paste(rep("=", 60), collapse = ""), "\n")
cat("【パート2: 検出力シミュレーション】\n")
cat(paste(rep("-", 40), collapse = ""), "\n")
# 中効果
cat("\n◆ 中効果 (p10=0.20, p01=0.10, 効果量=0.10)\n")
power_medium <- power_simulation(100, 0.20, 0.10, 0.30, 5000)
# 大効果
cat("\n◆ 大効果 (p10=0.25, p01=0.10, 効果量=0.15)\n")
power_large <- power_simulation(100, 0.25, 0.10, 0.30, 5000)
cat("\nシミュレーション完了！\n")
マクネマー検定 各方法の性能比較シミュレーション
============================================================
【パート1: Type I Error シミュレーション】
----------------------------------------
◆ 中確率パターン (p10=p01=0.15, p11=0.30)
小サンプル (n=30):
Type I Error シミュレーション (n=30, p10=p01=0.15, p11=0.30)
Method Type1_Error_Percent NA_Rate Expected
r_standard R標準(correct=FALSE) 4.18 0 5
r_corrected R標準(correct=TRUE) 1.96 0 5
exact Exact 1.96 0 5
midp Mid-p 3.88 0 5
大サンプル (n=200):
Type I Error シミュレーション (n=200, p10=p01=0.15, p11=0.30)
Method Type1_Error_Percent NA_Rate Expected
r_standard R標準(correct=FALSE) 5.40 0 5
r_corrected R標準(correct=TRUE) 4.18 0 5
exact Exact 4.22 0 5
midp Mid-p 5.36 0 5
============================================================
【パート2: 検出力シミュレーション】
----------------------------------------
◆ 中効果 (p10=0.20, p01=0.10, 効果量=0.10)
検出力シミュレーション (n=100, p10=0.20, p01=0.10, 効果量=0.10)
Method Power_Percent Rank
r_standard R標準(correct=FALSE) 45.5 1
r_corrected R標準(correct=TRUE) 37.5 4
exact Exact 37.6 3
midp Mid-p 43.0 2
◆ 大効果 (p10=0.25, p01=0.10, 効果量=0.15)
検出力シミュレーション (n=100, p10=0.25, p01=0.10, 効果量=0.15)
Method Power_Percent Rank
r_standard R標準(correct=FALSE) 73.6 1
r_corrected R標準(correct=TRUE) 66.6 4
exact Exact 66.8 3
midp Mid-p 72.3 2
シミュレーション完了！
--------------------------------------------------------------
Analysis is conducted using R version 4.5.1 (2025-06-13)
The script uses the following packages and versions:
compiler 4.5.1
| 条件 | R標準(FALSE) | R標準(TRUE) | Exact | Mid-p |
|---|---|---|---|---|
| 小サンプル(n=30) | 4.18% | 1.96% | 1.96% | 3.88% |
| 大サンプル(n=200) | 5.40% | 4.18% | 4.22% | 5.36% |
Mid-p testが名目水準(5%)に最も近いことが確認できます。
| 効果量 | R標準(FALSE) | R標準(TRUE) | Exact | Mid-p |
|---|---|---|---|---|
| 中(0.10) | 45.5% | 37.5% | 37.6% | 43.0% |
| 大(0.15) | 73.6% | 66.6% | 66.8% | 72.3% |
一貫したパターン: R標準(correct=FALSE) > Mid-p > Exact ≈ R標準(correct=TRUE)
| 手法 | Type I Error制御 | 検出力 | 総合評価 |
|---|---|---|---|
| Mid-p test | ◎最適 | ◎高 | 🥇推奨 |
| R標準(FALSE) | ○良好 | ◎最高 | 🥈実用的 |
| R標準(TRUE) | △保守的 | △低 | ❌推奨しない |
| Exact test | △保守的 | △低 | ❌推奨しない |
# 推奨される判断フロー
discord_pairs <- table[1,2] + table[2,1]
if (discord_pairs < 25) {
# Mid-p推奨
result <- mcnemar_midp(table)
cat("Mid-p test を推奨 (小〜中サンプル)\n")
} else {
# どちらでも可 (Mid-pがより良い)
result_midp <- mcnemar_midp(table) # より良い
result_r <- mcnemar.test(table, correct=FALSE) # 実用的
cat("Mid-p test または R標準 どちらでも可\n")
}
| 研究タイプ | 不一致ペア<25 | 不一致ペア≥25 |
|---|---|---|
| 確認的臨床試験 | Mid-p | Mid-p |
| 探索的研究 | Mid-p推奨 | Mid-p/R標準 |
| 日常的解析 | Mid-p推奨 | R標準も可 |
# 実用的なMid-p McNemar test実装
mcnemar_midp <- function(table) {
# 入力検証
if (!is.matrix(table) || !all(dim(table) == c(2, 2))) {
stop("2x2の行列を入力してください")
}
b <- table[1, 2] # 1回目陽性→2回目陰性
c <- table[2, 1] # 1回目陰性→2回目陽性
n <- b + c
# 不一致ペアが0の場合
if (n == 0) {
return(list(
p.value = 1,
method = "Mid-p McNemar test",
statistic = c(b = b, c = c),
note = "不一致ペアが0のため、効果なしと判定"
))
}
# Mid-p値の計算
min_bc <- min(b, c)
p_midp <- 2 * (pbinom(min_bc, n, 0.5) - 0.5 * dbinom(min_bc, n, 0.5))
# 結果の返却
return(list(
p.value = min(p_midp, 1),
method = "Mid-p McNemar test",
statistic = c(b = b, c = c),
discord_pairs = n,
effect_direction = ifelse(b > c, "1→2で減少", ifelse(c > b, "1→2で増加", "変化なし"))
))
}
# 使用例
example_table <- matrix(c(20, 5, 15, 30), nrow = 2,
dimnames = list(
"治療前" = c("陽性", "陰性"),
"治療後" = c("陽性", "陰性")
))
result <- mcnemar_midp(example_table)
print(result)
# 包括的なマクネマー検定解析関数
comprehensive_mcnemar <- function(table, alpha = 0.05, method = "auto") {
# 基本情報の抽出
a <- table[1, 1]; b <- table[1, 2]
c <- table[2, 1]; d <- table[2, 2]
n_total <- sum(table)
discord_pairs <- b + c
cat("=== マクネマー検定 包括的解析 ===\n")
cat("分割表:\n")
print(table)
cat(sprintf("\n総ペア数: %d, 不一致ペア数: %d (%.1f%%)\n",
n_total, discord_pairs, discord_pairs/n_total*100))
# 方法の自動選択
if (method == "auto") {
if (discord_pairs < 25) {
selected_method <- "midp"
cat("自動選択: Mid-p test (推奨)\n")
} else {
selected_method <- "midp"
cat("自動選択: Mid-p test (最良の選択)\n")
}
} else {
selected_method <- method
}
# 複数方法での計算
results <- list()
# Mid-p test
results$midp <- mcnemar_midp(table)
# R標準
tryCatch({
r_result <- mcnemar.test(table, correct = FALSE)
results$r_standard <- list(
p.value = r_result$p.value,
method = "Standard McNemar test",
statistic = r_result$statistic
)
}, error = function(e) {
results$r_standard <- list(
p.value = NA,
method = "Standard McNemar test (failed)",
note = "計算エラー"
)
})
# 結果の表示
cat("\n=== 各方法の結果 ===\n")
methods_df <- data.frame(
Method = c("Mid-p", "R標準"),
P_value = c(
results$midp$p.value,
results$r_standard$p.value
),
Significant = c(
results$midp$p.value < alpha,
ifelse(is.na(results$r_standard$p.value), NA, results$r_standard$p.value < alpha)
)
)
print(methods_df)
# 効果量の計算
odds_ratio <- ifelse(c == 0, Inf, b/c)
# 治療前後の周辺確率
p1_before <- (a + b) / n_total
p1_after <- (a + c) / n_total
diff_proportion <- p1_after - p1_before
cat(sprintf("\n=== 効果量 ===\n"))
cat(sprintf("オッズ比: %.3f\n", odds_ratio))
cat(sprintf("確率差: %.3f\n", diff_proportion))
# 推奨事項
cat(sprintf("\n=== 推奨 ===\n"))
cat(sprintf("Mid-p test: p = %.4f %s\n",
results$midp$p.value,
ifelse(results$midp$p.value < alpha, "(有意)", "(非有意)")))
cat("Mid-p testが最も適切な方法です。\n")
return(invisible(results))
}
問題: b = c = 0の時の解釈
no_discord_table <- matrix(c(20, 0, 0, 30), nrow = 2)
解決策: 効果なしと判定
result <- mcnemar_midp(no_discord_table) # p.value = 1.0, 効果なしと判定
従来の指針: b + c < 25でExact testを使用
現在の推奨: 全ての場合でMid-p testを使用
# 新しい指針 (推奨) # 常にMid-p testを使用 result <- mcnemar_midp(table)
Fagerland, M. W., Lydersen, S., & Laake, P. (2013). The McNemar test for binary matched-pairs data: mid-p and asymptotic are better than exact conditional. BMC Medical Research Methodology, 13, 91.
McNemar, Q. (1947). Note on the sampling error of the difference between correlated proportions or percentages. Psychometrika, 12(2), 153-157.
Fagerland, M. W., Lydersen, S., & Laake, P. (2014). Recommended tests and confidence intervals for paired binomial proportions. Statistics in Medicine, 33(16), 2850-2875.
Wikipedia - McNemar's test
R Documentation - mcnemar.test
exact2x2 Package
# 最終的な推奨実装
final_mcnemar_test <- function(table, alpha = 0.05) {
# Mid-p testの実行
result <- mcnemar_midp(table)
# 結果の解釈
cat(sprintf("Mid-p McNemar test: p = %.4f %s\n",
result$p.value,
ifelse(result$p.value < alpha, "(有意)", "(非有意)")))
return(result)
}
# 実用例
example_table <- matrix(c(20, 5, 15, 30), nrow = 2)
final_mcnemar_test(example_table)
| 状況 | 推奨手法 | 理由 |
|---|---|---|
| 確認的研究 | Mid-p test | 最適なType I error制御 |
| 探索的研究 | Mid-p test または R標準 | どちらも実用可能 |
| 小サンプル | Mid-p test | より安定した性能 |
| 大サンプル | Mid-p test または R標準 | どちらも良好 |
| 日常的解析 | R標準も許容 | 実用上十分 |
統計解析にはMid-p McNemar testを使用した。有意水準は5%とした。 Mid-p testは従来のExact testの過度の保守性を改善し、Type I error 制御と検出力のバランスに優れることが知られている (Fagerland et al., 2013)。
統計解析にはMcNemar検定 (連続性補正なし) を使用した。有意水準は5%とした。 不一致ペア数は○○であり、漸近的近似が適用可能と判断した。
注意: 本レポートのコードは全て再現可能です。Rの基本パッケージのみを使用しているため、追加のパッケージインストールは不要です。
最終更新: 2025年8月
推奨引用: このレポートを引用する場合は、主要参考文献のFagerland et al. (2013)を併せて引用してください。