中野智文のブログ

データ・マエショリストのメモ

毎回ある確率未満であることを確認する必要なサンプルサイズ

背景

毎回ある確率未満であることを確認するのに、どれだけサンプルサイズが必要だろうか。 例えば、広告枠にて、0.05% 未満のクリック率の広告枠は興味がない場合、0.05% 未満であること自体は信頼上限の計算にて容易に確認できるが、毎日確率が変動する場合、ある日突然 0.05% 以上になっていたりすることがある。これはその広告枠に広告を出稿していなければわからない話だが、本来は興味が無いので、広告は出稿したくない。この場合最低どれほどサンプルを取得するべきなのだろうか。

導出

信頼上限の計算から逆算すれば、信頼上限がある確率未満になるのに必要なサンプルサイズが求められそうである。 Wilson score interval にて信頼上限を求めるとした場合を考える。 Nを全体の試行回数(求めたいサンプルサイズ)、 X を成功回数、Z正規分布の分位数、P をこれ未満にしたい信頼上限とする。このとき、Wilson socre interval の信頼上限は次のようになる。

P = (X+Z*Z/2+Z*sqrt(X*(N-X)/N+Z*Z/4))/(N+Z*Z)

これを maxima で解くことを考える。 以下やってみた。

(%i5) solve((X+Z*Z/2+Z*sqrt(X*(N-X)/N+Z*Z/4))/(N+Z*Z)=P,[N]);
                            2      2
                         N Z  - 4 X  + 4 N X         2    2
                  Z sqrt(-------------------) - 2 P Z  + Z  + 2 X
                                  N
(%o5)        [N = -----------------------------------------------]
                                        2 P

解こうとしたが、sqrt が外せないらしい。そこで危険ではあるが私が式を最小限変形して(過去に何度も間違えた経験あり)、解ける形にする。

{ \displaystyle
\begin{align}
\frac{X+\frac{Z^{2}}{2}+Z\cdot \sqrt{\frac{X\cdot (N-X)}{N}+\frac{Z^{2}}{4} }}{N+Z^{2}}  =  P
\end{align}
}

左辺の分母でかけて、

{ \displaystyle
\begin{align}
X+\frac{Z^{2}}{2}+Z\cdot\sqrt{\frac{X\cdot(N-X)}{N}+\frac{Z^{2}}{4} } =  P\cdot (N+Z^{2})
\end{align}
}

左辺の左の2つの項を右辺へ移動して、二乗根を取り除ける体制に、

{ \displaystyle
\begin{align}
Z\cdot\sqrt{\frac{X\cdot(N-X)}{N}+\frac{Z^{2}}{4} }  =  P\cdot(N+Z^{2}) - (X + \frac{Z^{2}}{2})
\end{align}
}

慎重にやったので問題ないだろう。 左辺の二乗

Z^2*(X*(N-X)/N+Z*Z/4)

から、右辺

P*(N+Z*Z)-(X+Z*Z/2)

を引いたら0という形で、maximasolve() に再びかける。

(%i13) display2d:false;
(%o13) false
(%i14) solve(%o6-%o7^2,[N]);
(%o14) [N = -Z^2,N = -(Z*sqrt((P^2-2*P+1)*Z^2+(4-4*P)*X)+(P-1)*Z^2-2*X)/(2*P),
        N = (Z*sqrt((P^2-2*P+1)*Z^2+(4-4*P)*X)+(1-P)*Z^2+2*X)/(2*P)]

2つの解が出たが、サンプルサイズは正の数であり、後者だろう。よって以下が求める式である。

{ \displaystyle N = \frac{Z\sqrt{(P^{2}-2P+1)Z^{2}+(4-4P)X}+ (1-P)Z^{2}+2X}{ 2P } }

ところで、X はどうやって決めればよいのだろうか。普通の解釈をすれば、X は動的な値である。すなわち、サンプルの取得中に、X は増えていくので、それに応じて N を増やすような形になる。 ところが、多くの場合は動的に取得できることは少ないので、サンプル数は増えるかもしれないが前回の値を直接使えばいいだろう。これは、信頼上限をこえた場合は一部のサンプルではなく、全取得するという前提である。

python で簡単なものを書いてみた。CTRは 0.5% 未満には興味がないという例である。

import math
def sample_size(p, x, z=1.96):
    return (z*math.sqrt((p**2-2*p+1)*z**2+(4-4*p)*x)+(1-p)*z**2+2*x)/(2*p)
def wilson_ucb(n, x, z=1.96):
    return (x+z*z/2+z*math.sqrt(x*(n-x)/n+z*z/4))/(n+z*z)

0クリックの枠

# 1 クリックもなかった広告枠に必要なサンプルサイズは
sample_size(0.005, 0)
764.4784
# そのサンプルサイズで0クリックだった場合のCTRの信頼上限は
wilson_ucb(765, 0)
0.004996607883860603

信頼上限が 0.5% を超えることはなかった。再び配信数(インプレッション)は制限される。

CTR が0.5% よりずっと低い枠

# 1万インプレッションで、5 クリックあった場合、CTR=0.05%だが、この枠のサンプルサイズは?
sample_size(0.005, 5)
2336.485083416135
# 2337 インプレッション取得したところ、1クリックあった。CTRの信頼上限は、
wilson_ucb(2337, 1)
0.002068315942437113

信頼上限が、0.5% を超えることはなかった。再び配信数(インプレッション)は制限される。

# 1クリックだったのでサンプルサイズは?
sample_size(0.005, 1)
1129.0503977791236
# 0 クリックだったのでサンプルサイズは?
sample_size(0.005, 0)
764.4784
# 以上、クリックがなければ上記か、クリックがあれば、2つ上を繰り返す。

CTR が 0.5% 付近の場合

# 1000インプレッションで、5 クリックあった場合。
sample_size(0.005, 5)
2336.485083416135
wilson_ucb(1000, 5)
0.011242912497455987
# 信頼上限は1.1%。この時点で、全取得と思われるが、あえてそれをせずシミュレーション
# 2337 インプレッションで、11 クリックあった場合。
sample_size(0.005, 11)
3934.259645857473

徐々に増えていくことがわかる。いずれ、その広告枠のインプレッション数の最大値に到達することだろう。

結論

もはや適切なサンプルサイズを調べるために、当選確実の計算方法や、統計のページを検索する必要はない。

N = (Z*sqrt((P^2-2*P+1)*Z^2+(4-4*P)*X)+(1-P)*Z^2+2*X)/(2*P)

が、必要なサンプルサイズだ。ここで、Nを全体の試行回数(求めたいサンプルサイズ)、 X を成功回数、Z正規分布の分位数、P をこれ未満にしたい信頼上限とする。