'Ntile()' の問題点
The Trouble with 'Ntile()' (blog.djnavarro.net)
要約
R言語のdplyrパッケージに含まれるntile()関数は、データを指定された数のグループに分割する際に、統計的な精度が重要な場合には使用すべきではないと警告しています。この関数は、SQLのNTILE()関数とは異なり、データの正確な分位点に基づいたビン分割を行わず、予期しない動作をすることがあります。記事では、薬剤開発の曝露-反応分析を模倣したデータセットを用いて、ntile()がどのように誤解を招く結果を生むかを具体的に解説しています。
全文翻訳
Ntile() の問題点
You keep using that function. I do not think it does what you think it does
このブログの長年の読者であれば、私がtidyverseのファンであることはご存知でしょう。Rプログラミング言語を深く愛していますが、ベースRの挙動にはいくつかの奇妙な点があり、tidyverseが非常にうまくこなしていることの多くは、それらの粗さを滑らかにすることです。ggplot2でのデータ可視化、dplyrでのデータ操作、あるいは一貫性のないベースRの正規表現ツールと格闘する代わりにstringrを使った文字列操作の地獄をナビゲートすることであれ、それは恵みでした。しかし、当たり前のことをわざわざ記事にするようなことはしないので、tidyverseの基盤となっているfsやcliのような、あまり注目されないパッケージについて書くことが多いです。これらは皆の生活を少し楽にしてくれます。
要するに、私はファンなのです。
そして、tidyverseパッケージがうまく機能しないことについて記事を書くことに、ある種の不安と後悔を感じています。結局のところ、私が批判できる立場にあるでしょうか？私は長年、ひどいコードを書いてきましたし、人々が実際に使うパッケージに、非常に賢明でない選択肢が含まれてしまったこともあります。しかし、それがポイントでもあります…どんなに優れたツールにも問題はありますし、どんなに優れたプログラマーでも間違いを犯します。それを否定して平然を装うのは愚かです。
それを念頭に置いて、この記事はdplyr::ntile()についてのものです。統計的な正確さが重要な目的には、これを使用すべきではないという警告です。もし、関数の名前に基づいてそれが何をするのかを推測した場合（そしてSQLを使った経験がない幸運な人の一人であれば）、それが実際に何をするのかを発見したときに痛い目を見ることになります。これは、データを分位数ベースのグループにビン分割するためのツールでは断じてなく、データ分析に使用すると間違いなく誤動作します。どうか注意してください。
library(dplyr)
library(ggplot2)
［画像: Rのライブラリ読み込みを示すコードブロックの後に表示される画像］
便利なデータセット
まず、ntile()を使用する際に発生する問題をわずかに誇張するように意図的に設計された架空のデータセットer_dataを紹介します。このデータセットは、薬物動態学の分野でよく行われる曝露-反応（ER）分析のような状況を模倣しています。設定は以下の通りです。ブログの読者のほとんどが薬物開発に携わっていないため、詳しく説明します。3つの研究からのデータがあります。
- 研究
S01は、「単回漸増投与」デザインの第1相用量漸増試験です。この詳細はこのおもちゃの例ではそれほど重要ではありませんが、重要なのはこれが逐次設計であるということです。最初の被験者群は非常に低い用量を受け取ります。副作用がなければ、次の群はより高い用量を受け取ります。そして、その逆も同様です。このようなデザインの研究を含めることは、安全性の目的のために非常に重要です。したがって、アウトカム変数が安全性エンドポイントである場合、データセットにはこのようなものが含まれている可能性が高いです。通常、この研究のすべての被験者は男性です。1 この研究で妊娠している可能性のある人がいるリスクを冒したくないからです。この研究で妊娠している可能性のある人がいるリスクを冒したくないからです。 - 研究
S02は、第1相研究の一部である「薬物間相互作用」研究です。これは特に経口避妊薬に関連しています。新しい薬がピル剤の効果を変える（またはその逆）リスクを冒したくありません。そのため、このような研究を目にする可能性も十分にあります。驚くことではありませんが、この研究のすべての参加者は女性であり、用量レベルは固定されています。 - 研究
S03は、通常、薬効を特定するために行われる第2相「用量設定」研究です。ここでは、より幅広い参加者（男性と女性の両方が含まれる）がおり、用量レベルにもばらつきがあります。
この状況では、研究と性別ごとの被験者数の集計は、次のようになります。
er_data |> count(study_id, phase, description, sex)
# A tibble: 4 × 5
study_id phase description sex n
<chr> <dbl> <chr> <chr> <int>
1 S01 1 SAD dose-escalation M 31
2 S02 1 DDI oral contraceptive F 24
3 S03 2 Dose-finding F 60
4 S03 2 Dose-finding M 60
架空のデータセットが曝露-反応分析シナリオを大まかに模倣しているため、データセットには薬物曝露の典型的な測定値（例: cmaxはピーク薬物濃度、aucは一定期間の総薬物曝露の「面積」を測定）と、関心のあるresponse測定値（例: 安全性指標、有効性指標など）が含まれています。現在の投稿ではそれほど重要ではありませんが、感覚をつかむために、このデータにおける曝露-反応関係は次のようになります。
er_data |>
ggplot(aes(auc, response)) +
geom_point(aes(color = study_id)) +
geom_smooth(formula = y ~ x, method = "lm", color = "#222")
［画像: 曝露-反応関係を示す散布図と回帰直線］
薬物曝露が高いほど、反応は強くなります。これは通常、曝露-反応分析を行う上で関心のあることですが、現在の投稿の中心ではないため、先に進みます。
しかし、現在の投稿の中心となるのは、データセットが体系的かつ合理的な方法で整理されていることです。er_dataの各行は特定の被験者に対応し、各列は特定の測定値に対応します。このデータセットを準備したデータプログラマーは無意味に意地悪ではないため、行はランダムに並べられていません。代わりに、分析者が理解しやすいように並べられています。行はstudy_id、次にdose_mg、そしてsubject_idで並べられています。
er_data
# A tibble: 175 × 11
row_id study_id phase description subject_id sex dose_mg wt_kg auc cmax response
<dbl> <chr> <dbl> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
1 1 S01 1 SAD dose-escalation S01-001 M 10 88.5 1.87 0.265 16.8
2 2 S01 1 SAD dose-escalation S01-002 M 10 82.3 1.71 0.262 29.1
3 3 S01 1 SAD dose-escalation S01-003 M 10 117 1.21 0.129 16.9
4 4 S01 1 SAD dose-escalation S01-004 M 10 85.7 0.871 0.0828 19.0
5 5 S01 1 SAD dose-escalation S01-005 M 10 83.9 1.54 0.147 16.6
6 6 S01 1 SAD dose-escalation S01-006 M 10 77.2 2.21 0.136 29.3
7 7 S01 1 SAD dose-escalation S01-007 M 10 77.2 1.31 0.106 19.8
8 8 S01 1 SAD dose-escalation S01-008 M 30 87.9 4.93 0.426 26.1
9 9 S01 1 SAD dose-escalation S01-009 M 30 84.1 4.46 0.494 20.3
10 10 S01 1 SAD dose-escalation S01-010 M 30 77.2 4.93 0.411 19.7
# ℹ 165 more rows
列名は、各変数が何を表すかを正確に示しています。
row_idはブックキーピング目的で存在し、元の行番号を含みます。study_idはデータがどの研究から来たかを示します。phaseは第1相研究か第2相研究かを示します。descriptionは研究の簡単な説明です。subject_idは各被験者のユニークな識別子を提供します。sexは被験者が男性か女性かを示します。dose_mgは投与された用量（ミリグラム単位）を指定します。wt_kgは体重（キログラム単位）を指定します。aucとcmaxは2つの曝露メトリックです。responseは、最も想像力のない方法で名付けられた応答変数です。
ここでも、ほとんどのことは本題とは関係ありません。最も重要なのはwt_kg変数です。これは、人の体重を小数点以下1桁までしか記録していません。実際の体重計は通常、そのレベルの精度で体重を表示するため、小数点以下1桁に丸められています。2 その結果、体重は「理論的には」連続的に変化する量ですが、実際のデータセットではそうではありません。すべての実生活のデータ分析では、「同じ」体重を持つ被験者がかなり多く存在することになります。これは、分布の中央付近で最も頻繁に発生します。
er_data |> filter(wt_kg == median(wt_kg))
# A tibble: 14 × 11
row_id study_id phase description subject_id sex dose_mg wt_kg auc cmax response
<dbl> <chr> <dbl> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
1 6 S01 1 SAD dose-escalation S01-006 M 10 77.2 2.21 0.136 29.3
2 7 S01 1 SAD dose-escalation S01-007 M 10 77.2 1.31 0.106 19.8
3 10 S01 1 SAD dose-escalation S01-010 M 30 77.2 4.93 0.411 19.7
4 12 S01 1 SAD dose-escalation S01-012 M 30 77.2 4.14 0.252 11.6
5 23 S01 1 SAD dose-escalation S01-023 M 100 77.2 16.1 1.34 42.6
6 38 S02 1 DDI oral contraceptive S02-038 F 100 77.2 14.2 1.65 36.8
7 40 S02 1 DDI oral contraceptive S02-040 F 100 77.2 26.1 2.62 41.8
8 41 S02 1 DDI oral contraceptive S02-041 F 100 77.2 13.3 1.42 36.5
9 48 S02 1 DDI oral contraceptive S02-048 F 100 77.2 19.7 1.68 38.5
10 71 S03 2 Dose-finding S03-071 F 50 77.2 8.99 0.731 32.1
11 86 S03 2 Dose-finding S03-086 M 50 77.2 17.3 1.29 34.2
12 146 S03 2 Dose-finding S03-146 M 200 77.2 31.6 2.79 74.9
13 155 S03 2 Dose-finding S03-155 F 200 77.2 39.1 2.57 54.5
14 156 S03 2 Dose-finding S03-156 M 200 77.2 34.0 4.19 57.6
はい、タイがあります。
しかし待ってください…175行のデータセットで、中央値で14個の値がタイになっている？データセットの8%です。これは実生活では驚きですが、投稿の著者が意図的に問題を誇張するためにデータセットを準備した場合には、完全に予想されることです。これらのタイは、これから書く問題を引き起こすまさにそのものであり、私は問題をより簡単に見えるようにするためにデータセットを準備しました。とはいえ、私の例のデータは少し作り物ですが、ひどく非現実的というわけではありません。実生活でデータセットの8%が中央値でタイの値を持つということはあまり考えられませんが、私は実世界のデータでサンプルの中央値体重に2%のデータが集中している状況に遭遇したことが絶対にあります。この「中央値でのタイ値の塊」状況は、実際には非常によく見られます。
この非常に長い前置きを踏まえて、データ分析を行い、それを実行する際にあまりにも残念な間違いとしてntile()関数を適用してみましょう。
［画像: データセットの概要を示す表の後に表示される画像］
Ntile() の問題点
さて、まず非常に普通のことをしてみましょう。薬物動態学のデータ分析では、被験者を体重ベースのビンにグループ化し、体重ビンごとに曝露量を別々に調べることはよくあります。例えば、このようなものでしょうか？
er_data |>
mutate(wt_bin = factor(ntile(wt_kg, n = 4))) |>
ggplot(aes(wt_bin, auc)) +
geom_boxplot() +
facet_wrap(~study_id)
［画像: 体重ビンごとのAUCの箱ひげ図］
ここでも、詳細はあまり重要ではなく、これを非常に現実的にしようともしていません。「連続的」な変数をビンにグループ化し、次にビンごとにデータを要約するという、これに似たことを行うのは非常に一般的です。
さて、不快な部分に移りましょう。「連続的」な変数wt_kgは実際には連続的ではなく、精度が低いレベルでしか測定されていないため、データセットにはそれらのタイ値が含まれています。先ほど見たように、中央値体重と「完全に」同じ体重を持つ被験者が14人います。それらの人々はどうなるのでしょうか？彼らの体重は正確に中央値であるため、ビン2に割り当てることも、ビン3に割り当てることも同様に合理的です。「連続的」な変数を離散的なカテゴリにビン分割する場合、ブレークポイントに偶然存在するタイ値について、どうするかを決定する必要があります。
さて…この決定に対する「正しい」答えは何でしょうか？まあ…誰に聞くかによります。例えば、統計学者の誰かに尋ねると、順序統計量から分布の分位数を推定することについての非常に長く痛みを伴う議論に巻き込まれるでしょう。そして、泣き出す可能性が高いです。そこで、統計学者に尋ねるのはやめようと決め、代わりにntile()がそのトピックについて何を言っているのか尋ねることにします。
er_data |>
mutate(wt_bin = ntile(wt_kg, n = 4)) |>
filter(wt_kg == median(wt_kg)) |>
split(~wt_bin)
$`2`
# A tibble: 4 × 12
row_id study_id phase description subject_id sex dose_mg wt_kg auc cmax response wt_bin
<dbl> <chr> <dbl> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <int>
1 6 S01 1 SAD dose-escalation S01-006 M 10 77.2 2.21 0.136 29.3 2
2 7 S01 1 SAD dose-escalation S01-007 M 10 77.2 1.31 0.106 19.8 2
3 10 S01 1 SAD dose-escalation S01-010 M 30 77.2 4.93 0.411 19.7 2
4 12 S01 1 SAD dose-escalation S01-012 M 30 77.2 4.14 0.252 11.6 2
$`3`
# A tibble: 10 × 12
row_id study_id phase description subject_id sex dose_mg wt_kg auc cmax response wt_bin
<dbl> <chr> <dbl> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <int>
1 23 S01 1 SAD dose-escalation S01-023 M 100 77.2 16.1 1.34 42.6 3
2 38 S02 1 DDI oral contraceptive S02-038 F 100 77.2 14.2 1.65 36.8 3
3 40 S02 1 DDI oral contraceptive S02-040 F 100 77.2 26.1 2.62 41.8 3
4 41 S02 1 DDI oral contraceptive S02-041 F 100 77.2 13.3 1.42 36.5 3
5 48 S02 1 DDI oral contraceptive S02-048 F 100 77.2 19.7 1.68 38.5 3
6 71 S03 2 Dose-finding S03-071 F 50 77.2 8.99 0.731 32.1 3
7 86 S03 2 Dose-finding S03-086 M 50 77.2 17.3 1.29 34.2 3
8 146 S03 2 Dose-finding S03-146 M 200 77.2 31.6 2.79 74.9 3
9 155 S03 2 Dose-finding S03-155 F 200 77.2 39.1 2.57 54.5 3
10 156 S03 2 Dose-finding S03-156 M 200 77.2 34.0 4.19 57.6 3
ああ。
ああ、いやだ。
それは全く正しく見えませんね？
パッケージのウェブサイトにはあまり詳しく書かれていませんが、ntile()関数は行番号に基づいてタイをブレークします。計算上は便利なことかもしれませんが、残念ながら実世界では、あなたの「分位数グループ」は、認識されているどの分布分位数推定値にも基づいていません。現実的なデータセットで、観測値が構造的かつ意味のある方法で配置されている場合（つまり、適切にフォーマットされたすべてのデータセット）、ntile()は非常に非ランダムな方法でタイをブレークします。
このシナリオでは、中央値でタイになっている14個の値のうち4個がQ2に割り当てられ、その4個はすべて第1相用量漸増試験から来ています。すべて男性です。そして、用量漸増試験の中でも、Q2に割り当てられなかった観測値を見てみましょうか？それは100mg用量群からのものです。研究S01は用量漸増試験であり、逐次設計であることを思い出してください。データセットは、この試験内の被験者を合理的に用量レベルで並べているため、「中央値体重以上」とラベル付けされるS01被験者の決定は、薬の最高用量を与えられた被験者であることが保証されています。やばい。
言い換えれば、ntile()を使用したことでデータセットが分位数ベースの体重ビンに分割されたと信じたくなるかもしれませんが、そうではありません。私たちのwt_bin変数は、体重四分位グループのように見えるものですが、それは全く別のものです。データセットの中央値体重でタイ値があり、データが合理的に整理されているため、米国FDAに提出するレポートで、私たちが何をしたかを最も正確に説明する方法は次のようになります。
ベースライン体重に主に、研究IDにも少し、そして投与量にもわずかに基づいた「体重」四分位に被験者を割り当てました。
FDAまたは世界の他のどの規制当局も、これを分析の許容できる方法とは見なさないと確信しています。確かに、もしあなたが単に詳細をあまり気にしないCスイートエグゼクティブ向けの迅速で簡単なダッシュボードを準備しているだけなら、どうぞ…もし好きならntile()を使用してください。私はあなたを止めませんし、おそらくあなたのボスも止めないでしょう。しかし、新しい薬が人間に承認されるかどうかを決定するために使用される分析を行っている場合は？別のツールを選んでください。
［画像: FDAへの提出を想定した分析の重要性を示す画像］
代わりに何をすべきか？
さて、ダニエル、ntile()がこの目的に合わないことに同意しますが、代わりに何をすべきでしょうか？さて、他のダニエル、まず素晴らしい質問をありがとう、そして第二に、ntile()がビン分割を行う際に生じる2つの関連するが本質的に異なる懸念を認識することが役立つと思います。
-
タイブレーク規則は統計的に擁護できない。 上記の私の異議の核心は、
ntile()がブレークポイントでタイ値に遭遇したときに、実生活では統計的性質がひどいタイブレーク規則（行番号でタイをブレークする）を適用することです。3分類規則にタイブレーク手順を組み込む場合、それは（a）ほとんど文書化されておらず、（b）非常に非ランダムであるという両方の条件を満たすものであってはなりません。それは本当に、本当に悪いことです。 -
グループを分位数に基づかせたいのであれば、タイをブレークすべきではない。 2番目の点はもう少し微妙です。ほとんどの実世界の分析では、私たちの目標は実際に「等しいサイズのグループ」を作成することではありません。目標は、サンプルから計算された分位数に基づいてデータを分割することです。サンプルサイズが大きく、真に連続的なデータがある場合、これは違いのない区別です。しかし、このブログ記事で使用している
er_dataデータセットのような準連続データがある場合、これら2つは乖離します。分位数ベースの分割は、等サイズ分割と同じではありません。
これらはいずれもntile()を使用しない良い理由ですが、それらは引き起こす問題の深刻度において異なります。「非ランダムなタイブレーク」問題は非常に危険であり、統計分析のツールとしてntile()を自動的に失格させます。しかし、もしntile()が代わりにランダムにタイをブレークするとしたらどうでしょうか？タイブレークが適切に処理される限り、「等サイズ」分類規則を持つことは本当に悪いことでしょうか？
正直に言って、私は思います…おそらくそうではないでしょう？「分位数分割をFDAに送信し、等サイズ分割を送信しない」という統計学の根本的な法則はありません。連続変数を離散カテゴリに分割するたびに、必然的にある程度の忠実度を失うことになるため、分位数ベースの分割と（適切に構築された）等サイズ分割の違いについて、あまり心配するのは愚かに思えます。4等サイズ分割ではなく分位数分割を使用するのは、主に慣習によるものです。それにもかかわらず、慣習に従うつもりであれば、実際に慣習に従う必要があります。等サイズ分割を実行してから、分位数分割を実行したとFDAに報告することはできません。それは許容できません。
［画像: 分位数分割と等サイズ分割の概念的な違いを示す図］
分位数分割の実装
さて、分位数ベースのグループを実装することにしたと仮定しましょう。そのため、まずタイを全くブレークしないという原則から始めます。このため、エッジケースが「上に上がる」か「下に下がる」かについて、最初に決定する必要があります。5どちらの方向を選んでも構いませんが、選択する必要があります。2人の人が同じ体重であると記録されている場合、それらを異なる体重ビンに割り当てるべきではありません。はい、それは私たちの四分位分割に同じ数の被験者が含まれないことを意味し、場合によってはグループサイズが非常に不均一になる可能性がありますが、それは全く問題ありません。
実際には、次のようなヘルパー関数を作成するかもしれません。
cut_quantile <- function(x, n = 4, right_closed = TRUE) {
breaks <- quantile(x, probs = (0:n)/n, na.rm = TRUE)
bin_num <- as.numeric(cut(x, breaks, labels = 1:n, include.lowest = TRUE, right = right_closed))
bin_fct <- factor(bin_num, levels = 1:n, labels = paste0("Q", 1:n))
return(bin_fct)
}
上記のコードでは、責任の明確な分割があります。quantile()はデータセットの経験的分位数を計算する責任を負い、cut()は各観測値をその値に基づいてビンに割り当てる責任を負います。この役割分担は重要です。
quantile()は、分析者が信頼する推定値に基づいて、ビンを定義する明確なルールが常に存在することを保証します。cut()は、そのルールが、観測値のその値のみに基づいて、データセットのすべての観測値に全く同じ方法で適用されることを保証します。
175人の被験者がいるer_dataに適用すると、分位数ビン分割は、ほぼ均等なサイズの4つのグループを絶対に生成しません。タイ値の数が多いため、44/44/44/43の分割になる可能性はほとんどありません。分位数ビン分割の目的は、サンプルで等しいサイズのグループを作成することではなく、母集団で等しい分割がどのようになるかについての最良の推測に基づいてグループを定義することです。それは基本的に推論統計です。6その結果、cut_quantile()関数を使用して体重四分位数でデータをビン分割すると、次のようになります。
er_data |>
mutate(wt_bin = cut_quantile(wt_kg)) |>
summarise(.by = wt_bin, n_subj = n(), min_wt = min(wt_kg), max_wt = max(wt_kg))
# A tibble: 4 × 4
wt_bin n_subj min_wt max_wt
<fct> <int> <dbl> <dbl>
1 Q4 44 85.6 117
2 Q3 33 77.5 85.5
3 Q2 54 67.6 77.2
4 Q1 44 52.3 67.4
ここで重要なのは、各ビン内の最小値と最大値に注目することです。すべての体重ビンの最大値は、次の最も高い体重ビンの最小値よりも常に厳密に低いです。それはそうでなければなりません。データセットには同じ体重で記録されている被験者が多数いますが、タイをブレークしないため、ブレークポイントにあるタイ値は常に互いに同じ方向に移動します。この例では、cut()がデフォルトで「左開右閉」区間を使用するため、常に下に移動します。体重77.2kgの被験者14人全員が「Q2」に割り当てられ、そのうちの誰も「Q3」に割り当てられません。その結果、「Q3」グループの最小体重は77.5kgとわずかに高くなります。
分位数ベースのビンでは常にこのパターンを期待すべきです。観測値の値を割り当てられるビンにマッピングする決定論的なルールが存在することは、論理的な要件として従います。分位数ベースのビン分割を行っているため、タイをブレークすることは許可されておらず、そのテーブルに同じ数値が2回現れることはありません。
対照的に、ntile()が何をするかを見てみましょう。
er_data |>
mutate(wt_bin = factor(paste0("Q", ntile(wt_kg, 4)))) |>
summarise(.by = wt_bin, n_subj = n(), min_wt = min(wt_kg), max_wt = max(wt_kg))
# A tibble: 4 × 4
wt_bin n_subj min_wt max_wt
<fct> <int> <dbl> <dbl>
1 Q4 43 85.7 117
2 Q3 44 77.2 85.6
3 Q2 44 67.6 77.2
4 Q1 44 52.3 67.4
はい、そうです。今回は、私たちが単純に期待していた44/44/44/43の分割が得られ、私自身もこのきれいな集計表を見ることに満足するかもしれませんが、分位数ビン分割を実装するという目標にとっては、これは実際には間違ったことです。n_subjのカウントを見て喜ぶ代わりに、各「分位数ビン」内のmin_wtとmax_wtの値を見て、恐怖に身をすくめるべきです。なぜなら、これらの列を見ると、最初に気づくのは、「Q2」の最も重い被験者が…「Q3」の最も軽い被験者と全く同じ体重であるということです。これは丸め誤差ではなく、本物の衝突です。正確に同じ測定体重を持つ2人の被験者が異なる「体重ビン」に割り当てられました。それは無意味です。つまり…それは単に間違っています。ntile()を使用することで得られる、あの見栄えの良いビンカウントは、統計的なエラーを隠蔽しています。これは分位数分割ではありません。7
［画像: ntile()による分割結果の表］
慣習について
この記事のこの時点では、いくつかの小さな詳細を説明し、最後に非常に皮肉な結論を述べる以外に、あまり言うことはありません。まず小さな詳細から始めましょう。
この時点で、分位数ビンを定義するという意図を明確にし、それを行うためのツールも手に入れたので、残る実質的な決定は1つだけです。ブレークポイントでのタイ値は「上に上がる」（左閉右開区間）べきか、それとも「下に下がる」（左開右閉区間）べきか？想像できるかもしれませんが、どちらを選んでも本当に本当に構いません。しかし、科学的には全く違いがないにもかかわらず、慣習を意味のあるものと誤解すると、あなたを迷わせる可能性のある隠れた「慣習」を解きほぐすのに役立ちます。
私たちの徹底的に退屈な物語は、「分布」を定義するための慣習から始まります。慣習として、確率変数 (X) の累積分布関数（CDF）を、その変数 (X) の観測値が (x) を超えない確率 (\mbox{Pr}(X \leq x)) として定義します。er_dataデータセットの体重の経験的分布に適用すると、非常に便利なベースRの関数ファクトリecdf()の助けを借りて、経験的CDFを簡単にプロットできます。以下に示します。
# in-built ecdf()関数を使用して、データセットの体重のカスタム関数を構築します
ecdf_wt <- ecdf(er_data$wt_kg)
# 見やすいプロットデータを作成するための簡単なデータ操作
plot_df <- bind_rows(
in_data = distinct(er_data, wt_kg),
extends = tibble(wt_kg = c(45, 125)),
.id = "type"
) |>
arrange(wt_kg) |>
mutate(
wt_cdf = ecdf_wt(wt_kg),
label = wt_kg %in% c(77.1, 77.2, 77.5),
label_x = case_when(
wt_kg == 77.1 ~ 77.7,
wt_kg == 77.2 ~ 76.8,
wt_kg == 77.5 ~ 78.0,
.default = NA
)
)
# 分類に関連する確率の位置を強調表示した、体重データの経験的CDFを示します
base_plt <- plot_df |>
ggplot(aes(wt_kg, wt_cdf)) +
geom_hline(yintercept = 0:4/4, linetype = "dashed", color = "red") +
geom_point(data = plot_df |> filter(type == "in_data"), size = 2) +
labs(
x = "Weight (kg)",
y = "Empirical CDF",
)
base_plt +
geom_step(direction = "hv") +
annotate(
geom = "label",
x = 50,
y = c(.125, .375, .625, .875),
label = paste0("Q", 1:4)
)
［画像: 体重の経験的CDFと分位数線］
そのプロットを見るだけで、著者がおもちゃのデータセットを作成する際に「スケールに親指を置いた」場所が簡単にわかります。77.2kg（偶然にもサンプル中央値）で発生するECDF関数の急激なジャンプは自然ではありません。それは、ntile()の問題を誇張するために私が仕込んだ部分です。実生活では、その場所でECDFに比較的大きなジャンプが発生することは確かにありますが…ええ、それほど大きくはありません。データセットの人工性は、このプロットではっきりと現れています。しかし、それは重要なことではありません。
ECDFプロットの重要な部分、ECDFが中央値77.2kgを超える領域にズームインしましょう。私たちのプロットは次のようになります。このプロットの各点は、データセットで実際に観測された体重値を示しています。CDFの定義が specifically (\mbox{Pr}(X \leq x)) を参照していることに注意してください。(77.2 \leq 77.2) であることは明らかであり、175個の観測値のうち95個が77.2kg以下であることが観測されているため、その重要な値に関するより正確な記述は、56%の人が77.2kg以下の体重であるということです。つまり、私たちのデータセットの77.2kgの体重を持つすべての人は、体重の56パーセンタイルにいます。8そのため、77.2kgに対応する小さな黒い点は、中央値の上に位置します。
zoom_plt <- base_plt +
geom_label(
data = plot_df |> filter(label == TRUE),
mapping = aes(x = label_x, label = wt_kg)
) +
coord_cartesian(xlim = c(74, 81), ylim = c(.4, .6))
zoom_plt + geom_step(direction = "hv")
［画像: 経験的CDFの中央値付近の拡大図］
この議論に同意する場合、ブレークポイントでのタイ値はすべて「上に上がる」べきです。つまり、77.2kgの体重を持つすべての人は「Q3」に割り当てられるべきであり、「Q2」ではありません。これは、分位数ビン分割に左閉右開区間を使用した場合に期待される動作です。残念ながら、これはRのcut()関数のデフォルトの動作ではありません。これはデフォルトで右閉区間を使用します。この違いに注意してください。
# 右閉区間を使用（cut()のデフォルトに一致）
er_data |>
mutate(wt_bin = cut_quantile(wt_kg, right_closed = TRUE)) |>
summarise(.by = wt_bin, n_subj = n(), min_wt = min(wt_kg), max_wt = max(wt_kg))
# A tibble: 4 × 4
wt_bin n_subj min_wt max_wt
<fct> <int> <dbl> <dbl>
1 Q4 44 85.6 117
2 Q3 33 77.5 85.5
3 Q2 54 67.6 77.2
4 Q1 44 52.3 67.4
そして、これとは異なります。
# 右開区間を使用（ECDFの暗示的な動作、およびsantoku::chop_quantiles()のデフォルトに一致）
er_data |>
mutate(wt_bin = cut_quantile(wt_kg, right_closed = FALSE)) |>
summarise(.by = wt_bin, n_subj = n(), min_wt = min(wt_kg), max_wt = max(wt_kg))
# A tibble: 4 × 4
wt_bin n_subj min_wt max_wt
<fct> <int> <dbl> <dbl>
1 Q4 44 85.6 117
2 Q3 47 77.2 85.5
3 Q2 40 67.6 77.1
4 Q1 44 52.3 67.4
上記の図で黒い点がプロットされている方法に基づくと、2番目のバージョンが取得される分位数ビン分割です。
さて…それは、客観的に正しい答えがあるということでしょうか？規制当局に提出されるデータ分析のために分位数ビンを定義する際には、左開右閉区間を使用すべきであり、話は終わりですか？
LOL。
いいえ。
LMAO、さらに。
これらすべては慣習です。例えば、歴史的な偶然によって、CDFの定義を (\mbox{Pr}(X \leq x)) ではなく (\mbox{Pr}(X < x)) と決定していた場合、前の議論のすべての側面が逆転します。その宇宙では、CDFの意味は微妙に変わりますが、それは依然として同じ数学的オブジェクトを記述しています。どちらの宇宙も道徳的または統計的に優れているわけではありません。分位数ビンを定義するために右開区間または右閉区間を使用しても、本当に構いません。
すべてがどれほど任意であるかを強調するために、上記のプロットコードでgeom_step()にdirection = "hv"という引数を設定したことに注意してください。これは、ステップ関数の描画方法が、まず水平に移動し、次に垂直に移動することを意味します。その慣習を採用すると、体重の累積分布を表すステップ関数は、50パーセンタイルを77.2kgで横切ります。しかし、神または人間の法律は、direction = "hv"が世界についての深い真実を表現しているとは言っていません。
zoom_plt + geom_step(direction = "vh")
［画像: direction = "vh"を使用したステップ関数のプロット］
見よ！体重分布関数は今や77.1kgで中央値を横切ります。(\leq)を(<)に、またはhvをvhに変更することの力はこのようなものです。つまり…あまり大したことではありません。1つの慣習を選んで、それに固執してください。大丈夫でしょう。
［画像: 慣習の任意性を示す図］
エピローグ：SQLのせいにする
この小さな演習を振り返り、dplyr::ntile()が実際に行うことの発見に対する最初の恐怖を考えると、私の不満は、その挙動そのものよりも、その名前にあるように思えます。統計的訓練を受けた人なら誰でも、その名前はそれが分位数ベースのビンを構築するための関数であることを強く示唆します。しかし、それが実際に行うことは、SQLの同様に愚かな挙動を模倣することです。9 SQLの対応物と同様に、ntile()は、より「均等な」分割を生成する場合、同じ値を異なるグループに割り当てることを強制します。しかし、SQLでは、この挙動はさらに奇妙な結果につながります。データベースには行順序という概念が本質的にないため10、NTILEのタイブレーク挙動は単に「非決定論的」と宣言されています。その挙動について本当に知らされることは次のとおりです。
パーティション内の行数がinteger_expressionで割り切れない場合、これは1つのメンバーで異なる2つのサイズのグループを引き起こします。OVER句で指定された順序では、大きいグループが小さいグループの前に来ます。例えば、総行数が53でグループ数が5の場合、最初の3つのグループは11行になり、残りの2つのグループはそれぞれ10行になります。一方、総行数がグループ数で割り切れる場合、行はグループ間で均等に分散されます。例えば、総行数が50でグループが5つある場合、各バケットには10行が含まれます。
NTILEはグループのサイズを指定しますが、そのメンバーシップや値をカテゴリにマッピングするルールは指定しません。不適切に似た名前にもかかわらず、NTILEは意味のある意味で分位数ビン分割のためのツールではありません。11実際、NTILE式のORDER BY句がブレークでタイを生成する場合…まあ、何をするかは誰にもわかりません。それは「非決定論的」であり、この文脈では「ランダム」とは全く同じではありません。ランダムにタイをブレークすることは、ランダムな割り当ての統計的性質はよく理解されているため、はるかに問題が少ないでしょう。しかし、未知であるが分析にとって重要な変数と体系的に関連している可能性のあるルールを使用してタイをブレークすることは？それは恐ろしく悪いです。NTILEとntile()が提供する「迅速で簡単な分割」が有用な文脈があるかもしれませんが、真剣なデータ分析はそうではありません。
そして間違いなく、製薬業界での規制作業にはどちらも適していません。
付録
- おもちゃのデータセットは、exposure_response_data.csvとしてダウンロード可能であり、それを生成するスクリプトもsimulate_data.Rとして利用可能です。
- また、ブレークポイントの分位数推定方法をユーザーが制御でき、タイ処理の3つのオプション（タイが常に上に上がる、タイが常に下に下がるという2つの従来の分位数/カットアプローチ、およびタイの分割をより均等なビンにするが、タイはランダムにのみブレークできる（シード制御付き）という3番目の「ntileライク」アプローチ）を備えた、より慎重に作業された
cut_quantile()関数のバージョンも作成しました。これは、データセットが要求されたビン数をサポートできない場合に自動的に検出し、必要な情報メッセージとエラーメッセージを提供します。スクリプト自体はminicuts.Rとして利用可能であり、簡単な付随するvignetteでその内容を文書化しています。
フットノート
- blah blah blah、性別とジェンダー、性別とジェンダーの分離は統計的過誤である。これは私がそのトピックを再訪する記事ではありません。↩︎
- それに、データセットがより高い精度で体重を報告していたとしても、それは意味がありません。体重は時間とともに変化し、数時間のうちにわずかに変化する可能性があります。被験者が空腹時で裸で、体重計が本当にその精度であったとしても、72.23482キログラムの体重を持つ被験者を報告することは不条理です。↩︎
- エピローグで、これを真の原因であるSQL NTILE関数に結びつけると、SQLが行うことはさらに悪いことが判明します。↩︎
- この記事のどこかで、分位数分割はすべてのデータ値を正確に1つのカテゴリにマッピングする明確に定義されたルールを持っているため、本質的に優れているという考えを検討しました。そして、私はそれを違いの点としてまだ気に入っています。しかし、本当にそれほど違うのでしょうか？ランダムなタイブレークを持つ等サイズルールは、内部のブレークポイントを除くすべてのデータ値に対して決定論的なルールを持っています。そしてそれらのために、観測値を割り当てるための適切に指定された確率的手順があります。したがって、違いは、等サイズ分割ではRNGの状態に注意を払う必要があるということです。全体として…これは分位数分割を等サイズ分割よりも好む理由であるように思えますが、それほど強い理由ではないように感じます。↩︎
- あるいは、洗練された言葉で言うと、ビンを定義する区間が「左閉右開」区間（エッジケースは上に移動）か、「左開右閉」区間（エッジケースは下に移動）かを事前に決定することです。↩︎
- この記事を書いている途中、ランク順統計量から分位数を推定するために使用される統計理論を説明するセクションを含めることを検討しました。なぜなら、
quantile()のドキュメントを読むと、サンプルから分位数推定値を構築する方法が実際にはたくさんあることがわかるからです。しかし、さまざまな異なる方法は、この投稿で説明している問題に対して異なる答えを与えません。本当に違いを生むのは、「上に上がる」か「下に下がる」かの分類ルールの側面です。それはcut()に関連しており、quantile()ではありません。したがって、その理論は重要であり、「タイブレークなし」の原則の論理の一部ですが、低レベルの詳細は投稿には重要ではありません。↩︎ - エピローグを予告すると、
ntile()がこれほど悪いのには歴史的な理由があることが判明します。dplyrパッケージはSQL構文に大きく影響されており、SQLの良い点と悪い点の両方を引き継いでいます。SQLはデータ操作操作を定義するための非常にエレガントな言語ですが、統計的基盤が全くありません。SQLはデータ分析ツールとしてよく参照されますが、データ分析の分析部分ではしばしば非常に劣っています。これは、通常、その挙動に対して非常に堅固な統計的正当性を持つが、プログラミングの基本がしばしば非常に shaky であるRとは正反対です。dplyrパッケージは、SQLへのその基盤がこの特定のタスクには不適切であるため、ntile()に関してこの奇妙な状況に陥っています。↩︎ - えーと…実際、最後の「56パーセンタイル」の部分は完全に真実ではありません。それは、使用する推定値によって重要になる値の1つであることが判明しました。
quantile()に渡すtype引数の値によっては、77.2という低い値から77.5という高い値まで得られます。この投稿で議論している特定の問題については、すべての推定値が中央値体重が77.2kgであると合意しているため、これは無関係です。そのため、タイ処理の問題（私が議論しようとしていること）はこの決定の影響を受けません。そのため、この投稿ではその問題について議論していません。しかし、推定値の選択によって中央値の値が決まる状況に遭遇することが時々あることに注意する価値はあります。例えば、25パーセンタイル（これはビン分割に使用されるブレークポイントの1つです）の位置は、使用する推定値によって67.375kgという低さから67.5kgという高さまであります。統計学は奇妙な学問です。↩︎ dplyr::ntile()のソースコードで内部的にこれを見ることができます。これは、SQL NTILE計算が何をするかの定義を参照しています。設計上、分位数ビンを構築するものではありません。似たようなことをしますが、構築上、NTILEスタイルの等サイズビン分割は分位数分割ビン分割ではありません。↩︎- はい、インデックスはありますが、それは別のことです。↩︎
- 面白がって、プログラマーが統計学者がプログラミングが得意でないために設計したプログラミング言語であるRについて絶えず不平を言うことに気づきます。さて、逆もまた真実であることが判明しました…SQLは、データ分析が得意でないプログラマーによって設計されたデータ分析フレームワークを提供します。
NTILEまたはntile()と呼んでも構いません。それは統計的にばかげており、真剣なデータ分析には適していません。↩︎
再利用
引用
BibTeX引用:
@online{navarro2026,
author = {Navarro, Danielle},
title = {The Trouble with “Ntile()”},
date = {2026-09-20},
url = {https://blog.djnavarro.net/posts/2026-09-20_the-trouble-with-ntile/},
langid = {en}
}
帰属表示については、この作品を次のように引用してください。
Navarro, Danielle. 2026. “The Trouble with ‘Ntile()’.” September 20. https://blog.djnavarro.net/posts/2026-09-20_the-trouble-with-ntile/.