Rでサンプルサイズの計算　②仮説検定

3.  臨床試験における2群検定

新薬の血圧低下への効果を調べる臨床試験を設計したいとします。コントロール群と治療群の平均血圧の比較を行う際のサンプルサイズを設定する例を紹介します。
治療によって平均血圧が5 mmHg低下することを期待し、標準偏差は10 mmHgであると仮定し検定力を80%、有意水準を5%に設定し、必要なサンプルサイズを計算します。

平均血圧の差:　5 mmHg
標準偏差: 10 mmHg
検定力（1-β）: 80％
有意水準: 5%

上記のような場合でも、一般的に使用される”両側検定”での計算を進めます。
“平均血圧が低下することを期待し”、とありますが、臨床試験に至るまでに新薬によって血圧が低下することがある程度根拠がある場合でも、一般的には、低下する方向の検定、つまり片側検定を行うことは多くないようです。
低下が予測できたとしても、その根拠が科学的に、公に、認められていなければ、両側検定を行うことが推奨される場合が多いです。片側検定は、両側検定に比べて、有意でない差を有意と判断してしまうリスクが高まり、必要なサンプルサイズも小さく見積もられることがあるため、計画段階から十分な注意が必要です。

・帰無仮説 (H₀): 治療によって血圧に変化はない
・対立仮説 (H₁): 治療によって血圧に変化がある（低下または上昇）

# pwrライブラリのインストール・読み込み
# install.packages("pwr")
library(pwr)

# パラメータ設定
mean_diff <- 5     # 期待される平均の差（治療群と対照群の血圧差）
std_dev <- 10      # 標準偏差
alpha <- 0.05      # 有意水準
power <- 0.8       # 検出力

# 効果量（Cohen's d）の計算
effect_size <- mean_diff / std_dev

# サンプルサイズの計算 (両側検定)
result <- pwr.t.test(d = effect_size, sig.level = alpha, power = power, type = "two.sample")
# 結果の表示
print(result)
#     Two-sample t test power calculation 
#              n = 63.76561
#              d = 0.5
#      sig.level = 0.05
#          power = 0.8
#    alternative = two.sided

#NOTE: n is number in *each* group

5.　サンプルサイズを可視化ー検出力との関係ー

4.  製品質量の品質管理の条件で、異なるサンプルサイズに対する検出力を可視化するものです。
サンプルサイズが増えるにつれて、検出力がどのように変化するかがわかります。

library(pwr)
library(ggplot2)

# パラメータ設定
mean_diff <- 5     # 期待される質量の差（±5g）
std_dev <- 5       # 製品質量の標準偏差
alpha <- 0.05      # 有意水準
power <- 0.8       # 検出力

# 効果量（Cohen's d）の計算
effect_size <- mean_diff / std_dev

# 検出力分析のためのサンプルサイズのリスト
sample_sizes <- seq(5, 50, by = 1)  # 5から50までのサンプルサイズ

# 各サンプルサイズに対する検出力を計算
powers <- sapply(sample_sizes, function(n) {
  result <- pwr.t.test(n = n, d = effect_size, sig.level = alpha, type = "two.sample", alternative = "two.sided")
  return(result$power)
})

# データフレームに整理
data_power <- data.frame(Sample_Size = sample_sizes, Power = powers)

# サンプルサイズと検出力の関係を可視化
ggplot(data_power, aes(x = Sample_Size, y = Power)) +
  geom_line(color = "green", size = 1) +
  geom_point(color = "blue", size = 3) +
  theme_minimal() +
  labs(title = "Sample Size vs Power",
       x = "Sample Size",
       y = "Power") +
  theme(plot.title = element_text(hjust = 0.5)) +
  geom_hline(yintercept = 0.8, linetype = "dashed", color = "red") +
  annotate("text", x = 30, y = 0.82, label = "Power = 0.8", color = "red")