ここまでの三つの記事では、何を選び直すかを決めれば標準誤差と信頼区間を計算できる、という流れで進めてきました。ただし、抽出単位を正しく指定しても、ブートストラップが元の推定誤差を再現できない推定量があります。
分布の上限を標本最大値で推定する場合が、その代表です。観測値を選び直す限り、元標本の最大値を超える値は現れないため、真の上限は区間の外に置き去りになります。イベントが一度も観測されなかったデータや、再抽出標本で回帰が推定できなくなる場面も、回数を増やすだけでは扱えません。
ここでは、一様分布からの最大値で近似が崩れる様子を実測し、小標本とまれなイベント、推定に失敗した反復の扱い、そしてブートストラップと置換検定の違いまでを順に扱います。抽出単位の決め方は、回帰・多施設・時系列での抽出単位で扱いました。
最大値では近似が崩れる例
観測値を選び直す方法では、元標本の最大値を超える新しい観測値は現れません。これは平均の標準誤差を求める時には直ちに問題になりませんが、分布の上限を標本最大値で推定する場合には、近似の成立に直接関わります。
ここでは、0から1の一様分布から30個の値をランダムに生成したデータセットを使います。真の分布上限は1ですが、30個の最大値は1より小さくなります。そして、30個をそのまま復元抽出して最大値を計算しても、その最大値を超えられません。
set.seed(42)
measurements <- runif(30)
maxima <- boot(measurements, statistic = function(data, indices) max(data[indices]), R = 10000)
quantile(maxima$t, probs = c(0.025, 0.975))
png("bootstrap-maximum-failure.png", width = 1100, height = 688, res = 200, type = "cairo")
plot(ecdf(maxima$t[, 1]), xlim = c(0.88, 1.01), main = "", pch = 19,
xlab = "Bootstrap sample maximum", ylab = "Cumulative probability")
abline(v = c(max(measurements), 1), lty = c(2, 3), col = c("blue3", "red3"))
legend("topleft", c("Observed maximum", "True endpoint"), col = c("blue3", "red3"), lty = c(2, 3), bty = "n")
dev.off()#> 真の分布上限 1.0000000
#> 元標本の最大値 0.9888917
#> 再抽出最大値の2.5%点 0.9400145
#> 再抽出最大値の97.5%点 0.9888917パーセンタイル区間の上限が0.9888917なので、真の上限1を含みません。しかも、この問題は今回の乱数だけでたまたま起きた結果ではありません。連続一様分布からの有限標本では最大値が真の上限未満になるため、通常の再抽出最大値によるパーセンタイル区間は、この設定で真の上限を覆えません。
最大値に残る確率の集中
元標本の最大値が一意なら、その対象が少なくとも一度抽出された時に、再抽出最大値$M_n^*$は元の最大値$M_n$と一致します。
$$\Pr_*(M_n^*=M_n) = 1-\left(1-\frac1n\right)^n \longrightarrow 1-e^{-1}\approx0.632 \tag{1}$$
式(1)は、最大値を持つ観測が一度以上選ばれる確率を求めています。この確率は$n$を増やしても0へ消えず、再抽出最大値が元の最大値に一致する集中が残ります。
実際、$n=30$ではこの確率は約0.6383で、今回の10,000回での割合は0.6350でした。母集団から新しく30個を取った最大値の分布には、未知の真の上限に同じ種類の点の集中はありません。この違いが、通常のブートストラップが最大値の誤差分布を再現できないことにつながります。
もっとも、基本区間などでは、計算上は観測最大値より上へ区間を伸ばせる場合があります。しかし、区間が観測範囲を超えられることと、ブートストラップ分布が元の誤差分布と合っていることは、別の話です。
そこで、分布の端点や極値が研究の中心なら、その対象に合うパラメトリックなモデル、極値に関する推論、条件に応じた部分標本法などを考えます。図1の分布を滑らかに描き直したり、反復回数を100,000回へ増やしたりするだけでは、元の観測情報の制約は変わりません。
小標本とまれなイベント
「20人未満は使えない」「30人以上なら使える」という共通の人数基準はありません。平均、中央値、99%点、イベント割合では、必要な情報と推定の不安定さが異なるためです。
そして、今回の20人の平均変化は計算手順の教材で、20人ならBCaの95%区間が十分に正確という証明ではありません。また、十分な標本数かどうかは、独立な対象者の数だけでなく、推定量の分布、外れ値、カテゴリごとの人数、再抽出でのモデルの成立にも関わります。人数ではなくイベント数が推論を支えるという関係は、必要イベント数でも扱っています。
イベントが0件の場合
100人でイベントが一度も観測されなかった場合、通常の行の復元抽出では全ての標本がイベント0件になります。得られる標本割合は全て0で、その標準偏差も0になりますが、この0は観測された標本の性質であって、母集団のイベント確率の推定結果とは違います。
ここで、独立な二項モデルを仮定すると、イベント確率$p$の下で100人全員が非イベントになる確率は$(1-p)^{100}$です。これを0.05とする値を求めると、片側95%信頼上限は約2.95%になります。これは通常のノンパラメトリック・ブートストラップとは別の、二項モデルに基づく計算です。
また、イベントが100人中1件の場合も、100回の復元抽出でその1件が一度も選ばれない確率は$(99/100)^{100}$で約36.6%です。そのような標本では、単純な割合は計算できても、イベントに関する回帰係数が推定できない場合があります。
したがって、イベントの少ない解析で、推定できた標本だけを集めれば十分と考えるのは注意が必要です。推定できない標本を大量に除くと、安定して推定できる場合に条件を付けた分布へ変わってしまいます。
中央値と極端な分位点
連続分布の中央値は、中央値付近の密度が正となることなどの条件の下で、通常のブートストラップが成り立つ代表的な推定量の一つです。元のデータを少し動かした時の標本中央値の見た目が滑らかでないことだけで、ブートストラップ不可と判断しません。
ただし、人数が少ない場合、0が多い場合、丸められた同値が多数ある場合には、再抽出中央値が少数の値へ集中することがあります。また、99%点を推定するなら、100人ではその上側に期待される観測は約1人で、1,000人でも約10人です。
そして、この情報量でどのような区間を作れるかは、分布の裾に関する仮定と切り離せません。中央値で使えた手順を99%点へ移す時には、同じ「分位点の推定」として処理をコピーする前に、端点付近の観測数を確認します。
外れ値と裾の重い分布
元標本に一つだけ大きな測定値があると、その値が0回、1回、複数回現れることで、再抽出推定値が大きく変わる場合があります。ブートストラップは外れ値を自動的に無害化する方法ではなく、平均のように外れ値へ敏感な推定量は、その感度を持ったまま再計算されます。
しかも、その測定値が記録誤りなのか、実際に起き得る値なのかは、再抽出分布だけでは判断できません。また、都合のよい区間になるまで観測値を削除すると、解析対象が結果に応じて変わってしまいます。
さらに、分布の裾が極端に重く、平均や分散が存在しない設定では、通常の平均の標準誤差や平方根$n$による近似自体が成り立たない場合があります。有限標本から小さな数値の標準誤差が出ても、母集団で同じ量が定義できることの確認にはなりません。
推定できない標本の扱い
ロジスティック回帰でイベントが消える、カテゴリの一つが欠ける、共線性で係数が識別できなくなるなど、元標本では動いた解析が再抽出標本で止まることがあります。そして、その失敗はプログラム上の問題だけでなく、元の情報の少なさを表している場合があります。
ここで、10,000回のうち2,000回が失敗したなら、8,000回も成功したという見方だけでは不十分です。どのような標本で失敗したかを調べ、その2,000回が推定値の分布の特定の領域に対応していないかを確認します。ここでの反復回数と失敗件数は説明用の設定で、今回のcd4の実行で起きた件数ではありません。
そのため、非有限の推定値を黙って取り除いたり、失敗がなくなるまで条件を変えながら追加抽出したりすることは推奨しません。成功した解析に条件を付けた分布を、予定していたブートストラップ分布として報告することになるためです。
そこで、解析方法そのものに正則化や小標本向けの推定が必要なら、元標本と再抽出標本で同じ方法を使う方針へ組み直します。また、カテゴリを統合する場合には研究上の意味を先に定め、特定の再抽出標本だけで都合よく分類を変えないようにします。
検定と置換検定の違い
通常のブートストラップは、元の推定値の周辺で標本変動を近似します。そのため、推定した治療効果を残したまま生成した値を並べると、その効果を中心とする分布が得られるだけで、治療効果0という帰無仮説の下での検定分布とは別のものになります。
この点について、HallとWilsonは、ブートストラップ検定で、中心化や標準化を含めて何の分布を近似するかを論じています。ブートストラップ分布の中で元の推定値より大きい値の割合を、そのままp値とする方法では、検定したい帰無仮説との対応がありません。
たとえば、平均変化0を検定したい場合、元標本の平均を差し引いたブートストラップ誤差を使う方法や、帰無仮説を満たすデータ生成モデルを使う方法が候補になります。複合仮説や回帰モデルでは、残るパラメータや分散をどう推定するかも関わります。
また、ブートストラップ推定値が0を超える割合だけから両側p値を作る近道も、特定の区間の反転と近似的に対応する場合はありますが、任意の推定量で有効な検定になるわけではありません。使用した信頼区間と、p値の計算原理を別々に明らかにします。
ラベルを入れ替える置換検定
これに対して置換検定では、観測値を複製して標本を作る代わりに、群ラベルなどを入れ替えて統計量を計算します。無作為割付の手順に基づくランダム化検定や、帰無仮説の下での交換可能性に基づく検定と結び付いています。
ただし、2群の平均が等しいという帰無仮説だけでは、群間で分散や分布形が異なってもラベルを自由に入れ替えられる、とは限りません。対応のある研究なら対応を保つ入れ替え方、層別割付ならその割付構造に合う方法が必要です。群間比較で重み付けをどう決めるかという問題は、ログランク検定とも共通します。
したがって、ブートストラップと置換検定は、いずれも再標本化の計算を使いますが、近似している対象が違います。推定誤差を近似するのか、帰無仮説や割付の下での比較を行うのかを区別します。
複雑な推定手順と非正則性
傾向スコアを推定して重みを作り、アウトカムモデルを当てはめるなど、複数段階の推定では、元データから学習した各段階が標本変動に関わります。重みや変数変換を元標本の値に固定したまま最後の回帰だけをやり直すと、手順全体の揺れとは異なる量を評価する場合があります。
そのため、再抽出内で推定手順をやり直す範囲を明確にします。ただし、全ての段階を繰り返したことだけで、ブートストラップが常に妥当になるわけではありません。
その一例として、AbadieとImbensは、固定数の最近傍を使うマッチング推定量について、通常のブートストラップが一般には成立しないことを示しています。推定量が漸近的に正規となる設定でも問題が生じるため、「最終的には正規に近づく」だけでは十分な根拠になりません。
ただし、この論文が扱うのは、固定数の最近傍を使う特定の推定量です。別のマッチング法や因果効果推定へ結果を広げる時には、傾向スコアマッチングで扱った推定対象の違いも含めて、使っている方式に対応した理論や解析的分散、別の再標本化法の成立条件を確認します。
そして、変数選択後の係数、0の境界にある分散成分、鋭い閾値で選別する推定量などにも、通常の滑らかな推定量とは異なる問題が生じ得ます。推定手順の複雑さだけを理由にブートストラップへ任せるのではなく、その手順がどの変化に不連続に反応するかを調べます。
近似が崩れる条件を先に確かめる
この記事で並べた場面には、共通するものがあります。推定量が観測範囲の端や、標本ごとの成否に強く依存していることです。
実際、最大値では元標本の最大値に確率が集中し、イベント0件では全ての標本が0になり、推定できない標本を除くと成功した標本に条件を付けた分布へ変わります。どれも、10,000回を100,000回にしても同じことが起こります。
私なら、再抽出を回す前に、推定したい量が観測値の並べ替えでどれくらい動くかを手元で確かめます。1回除いただけで大きく動く量や、値の集中が見える量なら、その時点で別の推論方法を探した方が早いためです。
次の記事では、推定誤差ではなく開発手順の楽観性を対象にする使い方、多重代入との組み合わせ、生存時間解析での再抽出、そして論文へ何を書くかを扱います。
参考文献
Bickel PJ, Freedman DA. Some asymptotic theory for the bootstrap. The Annals of Statistics. 1981;9(6):1196-1217.
Hall P, Wilson SR. Two guidelines for bootstrap hypothesis testing. Biometrics. 1991;47(2):757-762.
Efron B. Bootstrap methods: Another look at the jackknife. The Annals of Statistics. 1979;7(1):1-26.
Davison AC, Hinkley DV. Bootstrap Methods and Their Application. Cambridge University Press; 1997.
