p値

数学いらずの医科統計学を踏まえて。
 
問1
新薬Xは臨床試験で有効(p=0.04)と判定された。
Xはすごくよく効く。
Xは効く。
Xは効かないこともない。
Xは効かない。
 
問2
新薬Xは臨床試験で有効(p=0.04)と判定された。
「Xが効く」ということが正しい確率は 0.96
「Xが効く」ということが正しい確率は 0.04
「Xが効かない」ということが正しい確率は 0.04
「Xが効かない」ということが正しい確率は 0.96
 
問3
新薬Xは臨床試験で有効(p=0.0004)と判定された。
Xはすごくよく効く。
Xは効く。
Xは効かないこともない。
Xは効かない。
 
問4
新薬Xは臨床試験で無効(p=0.06)と判定された。
Xはすごくよく効く。
Xは効く。
Xは効かないこともない。
Xは効かない。
 
p値の定義は、「帰無仮説が真のとき、対立仮説より大きな差を生じる確率」、つまり、新薬Xと既存薬が同じ効果だと帰無仮説を設定して、対立仮説ではそれより差が出る(0でない)とする。しかもその差は0.000…1でもいいわけだから、p値の大きさ(小ささ)は効果に相関しない。
こちら

(略)両群の統計計算上の差を表すp値が0.001未満でした。
このp値が低ければ低いほど、シロリムスが効くということになりますので、非常に有効であったということです。

と言ってしまっているのは、シロリムスがびっくりするほど効くという意味なら間違いで、この場合は「シロリムスがプロセボと同じと言ってしまうのはびっくりするほどない確率だが、よく効くかはどうなんだろうね」くらいな気がする。
というわけで問1,3,4は統計学的には「Xは効かないこともない。」としかいいようがない。問2はそもそもp値の定義が、帰無仮説・対立仮説の正当性を評価するものではないので答えなし。
 
問5
絶対に無効な(けど、我々にはわからない)新薬Yがある。これを臨床研究で調べたら
絶対に小さいp値が出る。
小さいp値が出やすい。
大きいp値が出やすい。
絶対に大きいp値が出る。
どんなp値が出るかはランダム。
 
問6
絶対に有効な(けど、我々にはわからない)新薬Yがある。これを臨床研究で調べたら
絶対に小さいp値が出る。
小さいp値が出やすい。
大きいp値が出やすい。
絶対に大きいp値が出る。
どんなp値が出るかはランダム。
 
帰無仮説ではたいてい、「2群の差はない」とする。なので、絶対に差がないとしても、対立仮説「2群の差はある」は、有意水準\alphaで生じる。2群に差がなくても、どんなp値が出るかはシミュレーション上、ランダムである。
実際には2群に差がないという状況はほぼありえないと思う。なので、シミュレーション上、差が1とか3とかある状況を設定したが、帰無仮説「2群に差がない」がそもそも破綻しているので小さいp値が出やすい。

n <- 10 # 1群のサンプル数
d0 <- c(0, 1, 3) # 真の差
niter <- 3000 # シミュレーション回数
res <- matrix(0, niter, length(d0))

for(j in seq(niter)){
	for(i in seq(d0)){
		a1 <- rnorm(n)
		a2 <- rnorm(n, d0[i])
		res[j, i] <- t.test(a1, a2)$p.value
	}
}
res <- apply(res, 2, sort)
matplot(res, pch=16, ylab="p value")
legend("topleft", legend=d0, lwd=3, col=seq(d0), bty="n")


 
こちらでは高々2%の違いを判断するためにどれほどのサンプルが必要かを論じている。2%というと大したことない差とおもわれるかもしれないが、それが重篤な副作用を起こすのであればそれは大変だし、しかもそれが50人の試験のときにわかるかというと、統計学的仮説検定では
 
50人で試験をしたら、2人に副作用が出た(2/50=0.04)
既存薬では、こういう場合 0.02 の確率で副作用が出ることが多い。
では、新薬の副作用頻度は既存薬に等しい 0.02 というのは、統計学的にはどうですか?
という流れになる。

binom.test(2, 50, 0.02)
        Exact binomial test

data:  2 and 50 
number of successes = 2, number of trials = 50, p-value = 0.2642
alternative hypothesis: true probability of success is not equal to 0.02 
95 percent confidence interval:
 0.004881433 0.137137626 
sample estimates:
probability of success 
                  0.04 

p値的には統計学的に有意ではない。
(本当はもっと考察があったのだろうが)50人規模では判断できないとして、数千人規模でやったら、くっきり分かれてしまった。
図を見れば、素人目に見ても違うことがわかる。

p <- c(0.02, 0.04)
n <- c(50, 3300)

p50 <- mapply(function(x) dbinom(0:50, 50, x), p)
p3300 <- mapply(function(x) dbinom(0:3300, 3300, x), p)

par(mfcol=c(2, 1))
matplot(p50, type="l", lwd=3, xlab="# side effect", ylab="probability")
legend("topright", legend=p, lwd=3, col=seq(p), bty="n", cex=2)

matplot(p3300, type="l", lwd=3, xlim=c(0, 250), xlab="# side effect", ylab="probability")
legend("topright", legend=p, lwd=3, col=seq(p), bty="n", cex=2)


 
差がない2群で検定したときのp値の分布や、確率 0.02 で50人もしくは3000人のうち副作用が出る人数の分布はどうか、という問で同級生に図を描いてもらったのだが、プロットのようにはいかなかった(撮り忘れた)。
数学力は(たぶん)日本でもかなり上位に入る人たち(だった)でも、感覚的にこういう事象をきちんと捉えることは難しいのだろうと思った。
 
バックナンバー