一度観測してからサンプルサイズを積むとp-hackingになるシミュレーション

一度観測してからサンプルサイズを積むとp-hackingになるシミュレーション

要旨

• twitterで話題になってたので自分でもやってみた
• pythonで実施
• マジで誤謬の確率が高くなる
  • よくよく考えると
    • 実験→検定→p値観測→分岐発生
    • [A]0.05より大きい→サンプルサイズ貯まるまで待つ
    • [B]0.05以下→有意差ありとしてここでストップ
  • 以下のような話が隠れている？
    • 観測0.08→待つ
    • →観測0.06→待つ
    • →観測0.04→本来ならここで止めるが続けると・・・
    • →観測0.06（閾値の範囲内に収まってしまった）
• 気になりごと
  • 今回のシミュレーションは帰無仮説が成立する場合の話なので棄却できるときはまた異なる話になるか？

本文

import numpy as np
from scipy import stats
import seaborn as sns

# https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.norm.html
# rvs(): 正規分布からランダムに値を取得する（Random variates）
# loc: 平均、scale: 標準偏差, size: 取得する数※,で次元を指定
# locとscaleを省略すると標準正規分布になる
# rvs = stats.norm.rvs(loc=0, scale=10, size=(50, 2))
rvs = stats.norm.rvs(size=(5, 3))
print(rvs)

rvs = stats.norm.rvs(size=(10000,2))
sns.displot(rvs)

def get_ttest_pvalue():
  rvs1 = stats.norm.rvs(loc=5, scale=3, size=500)
  rvs2 = stats.norm.rvs(loc=5, scale=3, size=500)
  stat, pvalue = stats.ttest_ind(rvs1,rvs2)
  return pvalue

pvalues = []
for i in range(10000):
  pvalues.append(get_ttest_pvalue())

print(len([x for x in pvalues if x <= 0.05]) / len(pvalues))
sns.displot(pvalues, bins=20)

• p-value ≤ 0.05の割合：0.0487

def get_ttest_pvalue_alternative():
  rvs1 = stats.norm.rvs(loc=5, scale=3, size=500)
  rvs2 = stats.norm.rvs(loc=5, scale=3, size=500)
  stat, pvalue = stats.ttest_ind(rvs1, rvs2)
  if (pvalue > 0.05):
    rvs1alt = stats.norm.rvs(loc=5, scale=3, size=500)
    rvs2alt = stats.norm.rvs(loc=5, scale=3, size=500)
    rvs1 = np.append(rvs1, rvs1alt)
    rvs2 = np.append(rvs2, rvs2alt)
    stat, pvalue = stats.ttest_ind(rvs1, rvs2)
  if (pvalue > 0.05):
    rvs1alt = stats.norm.rvs(loc=5, scale=3, size=500)
    rvs2alt = stats.norm.rvs(loc=5, scale=3, size=500)
    rvs1 = np.append(rvs1, rvs1alt)
    rvs2 = np.append(rvs2, rvs2alt)
    stat, pvalue = stats.ttest_ind(rvs1, rvs2)
  return pvalue

pvalues = []
for i in range(10000):
  pvalues.append(get_ttest_pvalue_alternative())
  
print(len([x for x in pvalues if x <= 0.05]) / len(pvalues))
sns.displot(pvalues, bins=20)

• p-value ≤ 0.05の割合：0.1115