Yao, W., Li, Z., Graubard, B.I. (2014) Estimation of ROC curve with complex survey data. Statistics in Medicine, 34(8), 1293-1303.

　仕事の都合で読んだやつ。複雑な標本抽出デザインのデータからROC下面積(AUC)とその分散を推定するにはどうしたらいいか、という論文である。
　うーん、マニアックな話だなあ。夢もへったくれもないなあ。俺もChatGPTのマーケティング活用だとかそういうの語りたいのに。なんかこう、カメラの前でろくろとか回したいのに…

いや、本物のろくろを回したいわけじゃなくて… なんだかビジネスっぽいことを云って虚栄心を満たしたいなあ、という話です。

1. イントロダクション
　ご存じのように、縦軸にtrue positive rate, 横軸にfalse positve rateをとり、あらゆる閾値について点を打ってできる曲線をROC曲線という。その要約のために使われているのが下面積AUCである。
　AUCを推定したり比較したりする方法はいろいろ提案されている。AUCはランダムに選んだ患者の値がランダムに選んだ健常者の値以上である確率だと解釈できるので、ノンパラメトリックな方法についていえば、経験的に推定されたROCのAUCはMann-WhitneyのU統計量と同じである。分散推定についてもいくつか提案がある。
　さて、complex標本デザインにおいては[なんと訳せばいいのかいつも困るけど、層別抽出とかクラスタ抽出とかね]、デザインをちゃんと考慮しないと分散はふつう過小評価される。ウェイティングがまずくてもバイアスが起きる。単純無作為抽出のときよりずっとややこしい。
　本論文では、層別単純無作為標本から層別多段クラスタ標本までの範囲について、AUCのノンパラ推定とその分散推定の方法を提案します。

2. complex samplingデザインにおけるAUC推定
2.1 標準的ノンパラメトリック法
　罹患者にテストして値\(X_i (i=1, \ldots, M)\)を、非罹患者にテストして値\(Y_j (j=1,\ldots, N)\)を得た。経験的ROCはこう書ける。\(\psi(X_i,Y_j)\)を[原文にはカンマがないがいれた] \(Y_j \lt X_i\)のとき1, 等しいとき1/2, 逆のとき0となる関数として、$$ \hat\theta = \frac{1}{MN} \sum_j \sum_i \psi(X_i, Y_j)$$ 値が連続的ならタイは無視できるので期待値は母AUC \(\theta\)となる。つまり不偏推定量である。

2.2 層別単純無作為標本のAUC推定
　有限母集団のサイズを\(T\)、層\(y\)のサイズを\(T_h = N_h + M_h\)とする。標本サイズ\(t_h = n_h + m_h\)は定数だが内訳は確率変数である。
　標本への同時包含確率を\(\pi_{(hi, h’j)}\)と書く[層\(h\)の\(i\)さんと層\(h’\)の\(j\)さんがともに標本に包含される確率、という意味であろう]。その逆数\(w_{(hi, h’j)}\)を同時標本ウェイトとする。
　層\(h\)の罹患者と層\(h’\)の非罹患者に関する無限母集団AUCを\(\theta_{(h, h’)}\)と書く。[んんん？ 基本的に有限母集団の話をしているはずなのに、なぜここでは無限母集団だと書いてあるのかしらん]

　有限母集団に関して、層をプールしたAUCの推定量はこうなる。同時包含のインジケータを\(I_{(hi,h’j)}\)、罹患者のインジケータを\(\delta_{hi}\)として$$ \hat{\theta}_U = \frac{1}{\hat{N}\hat{M}} \sum_{h=1}^H \sum_{h’=1}^H \sum_{j=1}^{T_{h’}} \sum_{i=1}^{T_h} w_{(hi, h’j)} I_{(hi,h’j)} \delta_{hi} (1-\delta_{h’j})\psi(X_{hi},Y_{h’j}) $$ 標本ウェイト\(w_{(hi, h’j)}\)は、\(h=h’\)なら\(\frac{T_h(T_h-1)}{t_h(t_h-1)}\)だし、そうでなければ\(\frac{T_h}{t_h} \frac{T_{h’}}{t_{h’}}\)である。\(M_h, N_h\)はわかんないけど、ふつう\(T_h\)は既知だから、\(\hat{M}_h = \frac{m_h}{t_h} T_h\)という風に推定して、\(\hat{p}_h = \hat{M}_h / \hat{M}, \hat{q}_h = \hat{N}_h/ \hat{N}\)とすると、$$ = \sum_{h=1}^{H} \sum_{h’=1}^{H} \hat{p}_h \hat{q}_{h’} \hat{\theta}_{(h, h’)} $$ と書ける。
　[ここでの\(i,j\)は標本の罹患者・非罹患者の番号ではなく、単に母集団成員の番号であることに注意。標本のなかの罹患者\(i\)さんと非罹患者\(j\)さんのペアを取り出すために\(I_{(hi, h’j)} \delta_{hi} (1-\delta_{h’j})\)を掛けているわけだ。つまり、$$ \hat{\theta}_{(h, h’)} = \frac{1}{m_h n_{h’}} \sum_j \sum_i I_{(hi, h’j)}\delta_{hi} (1-\delta_{h’j}) \psi_{(X_{hi}, Y_{h’j})} $$ だということであろう。これを代入すると、$$ \hat{\theta}_U = \frac{1}{\hat{N} \hat{M}} \sum_h \sum_{h’} w_{(hi, h’j)} m_h n_{h’} \hat{\theta}_{(h, h’)} $$ めんどくさいからここから先は追いかけないけど、\(\frac{1}{\hat{N} \hat{M}} w_{(hi, h’j)} m_h n_{h’} \)を整理すると \(\hat{p}_h \hat{q}_{h’}\)になるんじゃないかと思う]

　\(\hat{p}_h\)からなる長さ\(H\)のベクトルを\(\hat{p}\), \(\hat{q}_h\)の同様のベクトルを\(\hat{q}\)と書き、\(\hat{\theta}_{(h, h’)}\)からなる\(H \times H\)行列を\( \hat{\theta}_{(H \times H)} \)と書けば、$$ =\hat{p}^\top \hat{\theta}_{(H \times H)} \hat{q} $$と書ける。

　これ、比推定量なので、不偏推定量というわけではないことに注意。分子と分母は元の\(\hat{\theta}\)の分子と分母の不偏推定量なので、デザイン一致推定量であるとはいえる。

2.3 層別多段クラスタ抽出のAUC推定
　層別2段クラスタ抽出について考えます。最初の層別クラスタ抽出は各層で独立にやり、抽出率が層によって違うとする。
　たとえばですね、PSU(一次抽出単位)の母集団サイズも標本サイズも同じで二段目が無作為なら、ウェイトはいらんです。PSUの母集団サイズが等しくなくて、一段目の抽出確率も等しくないとしても、たとえばそれがPPS抽出で、かつ標本PSUからの2段目抽出がほぼ同じ標本サイズの単純無作為抽出なら、やっぱりウェイトはいらない。[ぶひー、ややこしい。層が市,区,町,村の4つでPSUが市区町村だとしますね。各層において市区町村が抽出される確率がその人口に比例しており、かつ抽出した各市区町村から10人ずつ単純無作為抽出するんなら、ウェイトはいらないでしょ、という話だと思う]

　有限母集団が\(H\)枚の層を持ち、層\(h\)が\(K_h\)のPSUを持ち、PSU\(hg\)の母集団サイズが\(T_{hg} = N_{hg} + M_{hg}\)だとする。抽出するPSUの数を\(k_h\)、PSU\(hg\)から抽出する人数を\(t_{hg} = n_{hg} + m_{hg}\)とする。
　AUCの推定量は、例によって$$ \hat{\theta}_U = \frac{1}{\hat{N}\hat{M}} \sum_{h=1}^H \sum_{h’=1}^H \sum_{g=1}^{K_h} \sum_{g’=1}^{K_{h’}} \sum_{j=1}^{T_{h’g’}} \sum_{i=1}^{T_{hg}} w_{(hgi, h’g’j)} I_{(hgi, h’g’j)} \delta_{hgi} (1-\delta_{h’g’j}) \psi(X_{hgi}, Y_{h’g’j}) $$ となる。
　同時標本ウェイト\(w_{(hgi,h’g’j)}\)は、層もPSUも同じ二人なら\( \frac{K_h}{k_h} \frac{T_{hg}(T_{hg}-1)}{t_{hg}(t_{hg}-1)}\)、層は同じだけどPSUが違う二人なら \(\frac{K_h(K_h -1)}{k_h(k_h -1)} \frac{T_{hg} T_{hg’}}{t_{hg} t_{hg’}}\)、層の違う二人なら\( \frac{K_h K_{h’} T_{hg} T_{h’g’}}{k_h k_{h’} t_{hg} t_{h’g’}} \)となる。さきほどのようにベクトルと行列で略記することもできよう。

2.4 分散推定
　2つの方法がある。
　その1、ジャックナイフLOO。\(hg\)さんを抜いたときの推定量を\(\hat{\theta}_{(hg)}\)として、$$ \hat{V}_{jk} (\hat{\theta}_U) = \sum_h^H \left[ \frac{k_h -1}{k_h} \sum_g^{k_h} (\hat{\theta}_{(hg)} – \hat{\theta}_U)^2 \right] $$ である。[うわあ、計算が大変そうだ]
　その2, BRR推定量。各層で抽出するPSUが2個のときに使える。[…略…]

2.5 ドメインのAUC推定
　たとえば男だけのAUCを求めたい、という場合も話は簡単で、男だけに絞って求めればよい。ウェイトは変わらない。分散推定も同様。

3. シミュレーション
[モンテカルロシミュレーションで推定性能を調べている模様。パス]

4. 適用例
[Hispanic Health and Nutrition Examination Surveyという調査で、測定したBMIで定義した肥満を自己申告のBMIとかtricep skinfold測定値とかsubscapular skinfold測定値とかで予測する、というのをやっている。えーと、上腕三頭筋皮下脂肪厚、肩甲骨下部皮下脂肪厚ってことかな？ なんか知らんけどパス]

5. 考察
　[…]
　ふつう同時抽出確率というのはわからないのでなんらかの仮定が必要になる。本論文の適用例では、抽出した家族のなかで対象者はポアソン抽出だと仮定した。[4章を読んでないからよくわからんが、たぶん二段抽出の二段目が、それぞれの人についての独立ベルヌーイ試行になっていたと仮定したという話であろう]
　層別抽出の場合の先行研究では、層別にAUCを求めその精度の逆数で重みづけして結合したりしていたが[そうそう、私も最初そうだと思った]、テスト値の分布が層によって違う場合にはバイアスが生じる。いっぽう提案した推定量はほぼ不偏である。層内の罹患/非罹患者ペアだけでなく層間のペアもプールしているからである。[なるほどね]
　今後の課題: もっと効率的な分散推定(テイラー近似を使うとか)。パラメトリックないしセミ・パラメトリックな疾患予測モデルへの拡張。
————
　このたび仕事の都合で、調査ウェイト付きデータからAUCを求めたいという場面があったんだけど、層別非比例抽出だったらMantel-Haenszelオッズ比みたいに層別にAUCを推定して重みづけ平均するのかな… なんだか妙な話だな… よくわからんからウェイトを抽出確率にしてブートストラップ推定するしかないか… やだなあ、めんどくさいなあ… などと思っていたのである。まさか、層の総当たりペアごとにAUCを推定して重みづけ平均するとはおもわなんだ。勉強になりましたですー。
　サーヴェイ・データ分析のあるあるだけど、ローデータに個体ウェイトがついている、抽出デザインまでは遡れない、という場合はどうするんだろう？ 抽出が個体間で独立だとみて(そうみるしかないよね)、二人の同時ウェイトを単に個体ウェイトの逆数の積の逆数とし、患者と健常者の総当たりペアについて\(\psi\)の重みづけ和をとってペア数で割る、でいいのかな？ うーん、総当たりペア数が大きくなるとそれなりに大変そうだな。ブートストラップよりはましだけどさ。