分割表

Cochran-Armitage検定の基礎と製薬統計への応用〜数理的導入からRによる実装まで〜

はじめに

製薬業界における臨床試験では、用量反応関係(dose-response relationship)を検証する場面が頻繁に登場します。特に第II相試験などでは、複数の用量群を設定し、ある二値アウトカム(有効/無効、副作用の有無など)が用量に応じて増加または減少するかを確認することが重要です。
このとき単純なカイ二乗検定を用いると、「群間に差があるかどうか」は分かりますが、差が単調に増加・減少しているかまでは評価できません。そこで登場するのが Cochran-Armitage検定(Cochran-Armitage trend test) です。本検定は、二値アウトカムと順序カテゴリ(例:用量群)の間に線形トレンドが存在するかを検証するための統計手法です。
本稿では、Cochran-Armitage検定について数理的導入、R言語での実装、製薬業界での応用背景について解説していこうと思います。

数理的導入

データの枠組み

以下のような\(2 \times k\)分割表を考えていきます。

用量
反応\(d_{0}\)\(d_{1}\)\(\cdots\)\(d_{k}\)
あり\(x_{0}\)\(x_{1}\)\(\cdots\)\(x_{k}\)\(x_{+}\)
なし\(n_{0} – x_{0}\)\(n_{1} – x_{1}\)\(\cdots\)\(n_{k} – x_{k}\)\(N – x_{+}\)
\(n_{0}\)\(n_{1}\)\(\cdots\)\(n_{k}\)\(N\)
表1: 用量、反応関係の\(2 \times k\)表 (データ)
用量
反応\(d_{0}\)\(d_{1}\)\(\cdots\)\(d_{k}\)
あり\(P_{0}\)\(P_{1}\)\(\cdots\)\(P_{k}\)
なし\(1 – P_{0}\)\(1 – P_{1}\)\(\cdots\)\(1 – P_{k}\)
11\(\cdots\)1
表2: 用量、反応関係の\(2 \times k\)表 (セル確率)

用量に関して順序があり、\(d_{0} \leq d_{0} \leq \cdots \leq d_{K}\)とする。

今表1のデータにk + 1コの独立な二項モデルを想定します。

\[L = \prod_{i=1}^{K} \binom{n_{i}}{x_{i}}P_{i}^{x_{i}}(1 – P_{i})^{n_{i} – x_{i}} — (1) \]

帰無仮説と対立仮説

帰無仮説 H₀: 成功確率 \(P_i\) はすべての群で等しい(トレンドなし)すなわち、
\[P_{0} = P_{1} = \cdots = P_{K}\]

対立仮説 H₁: 成功確率 \(P_i\) は \(x_i\) に沿って単調に増加または減少する(線形トレンドあり)すなわち、
\[P_{0} \leq P_{1} \leq \cdots \leq P_{K}\]
\[P_{0} \ge P_{1} \ge \cdots \ge P_{K}\]

2回微分可能な単調関数hに対して、以下のようにおく。

\[P_{i} = h(\alpha +\beta d_{i}) (i=1,2,…,K) —(2)\]

ただし、\(P_{0} = h(\alpha)\)。hは単調増加関数と呼び、反応関数・反応曲線とも呼ばれます。

この時、用量反応関係を調べる仮説は以下のように表現することもできます。

$$
\begin{eqnarray}
\begin{cases}
H_0: & P_{0} = P_{1} = \cdots = P_{K}\\
H_1: & P_{0} \leq P_{1} \leq \cdots \leq P_{K}
\end{cases}
\end{eqnarray}
$$

$$
\begin{eqnarray}
\leftrightarrow
\begin{cases}
H_0: & \beta = 0\\
H_1: & \beta > 0
\end{cases}
\end{eqnarray}
$$

検定統計量の導出

今回はスコア検定で導出していく。(1)式に(2)式を代入すると、

\[L = \prod_{i=0}^{K}\binom{n_{i}}{x_{i}}[h(\alpha +\beta d_{i})]^{x_{i}}[1 – h(\alpha +\beta d_{i})]^{(n_{i} – x_{i})}\]

上記式の対数尤度関数は

\[log L(x_{i} ;\alpha,\beta) = log C + \sum_{i=0}^{K}x_{i}log h(\alpha +\beta d_{i}) + \sum_{i=0}^{K}(n_{i} – x_{i})log[1 – h(\alpha +\beta d_{i})]\]

ただし、Cは定数であり、\(\alpha,\beta\)に無関係な定数となります。今、\(H_{0}\)の下で\(\beta = 0\)の最尤推定量\(\hat{\alpha}\)は

$$
\begin{eqnarray}
\frac{\partial log L}{\partial \alpha} =\sum_{i=0}^{K}x_{i}\frac{\frac{\partial}{\partial \alpha}h(\alpha)}{h(\alpha)} – \sum_{i=0}^{K}(n_{i} – x_{i})\frac{\frac{\partial}{\partial \alpha}h(\alpha)}{1 – h(\alpha)}(\equiv 0)\\
\leftrightarrow (1 – h(\hat{\alpha}))\sum_{i=0}^{K}x_{i} – h(\hat{\alpha})\sum_{i=0}^{K}(n_{i} – x_{i}) = 0\\
\leftrightarrow h(\hat{\alpha}) = \frac{\sum_{i=0}^{K}x_{i}}{\sum_{i=0}^{K}n_{i}} = \frac{x_{+}}{N} = \hat{h}(\alpha)
\end{eqnarray}
$$

ここで、\(H_0\)の下での\(\beta\)の値を\(\beta_0\)(実際には\(\beta_0 = 0\))とおく。
今、log L を\(\alpha , \beta\)で偏微分していく。

\[{\bf S(x_{i};\alpha , \beta)} =\left[ \begin{array}{c} \frac{\partial}{\partial \alpha}log L \\ \frac{\partial}{\partial \beta}log L \end{array} \right]\]

この時、\(\frac{1}{\sqrt{n}}S(x_{i};\alpha , \beta) \xrightarrow{d} N(0, \frac{1}{n} I(\alpha, \beta))\)となる。
ただし、 \(I(\alpha, \beta)\)はフィッシャー情報行列と呼ばれる。

これより、

\[S(x_{i};\alpha , \beta)^{t} I^{-1}(\alpha, \beta) S(x_{i};\alpha , \beta) \sim \chi^{2}_{(1)}\]

よって求めていく検定統計量Qは

\[Q=S(x_{i};\hat{\alpha}, \beta_{0})^{t} I^{-1}(\hat{\alpha}, \beta_{0}) S(x_{i};\hat{\alpha}, \beta_{0})\]

ここで、

\[{\bf I (\alpha, \beta)} =\left[ \begin{array}{cc} I_{11} (\alpha, \beta) & I_{12} (\alpha, \beta)\\ I_{21} (\alpha, \beta) & I_{22} (\alpha, \beta) \end{array} \right]\]

\[I_{12} (\alpha, \beta) = I_{21} (\alpha, \beta)\]

$$
\begin{eqnarray}
Q &=& \left[\frac{\partial}{\partial \alpha} log L , \frac{\partial}{\partial \beta} log L \right] \left[ \begin{array}{cc}
I_{11} (\hat{\alpha}, \beta_{0}) & I_{12} (\hat{\alpha}, \beta_{0})\\
I_{21} (\hat{\alpha}, \beta_{0}) & I_{22} (\hat{\alpha}, \beta_{0})
\end{array} \right] \left( \begin{array}{c} \frac{\partial}{\partial \alpha} log L \\ \frac{\partial}{\partial \beta} log L \end{array} \right)\\
&=& \frac{(\frac{\partial}{\partial \alpha} log L)^{2} I_{22} (\hat{\alpha}, \beta_{0}) – 2(\frac{\partial}{\partial \alpha} log L) (\frac{\partial}{\partial \beta} log L) + (\frac{\partial}{\partial \beta} log L)^{2}I_{11} (\hat{\alpha}, \beta_{0})}{I_{11} (\hat{\alpha}, \beta_{0}) I_{22} (\hat{\alpha}, \beta_{0}) – I^{2}_{12} (\hat{\alpha}, \beta_{0})}
\end{eqnarray}
$$

ここで\(\frac{\partial}{\partial \alpha} log L\)は\(\alpha = \hat{\alpha}\)の時0となります。

$$
\begin{eqnarray}
Q &=& \frac{(\frac{\partial}{\partial \beta} log L)^{2}I_{11} (\hat{\alpha}, \beta_{0})}{I_{11} (\hat{\alpha}, \beta_{0}) I_{22} (\hat{\alpha}, \beta_{0}) – I^{2}_{12} (\hat{\alpha}, \beta_{0})}\\
&=& \frac{(\frac{\partial}{\partial \beta} log L)^{2}}{I_{22} (\hat{\alpha}, \beta_{0}) – \frac{I^{2}_{12} (\hat{\alpha}, \beta_{0})}{I_{11} (\hat{\alpha}, \beta_{0})}}
\end{eqnarray}
$$

\[\frac{\partial}{\partial \beta} log L = \sum_{i=0}^{K}x_{i}\frac{d_{i}\frac{\partial}{\partial \beta}h(\alpha + \beta d_{i})}{h(\alpha + \beta d_{i})} + \sum_{i=0}^{K}(n_{i} – x_{i})\frac{d_{i}\frac{\partial}{\partial \beta}h(\alpha + \beta d_{i})}{(1 – h(\alpha + \beta d_{i}))}\]

\[\left. \frac{\partial \log L}{\partial \beta} \right|{\alpha = \hat{\alpha}, \, \beta = \beta_{0}}
= \left( \sum_{i=0}^{K} x_{i} \frac{d_{i} h'(\hat{\alpha})}{h(\hat{\alpha})} \right) +
\left( \sum_{i=0}^{K} (n_{i} – x_{i}) \frac{d_{i} h'(\hat{\alpha})}{1 – h(\hat{\alpha})} \right)\]

ただし、\(h^{‘}(\hat{\alpha}) = \frac{\partial}{\partial \beta}h(\alpha + \beta d_{i})\) , \(h(\hat{\alpha}) = \frac{x_{+}}{N}\)より、

$$
\begin{eqnarray}
\text{与式} &=& \frac{x_{+}}{N}h^{‘}(\hat{\alpha}) \sum_{i=0}^{K} x_{i} d_{i} + \frac{N h^{‘}(\hat{\alpha})}{N – x_{+}}\sum_{i=0}^{K}d_{i}(n_{i} – x_{i})\\
&=& \frac{N^{2} h^{‘}(\hat{\alpha})}{x_{+} (N -x_{+})}\left(\frac{N – x_{+}}{N}\sum_{i=0}^{K}x_{i} d_{i} + \frac{x_{+}}{N}d_{i}(n_{i} – x_{i}) \right)\\
&=& \frac{N^{2} h^{‘}(\hat{\alpha})}{x_{+} (N -x_{+})}\left(\sum_{i=0}^{K}x_{i} d_{i} – \frac{x_{+}}{N}d_{i} n_{i} \right) —(3)
\end{eqnarray}
$$

次に、\(I(\alpha, \beta)\)について考えていく。

\[I_{11}(\alpha, \beta) E\left[\frac{\partial^{2} log L}{\partial \alpha^{2}} \right]\]

$$
\begin{eqnarray}
\frac{\partial^{2} log L}{\partial \alpha^{2}} &=& \frac{\partial}{\partial \alpha}\left(\sum_{i=0}^{K} x_{i}\frac{h^{‘}(\alpha + \beta d_{i})}{h(\alpha + \beta d_{i})} – \sum_{i=0}^{K} (n_{i} – x_{i})\frac{h^{‘}(\alpha + \beta d_{i})}{1 – h(\alpha + \beta d_{i})} \right)\\
&=& \sum_{i=0}^{K} x_{i} \frac{h^{”}(\alpha + \beta d_{i})h(\alpha + \beta d_{i}) – h^{‘}(\alpha + \beta d_{i}) \times h^{‘}(\alpha + \beta d_{i})}{(h(\alpha + \beta d_{i}))^{2}}\\
& – & \sum_{i=0}^{K}(n_{i} – x_{i}) \frac{h^{”}(\alpha + \beta d_{i})h(\alpha + \beta d_{i}) – h^{‘}(\alpha + \beta d_{i}) \times h^{‘}(\alpha + \beta d_{i})}{(1 – h(\alpha + \beta d_{i}))^{2}}
\end{eqnarray}
$$

ただし、\(h^{‘}(\alpha + \beta d_{i}) = \frac{\partial}{\partial \alpha}log L\)、
\(h^{”}(\alpha + \beta d_{i}) = \frac{\partial}{\partial \alpha^{2}}log L\)と表記した。

今、\(x_{i} (i = 0,1,…,K) ; \mathop{\perp\!\!\!\!\perp} \sim Bi (n_{i}, P_{i})\)より、

\[E[x_{i}] = n_{i}P_{i} = n_{i} h(\alpha + \beta d_{i})\]

よって、

$$
\begin{eqnarray}
E \left[\frac{\partial^{2} log L}{\partial \alpha^{2}} \right]
&=& \sum_{i=0}^{K}n_{i}\frac{h^{”}(\alpha + \beta d_{i})h(\alpha + \beta d_{i}) – (h^{‘}(\alpha + \beta d_{i}))^{2}}{h(\alpha + \beta d_{i})}\\
&-& \sum_{i=0}^{K}n_{i}\frac{h^{”}(\alpha + \beta d_{i})(1 – h(\alpha + \beta d_{i})) + (h^{‘}(\alpha + \beta d_{i}))^{2}}{1 – h(\alpha + \beta d_{i})}\\
& = & – \sum_{i=0}^{K}n_{i}\frac{(h^{‘}(\alpha + \beta d_{i}))^{2}}{h(\alpha + \beta d_{i})(1 – h(\alpha + \beta d_{i}))}
\end{eqnarray}
$$

\[\therefore I_{11}(\alpha , \beta) = -\sum_{i=0}^{K}n_{i}A_{i}\]

同様にして、
\(I_{12}(\alpha , \beta) = I_{21}(\alpha , \beta) = E \left[\frac{\partial^{2} log L}{\partial \alpha \partial \beta} \right] = – \sum_{i=1}^{K}n_{i}A_{i}\)
\(I_{22}(\alpha , \beta) = E \left[\frac{\partial^{2} log L}{\partial \beta^{2}} \right] = \sum_{i=1}^{K}n_{i}d^{2}_{i}A_{i}\)
\(\because i=0 \text{の時、}d_{0} = 0\)

\(\alpha = \hat{\alpha} , \beta = \beta_{0}\)とすると、

\[A_{i} = \frac{\left(h^{‘}(\hat{\alpha}) \right)^{2}}{h(\hat{\alpha})(1 – h(\hat{\alpha})} = \frac{\left(N h^{‘}(\hat{\alpha}) \right)^{2}}{x_{+}(N – x_{+})}\]

\[I_{22}(\hat{\alpha} , \beta_{0}) – \frac{I_{12}^{2}(\hat{\alpha} , \beta_{0})}{I_{11}(\hat{\alpha} , \beta_{0})} = \frac{\left(N h^{‘}(\hat{\alpha}) \right)^{2}}{x_{+}(N – x_{+})} \left[\sum_{i=1}^{K}n_{i}d^{2}_{i} – \frac{1}{N}(\sum_{i=1}^{K}n_{i}d_{i})^{2} \right] —(4)\]

したがって、検定統計量Qは(3),(4)の式から、\(\frac{(3)^{2}}{(4)}\)であるから、

$$
\begin{eqnarray}
Q &=& \frac{\left( \dfrac{N^{2} h(\hat{\alpha})}{x_{+}(N – x_{+})} \right)^{2}}
{\dfrac{\left(N h'(\hat{\alpha}) \right)^{2}}{x_{+}(N – x_{+})}}
\times
\frac{\left( \sum_{i=1}^{K} x_{i} d_{i} – \frac{x_{+}}{N} d_{i} n_{i} \right)^{2}}
{\sum_{i=1}^{K} n_{i} d_{i}^{2} – \frac{1}{N} \left( \sum_{i=1}^{K} n_{i} d_{i} \right)^{2}} \\
&=& \frac{N^{2}}{x_{+}(N – x_{+})}
\times
\frac{\left( \sum_{i=1}^{K} x_{i} d_{i} – \frac{x_{+}}{N}\sum_{i=1}^{K} d_{i} n_{i} \right)^{2}}
{\sum_{i=1}^{K} n_{i} d_{i}^{2} – \frac{1}{N} \left( \sum_{i=1}^{K} n_{i} d_{i} \right)^{2}}
\end{eqnarray}
$$

\(\frac{x_{+}}{N} = \hat{p}\)とすると、

\[Q = \frac{\left(\sum_{i=0}^{K}x_{i} d_{i} – \hat{p}d_{i} n_{i} \right)^{2}}{\hat{p}(1 – \hat{p})\left[\sum_{i=1}^{K}n_{i}d^{2}_{i} – \frac{1}{N}(\sum_{i=1}^{K}n_{i}d_{i})^{2}\right]}\]

上記式からも分かるように、反応曲線hに関係しないことが重要であり、hがどのような関数でもQに基づく検定は漸近的に最良であると考えられます。

有意水準を\(\alpha\)とすると、片側検定の場合

$$
\begin{eqnarray}
\begin{cases}
P-value = Pr(Z \ge \sqrt{Q}) & < 0.05 \rightarrow H_0\text{を}reject\\
P-value = Pr(Z \ge \sqrt{Q}) & \ge 0.05 \rightarrow H_0\text{を}accept
\end{cases}
\end{eqnarray}
$$

ただし、ZはN(0,1)に従います。

製薬業界での応用背景

用量反応関係の検証

新薬開発において、薬効が用量に応じて増加するか、副作用が用量に応じて増加するかを確認することは極めて重要です。Cochran-Armitage検定は、「用量が上がるにつれて有効率が増加する」といった仮説を直接検証できるため、探索的試験や規制当局への申請資料においても有用です。

規制科学との関係

ICHガイドライン(E9など)では、統計的検証の透明性と妥当性が求められています。Cochran-Armitage検定は、順序性を持つカテゴリ変数と二値アウトカムの関係を評価する標準的手法として広く認知されており、規制当局も理解しやすい検定です。

Rによる実装例

仮想的な臨床試験データを考えます。

  • プラセボ群(0mg): 40例中 5例が有効
  • 低用量群(10mg): 40例中 10例が有効
  • 中用量群(20mg): 40例中 18例が有効
  • 高用量群(40mg): 40例中 25例が有効
データ

#データ作成
dose <- c(0, 10, 20, 40)
success <- c(5, 10, 18, 25)
failure <- c(35, 30, 22, 15)

# 2×k表
tab <- cbind(success, failure)
rownames(tab) <- paste0(dose, “mg”)
tab

Cochran-Armitage検定の実行

# パッケージ読み込み
library(DescTools)

# Cochran-Armitage検定
CochranArmitageTest(tab, scores = dose, alternative = “increasing”)

この結果から、p値が有意に小さいため「用量が増えるにつれて有効率が増加する」というトレンドが統計的に支持されます。

実務上の注意点

スコアの設定: 通常は用量をそのままスコアに用いるが、対数変換や等間隔スコアを使う場合もある。

片側検定か両側検定か: 製薬業界では「増加方向のみ」を仮説とする片側検定が多い。

サンプルサイズ: 小規模試験では漸近近似が不十分になるため、正確検定(exact test)を検討する。

他の手法との比較
  • カイ二乗検定: 群間差は検出できるが、順序性を利用しないため検出力が低下する。
  • ロジスティック回帰: Cochran-Armitage検定はロジスティック回帰のスコア検定と等価であり、より一般化した枠組みで理解できる。
  • Jonckheere-Terpstra検定: こちらは連続アウトカムや順序尺度に適用可能で、二値データに特化したCochran-Armitage検定と使い分ける。

まとめ

今回はCochran-Armitage検定について数理的導入、R言語での実装、製薬業界での応用背景について解説していきました。(数式が多くなり失礼いたしました。)
Cochran-Armitage検定は、二値アウトカムと順序カテゴリの間に線形トレンドが存在するかを評価する有効な手法です。製薬業界では用量反応関係の検証に広く用いられ、カイ二乗検定よりも検出力が高い点が特徴です。RではCochranArmitageTest関数で容易に実装でき、規制当局への説明責任にも適した透明性の高い解析法といえます。

ABOUT ME
tomokichi
外資系製薬会社で生物統計家として働ている1児のパパ。生物統計家とは何か、どのようなスキルが必要か、何を行っているのかを共有していきたいと思っております!生物統計に関する最新情報を皆様にお届けすべく、日々奮闘中です。趣味は筋トレ、温泉巡り、家族と散歩。