ある教育プログラムの効果を検証するために、2つのクラスにそれぞれ異なる教授法を適用し、学期末の試験スコアを比較したいとします。しかし、クラスの学生は学期前の学力が異なるかもしれません。もし一方のクラスに学力の高い学生が偏っていたら、試験スコアの差は教授法の効果なのか、元々の学力差なのか区別できません。
このような交絡変数（confounding variable）の影響を統計的に除去し、処理効果を正しく評価する手法が共分散分析（Analysis of Covariance, ANCOVA）です。ANCOVAは分散分析（ANOVA）と回帰分析を組み合わせた手法であり、群間比較を行いながら連続的な共変量の影響を調整します。
ANCOVAを理解すると、以下のような応用が可能です。
- 臨床試験: 治療効果を評価する際に、ベースラインの測定値や年齢などの影響を調整します
- 教育研究: 教授法の比較で、事前テストのスコア（学力差）を共変量として統制します
- 社会科学: 所得の地域差を分析する際に、教育水準や年齢構成の違いを調整します
- 農学: 品種間の収量比較で、土壌の肥沃度の差を共変量として除去します
本記事の内容
- ANCOVAの直感的な理解 — なぜ共変量の調整が必要か
- ANCOVAのモデルと数学的定式化
- パラメータの最小二乗推定とその導出
- 調整平均の導出と解釈
- F統計量の導出と検定の手順
- 回帰の平行性の仮定と検証
- ANOVAとの比較と検出力の理論的分析
- 複数共変量への拡張
- 実践的なガイドラインと注意点
- Pythonによる実装と可視化
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
共変量の調整 — なぜ必要か
交絡の問題
2つの教授法（A, B）を比較する例に戻りましょう。各群の事前テストスコア（共変量）と事後テストスコア（目的変数）の関係を散布図で描くと、次のような状況が起こり得ます。
- 教授法A群の学生は事前テストの平均が70点で、事後テストの平均が80点
- 教授法B群の学生は事前テストの平均が60点で、事後テストの平均が75点
単純にANOVAで事後テストスコアを比較すると、「教授法Aの方が5点優れている」と結論するかもしれません。しかし、教授法A群は元々学力が高い学生で構成されているため、この5点の差は教授法の効果ではなく、元の学力差を反映しているだけかもしれません。
この状況を数値で具体的に見てみましょう。事前テストと事後テストの間に正の相関（たとえば $r = 0.7$）があるとします。事前テストの群間差は $70 – 60 = 10$ 点です。回帰係数を $\beta = 0.6$ と推定すると、事前テストの差10点は事後テストに $0.6 \times 10 = 6$ 点の差を生み出すと予測されます。つまり、観測された5点の差のうち6点が学力差の反映だとすると、実際にはB群の教授法の方が $5 – 6 = -1$ 点、すなわち1点ほど優れている可能性があります。
ANCOVAは、事前テストスコアの影響を統計的に除去した上で事後テストスコアを比較します。つまり、「もし両群の事前テストの平均が同じだったら、事後テストの差はどうなるか」を推定します。
調整の仕組み — 回帰を用いた除去
ANCOVAの調整メカニズムは、回帰分析と分散分析の組み合わせとして理解できます。
- ステップ1: 目的変数 $Y$ の共変量 $X$ への回帰直線を当てはめる
- ステップ2: 回帰直線からの残差（共変量で説明できない部分）を取り出す
- ステップ3: この残差について群間の差を検定する
直感的には、「共変量の影響を取り除いた後に残る群間差」を見ています。これにより、共変量の群間差による交絡が除去されます。
たとえば、事前テストの点数が高い学生ほど事後テストの点数も高い傾向があるとします。回帰係数が $\beta = 0.6$ ならば、事前テストが10点高い群は事後テストでも $0.6 \times 10 = 6$ 点高くなると予測されます。ANCOVAはこの6点を「共変量の影響」として差し引き、残りの差を「処理効果」として評価します。
このプロセスは以下のように数式で表現できます。個体 $(i, j)$ の「共変量調整済みスコア」は、
$$ Y_{ij}^{\text{adj}} = Y_{ij} – \hat{\beta}(X_{ij} – \bar{X}_{\cdot\cdot}) $$
として計算されます。この調整済みスコアの群平均が調整平均であり、群間差の検定がANCOVAのF検定に対応します。ただし実際の推定はもう少し複雑で、群効果と回帰係数を同時に推定する必要があります。
ANCOVAが本質的にやっていること — 3つの見方
ANCOVAの本質を理解するために、3つの同値な見方を紹介します。
見方1: 残差の分散分析
目的変数 $Y$ と共変量 $X$ それぞれを群のダミー変数に回帰し、その残差 $e_Y$ と $e_X$ を求めます。次に $e_Y$ を $e_X$ に回帰し、その残差を群間で比較します。この「二重の残差化」が、共変量と群効果の両方を調整する仕組みです。この見方はフリッシュ・ウォーの定理（Frisch-Waugh-Lovell theorem）に基づいており、部分回帰の考え方と直結しています。
具体的には、まず $Y$ を群のダミー変数に回帰して残差 $e_Y = Y – \hat{Y}_{\text{group}}$ を得ます。同様に $X$ を群のダミー変数に回帰して残差 $e_X = X – \hat{X}_{\text{group}}$ を得ます。これらの残差はそれぞれ「群の影響を除いたY」「群の影響を除いたX」を表しています。そして $e_Y$ を $e_X$ に回帰したときの回帰係数は、フルモデルにおける共変量の回帰係数 $\hat{\beta}$ と一致します。
見方2: 共変量を追加した一般線形モデル
ANOVAのモデルに回帰項を追加し、一般線形モデルの枠組みでF検定を行います。群効果と共変量の効果を同時に推定し、共変量を制御した上での群効果を評価します。これは最も直接的な定式化であり、後述するモデルの定義はこの見方に対応しています。一般線形モデルの枠組みでは、群のダミー変数と共変量がともに説明変数として含まれ、通常の最小二乗法でパラメータが推定されます。
見方3: 条件付き比較
共変量の値を固定した場合の群間差を推定します。回帰の平行性の仮定（後述）が成り立てば、固定する共変量の値によらず群間差は一定です。全体平均で固定した値が調整平均です。この見方は最も直感的で、「もしすべての被験者が同じ事前テストスコアだったら」という反事実的な問いに対応しています。
これらの見方はすべて同じ結果を与えますが、場面によって直感的な理解の助けになります。第1の見方は偏相関の概念との関連を示し、第2の見方は統計ソフトウェアの実装に直結し、第3の見方は結果の解釈に最も役立ちます。
この直感を数学的に定式化しましょう。
ANCOVAのモデル
モデルの定義
$k$ 群の比較で、共変量が1つの場合のANCOVAモデルは次の通りです。
$$ \begin{equation} Y_{ij} = \mu + \tau_j + \beta(X_{ij} – \bar{X}_{\cdot\cdot}) + \varepsilon_{ij} \end{equation} $$
- $Y_{ij}$: 群 $j$ の $i$ 番目の観測値（目的変数）
- $X_{ij}$: 対応する共変量の値
- $\mu$: 全体平均
- $\tau_j$: 群 $j$ の効果（$\sum_j n_j \tau_j = 0$）
- $\beta$: 共変量の回帰係数（全群共通）
- $\bar{X}_{\cdot\cdot}$: 共変量の全体平均
- $\varepsilon_{ij} \overset{\text{i.i.d.}}{\sim} N(0, \sigma^2)$
このモデルを通常のANOVAモデル $Y_{ij} = \mu + \tau_j + \varepsilon_{ij}$ と比較すると、違いは $\beta(X_{ij} – \bar{X}_{\cdot\cdot})$ の項が追加されていることです。この項が共変量の影響を吸収し、残差の分散を小さくします。
共変量を全体平均からの偏差 $X_{ij} – \bar{X}_{\cdot\cdot}$ で表している理由は、$\mu$ と $\tau_j$ の解釈を明確にするためです。$X_{ij} = \bar{X}_{\cdot\cdot}$（共変量が全体平均に等しい）のとき、$Y_{ij}$ の期待値は $\mu + \tau_j$ となり、$\tau_j$ がそのまま群効果を表します。もし中心化しなければ $\mu$ の解釈が「$X=0$ における全体平均」となり、$X=0$ に意味がない場合（たとえばテストスコアが0点の場合）に解釈が困難になります。
このモデルの幾何学的な意味を考えましょう。$Y$ を $X$ の関数として描くと、各群は傾きが同じで切片が異なる直線に対応します。群 $j$ の回帰直線は、
$$ E[Y_{ij} \mid X_{ij}] = (\mu + \tau_j) + \beta(X_{ij} – \bar{X}_{\cdot\cdot}) $$
であり、$k$ 本の平行な直線が得られます。ANCOVAの帰無仮説 $\tau_1 = \tau_2 = \cdots = \tau_k = 0$ は、「これらの平行線がすべて一本の直線に重なる」ことを意味します。
行列表現
一般線形モデルの枠組みでは、ANCOVAモデルは次のように表現されます。
$$ \bm{Y} = \bm{X}_{\text{design}} \bm{\theta} + \bm{\varepsilon} $$
2群（$k=2$）の場合を具体的に書くと、セル平均モデルの形式では、
$$ \begin{pmatrix} Y_{11} \\ Y_{21} \\ \vdots \\ Y_{n_1,1} \\ Y_{12} \\ \vdots \\ Y_{n_2,2} \end{pmatrix} = \begin{pmatrix} 1 & 0 & X_{11} – \bar{X}_{\cdot\cdot} \\ 1 & 0 & X_{21} – \bar{X}_{\cdot\cdot} \\ \vdots & \vdots & \vdots \\ 1 & 0 & X_{n_1,1} – \bar{X}_{\cdot\cdot} \\ 0 & 1 & X_{12} – \bar{X}_{\cdot\cdot} \\ \vdots & \vdots & \vdots \\ 0 & 1 & X_{n_2,2} – \bar{X}_{\cdot\cdot} \end{pmatrix} \begin{pmatrix} \mu_1 \\ \mu_2 \\ \beta \end{pmatrix} + \bm{\varepsilon} $$
ここで $\mu_j = \mu + \tau_j$ はセル平均パラメータです。この行列表現を使えば、ANCOVAのパラメータ推定を通常の最小二乗法 $\hat{\bm{\theta}} = (\bm{X}^{\top}\bm{X})^{-1}\bm{X}^{\top}\bm{Y}$ で行うことができます。
一般の $k$ 群の場合には、効果コーディング（制約 $\sum_j n_j \tau_j = 0$）を使ったデザイン行列は次のようになります。
$$ \bm{X}_{\text{design}} = \begin{pmatrix} \bm{1}_N & \bm{D} & \bm{x} – \bar{X}_{\cdot\cdot}\bm{1}_N \end{pmatrix} $$
ここで $\bm{1}_N$ は $N$ 次元のすべての要素が1のベクトル、$\bm{D}$ は群のダミー変数からなる $N \times (k-1)$ 行列（基準群を省略）、$\bm{x}$ は共変量ベクトルです。パラメータベクトルは $\bm{\theta} = (\mu, \tau_1, \dots, \tau_{k-1}, \beta)^{\top}$ です。
帰無仮説と対立仮説
ANCOVAの帰無仮説は、共変量を調整した後に群効果がないということです。
$$ H_0: \tau_1 = \tau_2 = \cdots = \tau_k = 0 $$
$$ H_1: \text{少なくとも1つの } \tau_j \neq 0 $$
この検定は、2つのネストされたモデルの比較として行います。
フルモデル（制約なし）: $Y_{ij} = \mu + \tau_j + \beta(X_{ij} – \bar{X}_{\cdot\cdot}) + \varepsilon_{ij}$
制約モデル（$\tau_j = 0$）: $Y_{ij} = \mu + \beta(X_{ij} – \bar{X}_{\cdot\cdot}) + \varepsilon_{ij}$
フルモデルと制約モデルの残差平方和の差が十分大きければ、群効果が有意であると判断します。
ここで注目すべきは、帰無仮説が「共変量を調整した上での群間差なし」であり、通常のANOVAの帰無仮説とは微妙に異なる点です。ANOVAでは「群の平均が等しい」ですが、ANCOVAでは「共変量が等しい条件の下で群の平均が等しい」を検定しています。
パラメータ推定の具体的な手順を見ていきましょう。
パラメータの最小二乗推定
正規方程式の導出
ANCOVAモデルのパラメータ $\mu, \tau_1, \dots, \tau_k, \beta$ を最小二乗法で推定します。最小化する目的関数は、
$$ S(\mu, \tau_1, \dots, \tau_k, \beta) = \sum_{j=1}^{k} \sum_{i=1}^{n_j} \left[ Y_{ij} – \mu – \tau_j – \beta(X_{ij} – \bar{X}_{\cdot\cdot}) \right]^2 $$
各パラメータについて偏微分をゼロに置くと、正規方程式が得られます。
$$ \frac{\partial S}{\partial \mu} = -2 \sum_{j} \sum_{i} \left[ Y_{ij} – \mu – \tau_j – \beta(X_{ij} – \bar{X}_{\cdot\cdot}) \right] = 0 $$
$$ \frac{\partial S}{\partial \tau_j} = -2 \sum_{i} \left[ Y_{ij} – \mu – \tau_j – \beta(X_{ij} – \bar{X}_{\cdot\cdot}) \right] = 0 \quad (j = 1, \dots, k) $$
$$ \frac{\partial S}{\partial \beta} = -2 \sum_{j} \sum_{i} (X_{ij} – \bar{X}_{\cdot\cdot}) \left[ Y_{ij} – \mu – \tau_j – \beta(X_{ij} – \bar{X}_{\cdot\cdot}) \right] = 0 $$
推定量の解法
第1の方程式から、制約 $\sum_j n_j \tau_j = 0$ の下で、
$$ \hat{\mu} = \bar{Y}_{\cdot\cdot} $$
が得られます（$\sum_j \sum_i (X_{ij} – \bar{X}_{\cdot\cdot}) = 0$ を利用）。
第2の方程式から、
$$ \hat{\mu} + \hat{\tau}_j + \hat{\beta}(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot}) = \bar{Y}_{\cdot j} $$
したがって、
$$ \hat{\tau}_j = (\bar{Y}_{\cdot j} – \bar{Y}_{\cdot\cdot}) – \hat{\beta}(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot}) $$
この式は重要な意味を持ちます。群効果の推定量 $\hat{\tau}_j$ は、ANOVAにおける群効果 $(\bar{Y}_{\cdot j} – \bar{Y}_{\cdot\cdot})$ から、共変量の群間差による寄与 $\hat{\beta}(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot})$ を差し引いたものです。まさに「共変量の影響を除去した群効果」です。
第3の方程式から、$\hat{\beta}$ の推定に進みます。群効果を含むモデルにおける $\hat{\beta}$ は、群内平方和と共積和から求められます。具体的には、以下の量を定義します。
$$ SS_{XY, \text{within}} = \sum_{j} \sum_{i} (X_{ij} – \bar{X}_{\cdot j})(Y_{ij} – \bar{Y}_{\cdot j}) $$
$$ SS_{X, \text{within}} = \sum_{j} \sum_{i} (X_{ij} – \bar{X}_{\cdot j})^2 $$
これらを使って、共変量の回帰係数は次のように推定されます。
$$ \hat{\beta} = \frac{SS_{XY, \text{within}}}{SS_{X, \text{within}}} $$
ここで「群内」とは、各群の平均からの偏差で計算した量を意味します。これは重要な結果です。ANCOVAの回帰係数は、群内変動のみに基づいて推定されます。群間変動は含みません。つまり、$\hat{\beta}$ は「同じ群の中で共変量が1単位変化したとき、目的変数がどれだけ変化するか」を表しており、群の違いによる変動は排除されています。
この推定量が群間変動を含まない理由は、デザイン行列の構造に由来します。群のダミー変数が含まれているため、共変量の効果は「群効果で説明できない部分」から推定されます。これはフリッシュ・ウォーの定理の直接的な帰結です。
$\hat{\beta}$ の導出の詳細
第3の正規方程式を展開します。
$$ \sum_{j} \sum_{i} (X_{ij} – \bar{X}_{\cdot\cdot}) Y_{ij} = \hat{\mu} \sum_{j} \sum_{i} (X_{ij} – \bar{X}_{\cdot\cdot}) + \sum_{j} \hat{\tau}_j \sum_{i} (X_{ij} – \bar{X}_{\cdot\cdot}) + \hat{\beta} \sum_{j} \sum_{i} (X_{ij} – \bar{X}_{\cdot\cdot})^2 $$
左辺の第1項は $\sum_{j} \sum_{i} (X_{ij} – \bar{X}_{\cdot\cdot}) = 0$ なので消えます。$\hat{\tau}_j$ を代入し、$X_{ij} – \bar{X}_{\cdot\cdot} = (X_{ij} – \bar{X}_{\cdot j}) + (\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot})$ と分解すると、整理の結果として、
$$ \hat{\beta} = \frac{\sum_{j} \sum_{i} (X_{ij} – \bar{X}_{\cdot j})(Y_{ij} – \bar{Y}_{\cdot j})}{\sum_{j} \sum_{i} (X_{ij} – \bar{X}_{\cdot j})^2} = \frac{SS_{XY, \text{within}}}{SS_{X, \text{within}}} $$
を得ます。これは各群のデータについて個別に計算した回帰係数 $\hat{\beta}_j = SS_{XY,j} / SS_{X,j}$ のプール推定量です。各群の $SS_X$ で重みづけた加重平均として解釈できます。
$$ \hat{\beta} = \frac{\sum_j SS_{XY,j}}{\sum_j SS_{X,j}} = \frac{\sum_j SS_{X,j} \hat{\beta}_j}{\sum_j SS_{X,j}} $$
回帰係数が群間で共通であるという仮定（平行性の仮定）の下では、このプール推定量が最も効率的な推定を与えます。
それでは、推定されたパラメータに基づいて調整平均を定義しましょう。
調整平均の導出
調整平均の定義
ANCOVAの最も重要な出力の一つが調整平均（adjusted mean, least squares mean, estimated marginal mean）です。群 $j$ の調整平均は次のように定義されます。
$$ \begin{equation} \bar{Y}_{j, \text{adj}} = \hat{\mu} + \hat{\tau}_j = \bar{Y}_{\cdot j} – \hat{\beta}(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot}) \end{equation} $$
この式の意味は明快です。群 $j$ の観測された平均 $\bar{Y}_{\cdot j}$ から、共変量の群間差 $(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot})$ による寄与 $\hat{\beta}(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot})$ を引いています。
もし群 $j$ の共変量の平均 $\bar{X}_{\cdot j}$ が全体平均 $\bar{X}_{\cdot\cdot}$ より大きければ、$\hat{\beta} > 0$ のとき調整平均は観測平均よりも下方修正されます。逆に、共変量の平均が低い群は上方修正されます。
具体例で確認しましょう。2群の場合で、全体平均 $\bar{X}_{\cdot\cdot} = 65$、$\hat{\beta} = 0.6$ とします。
- 群A: $\bar{X}_{\cdot 1} = 70$, $\bar{Y}_{\cdot 1} = 80$ → $\bar{Y}_{1, \text{adj}} = 80 – 0.6 \times (70 – 65) = 80 – 3 = 77$
- 群B: $\bar{X}_{\cdot 2} = 60$, $\bar{Y}_{\cdot 2} = 75$ → $\bar{Y}_{2, \text{adj}} = 75 – 0.6 \times (60 – 65) = 75 + 3 = 78$
観測された平均では群Aが5点高かったのに、調整平均では群Bが1点高くなりました。群Aの高い平均は事前スコアの高さによるものであり、教授法の効果としては群Bの方が上回っていたのです。
この「逆転」はANCOVAの特徴的な結果であり、共変量の調整がいかに重要であるかを示しています。シンプソンのパラドクスとも関連する現象です。
調整平均の標準誤差
調整平均の標準誤差は、通常の群平均の標準誤差より複雑です。群 $j$ の調整平均の分散を導出しましょう。
$\bar{Y}_{j, \text{adj}} = \bar{Y}_{\cdot j} – \hat{\beta}(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot})$ について、$\bar{Y}_{\cdot j}$ と $\hat{\beta}$ はともに確率変数です。$\hat{\beta}$ の分散は $\text{Var}(\hat{\beta}) = \sigma^2 / SS_{X, \text{within}}$ です。$\bar{Y}_{\cdot j}$ と $\hat{\beta}$ の共分散も考慮すると、
$$ \text{Var}(\bar{Y}_{j, \text{adj}}) = \sigma^2 \left(\frac{1}{n_j} + \frac{(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot})^2}{SS_{X, \text{within}}}\right) $$
したがって標準誤差は、
$$ SE(\bar{Y}_{j, \text{adj}}) = \sqrt{MS_{\text{res}} \left(\frac{1}{n_j} + \frac{(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot})^2}{SS_{X, \text{within}}}\right)} $$
ここで $MS_{\text{res}} = SS_{\text{res, full}} / (N – k – 1)$ はフルモデルの残差平均平方であり、$\sigma^2$ の不偏推定量です。
この式の第2項 $(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot})^2 / SS_{X, \text{within}}$ は、共変量の群平均が全体平均から離れるほど調整平均の不確実性が増すことを反映しています。共変量の群間差が大きいほど、外挿的な調整が必要になり、推定精度が下がるのです。これは回帰分析において予測値が平均から離れるほど予測区間が広がることと本質的に同じ現象です。
特に、ランダム化実験では $\bar{X}_{\cdot j} \approx \bar{X}_{\cdot\cdot}$ が期待されるため、第2項はほぼゼロとなり、標準誤差は $\sqrt{MS_{\text{res}} / n_j}$ に近づきます。
調整平均間の差の検定
2群の調整平均の差 $\bar{Y}_{1, \text{adj}} – \bar{Y}_{2, \text{adj}}$ の分散は、
$$ \text{Var}(\bar{Y}_{1, \text{adj}} – \bar{Y}_{2, \text{adj}}) = \sigma^2 \left(\frac{1}{n_1} + \frac{1}{n_2} + \frac{(\bar{X}_{\cdot 1} – \bar{X}_{\cdot 2})^2}{SS_{X, \text{within}}}\right) $$
第3項 $(\bar{X}_{\cdot 1} – \bar{X}_{\cdot 2})^2 / SS_{X, \text{within}}$ は2群の共変量の平均差に起因する付加的な不確実性です。t統計量は、
$$ t = \frac{\bar{Y}_{1, \text{adj}} – \bar{Y}_{2, \text{adj}}}{\sqrt{MS_{\text{res}} \left(\frac{1}{n_1} + \frac{1}{n_2} + \frac{(\bar{X}_{\cdot 1} – \bar{X}_{\cdot 2})^2}{SS_{X, \text{within}}}\right)}} $$
帰無仮説のもとで $t \sim t(N – k – 1)$ です。分母の第3項は、共変量の群間差による調整の不確実性を反映した「ペナルティ」です。2群の共変量の平均が離れているほど分母が大きくなり、有意な差を検出しにくくなります。この性質は、共変量のバランスが良い（ランダム化がうまくいっている）ほどANCOVAの効率が良いことを意味しています。
F統計量の導出と平方和の分解
ネストされたモデルの比較
ANCOVAのF統計量は、共変量を調整した後の群間差を検定します。
フルモデルの残差平方和を $SS_{\text{res, full}}$、制約モデル（群効果なし、共変量のみ）の残差平方和を $SS_{\text{res, reduced}}$ とすると、
$$ \begin{equation} F = \frac{(SS_{\text{res, reduced}} – SS_{\text{res, full}}) / (k – 1)}{SS_{\text{res, full}} / (N – k – 1)} \end{equation} $$
帰無仮説のもとで $F \sim F(k-1, N-k-1)$ です。通常のANOVAの分母の自由度は $N – k$ ですが、ANCOVAでは共変量のパラメータ $\beta$ の推定に1自由度を使うため $N – k – 1$ になります。
各残差平方和の具体的な計算
フルモデルの残差平方和は、
$$ SS_{\text{res, full}} = SS_{Y, \text{within}} – \frac{(SS_{XY, \text{within}})^2}{SS_{X, \text{within}}} $$
ここで $SS_{Y, \text{within}} = \sum_j \sum_i (Y_{ij} – \bar{Y}_{\cdot j})^2$ は $Y$ の群内平方和です。右辺第2項は「群内回帰で説明される変動」であり、群内変動から回帰で説明できる部分を差し引いたものがフルモデルの残差です。
制約モデル（群効果なし）の残差平方和は、
$$ SS_{\text{res, reduced}} = SS_{Y, \text{total}} – \frac{(SS_{XY, \text{total}})^2}{SS_{X, \text{total}}} $$
ここで $SS_{Y, \text{total}} = \sum_j \sum_i (Y_{ij} – \bar{Y}_{\cdot\cdot})^2$、$SS_{XY, \text{total}} = \sum_j \sum_i (X_{ij} – \bar{X}_{\cdot\cdot})(Y_{ij} – \bar{Y}_{\cdot\cdot})$、$SS_{X, \text{total}} = \sum_j \sum_i (X_{ij} – \bar{X}_{\cdot\cdot})^2$ はそれぞれ全体の平方和と共積和です。
したがって、F統計量の分子（共変量調整済み群間平方和）は、
$$ SS_{\text{group, adj}} = SS_{\text{res, reduced}} – SS_{\text{res, full}} $$
$$ = \left( SS_{Y, \text{total}} – \frac{(SS_{XY, \text{total}})^2}{SS_{X, \text{total}}} \right) – \left( SS_{Y, \text{within}} – \frac{(SS_{XY, \text{within}})^2}{SS_{X, \text{within}}} \right) $$
平方和の分解とType I/II/III
ANCOVAにおける平方和の分解を詳しく見ましょう。全体の変動を以下のように分解します。
$$ SS_{\text{total, adjusted}} = SS_{\text{group, adj}} + SS_{\text{covariate, adj}} + SS_{\text{residual}} $$
ここで、
- $SS_{\text{group, adj}}$: 共変量を調整した後の群間平方和
- $SS_{\text{covariate, adj}}$: 群効果を調整した後の共変量による平方和
- $SS_{\text{residual}}$: 残差平方和
注意すべき点は、この分解における各項は一般に直交（独立）ではないということです。群効果と共変量は相関しうるため、平方和の分解はモデルに変数を投入する順序に依存します。これはType I/II/III の平方和の区別として知られています。
Type I（逐次型）: 変数を順番に投入し、各変数の追加的な貢献を評価します。共変量を先に投入するか群効果を先に投入するかで結果が異なります。順序 $X \to \tau$ とすると、$SS_X$ は共変量単独の効果、$SS_{\tau \mid X}$ は共変量を制御した後の群効果です。通常のANCOVAではこの順序が適切です。
Type II: 各変数の効果を、その変数を含まないが他のすべての主効果を含むモデルと比較して評価します。つまり群効果は共変量を含むモデルとの比較、共変量は群効果を含むモデルとの比較で評価されます。交互作用がない場合にType IIIと一致します。
Type III: 各変数の効果を、他のすべての変数を含むモデルと比較して評価します。不均衡データでも解釈が明確です。ANCOVAでは通常 Type III が推奨されます。
均等なサンプルサイズ（$n_1 = n_2 = \cdots = n_k$）のとき、ランダム化実験であれば $\bar{X}_{\cdot j} \approx \bar{X}_{\cdot\cdot}$ が成り立ち、群効果と共変量が（近似的に）直交するため、Type I, II, III はほぼ同じ結果を与えます。不均衡データや観察研究ではType IIIの使用が強く推奨されます。
ANCOVAの仮定
回帰の平行性（等傾斜）
ANCOVAの最も重要な仮定は、回帰直線が全群で平行であることです。つまり、共変量の回帰係数 $\beta$ が群によらず同じであることを要求します。
$$ \beta_1 = \beta_2 = \cdots = \beta_k = \beta \quad \text{（平行性の仮定）} $$
この仮定が成り立たない場合（群ごとに傾きが異なる場合）、共変量の値によって群間差の大きさが変わってしまい、「共変量の影響を除去した後の群間差」という概念が単一の値として意味を持たなくなります。
たとえば、教授法Aでは事前テストが高い学生ほど事後テストの改善が大きく（$\beta_A = 0.8$）、教授法Bではそのような傾向が弱い（$\beta_B = 0.3$）場合を考えましょう。事前テスト60点の学生では教授法Bが優れ、事前テスト80点の学生では教授法Aが優れるかもしれません。「どちらの教授法が優れているか」は事前テストの値によって答えが変わるのです。
平行性の検定
平行性の仮定を検定するには、群×共変量の交互作用項を含むモデルを当てはめます。
$$ Y_{ij} = \mu + \tau_j + \beta_j(X_{ij} – \bar{X}_{\cdot\cdot}) + \varepsilon_{ij} $$
交互作用の検定は以下の帰無仮説・対立仮説で行います。
$$ H_0: \beta_1 = \beta_2 = \cdots = \beta_k \quad \text{vs.} \quad H_1: \text{少なくとも2つの } \beta_j \text{ が異なる} $$
フルモデル（交互作用あり）: 各群が独自の傾き $\beta_j$ を持ちます。自由度は $N – 2k$ です（$k$ 個の切片と $k$ 個の傾きで $2k$ 個のパラメータ）。
制約モデル（平行性の仮定）: 共通の傾き $\beta$ を持ちます。自由度は $N – k – 1$ です。
交互作用のF検定統計量は、
$$ F_{\text{interaction}} = \frac{(SS_{\text{res, parallel}} – SS_{\text{res, interaction}}) / (k – 1)}{SS_{\text{res, interaction}} / (N – 2k)} $$
この検定が有意（$p < 0.05$）であれば、平行性の仮定は棄却され、通常のANCOVAを適用すべきではありません。ただし、この検定は検出力が比較的低いため、非有意であっても平行性が完全に保証されるわけではありません。散布図による視覚的確認も併用することが推奨されます。
平行性が成り立たない場合の対処法
平行性の仮定が棄却された場合、いくつかの代替手法があります。
1. Johnson-Neyman法: 「どの共変量の範囲で群間差が有意であるか」を特定します。具体的には、群間差のt統計量が $\pm t_{\alpha/2}$ に等しくなる共変量の値（交差点）を求め、その値を境に有意な領域と非有意な領域を特定します。
2群の場合、群間差を共変量 $X$ の関数として $d(X) = (\hat{\mu}_1 – \hat{\mu}_2) + (\hat{\beta}_1 – \hat{\beta}_2)(X – \bar{X}_{\cdot\cdot})$ と表すと、分散は $X$ の2次関数になり、$|d(X)| / SE(d(X)) = t_{\alpha/2, N-2k}$ を解いて2つの境界値を求めます。
2. 特定の共変量値での比較（pick-a-point法）: 共変量の特定の値（たとえば第1四分位数、中央値、第3四分位数）での群間差をそれぞれ推定・検定します。傾きが異なる場合、群間差は共変量の値によって変わるため、複数の点での比較が適切です。
3. ノンパラメトリック手法: 共変量と目的変数の関係が非線形である可能性がある場合、スプライン回帰やカーネル法を用いた一般化ANCOVAが利用できます。一般化加法モデル（GAM）を使えば、$E[Y] = \mu + \tau_j + f(X)$ のように非線形な共変量効果をモデル化できます。
その他の仮定と診断法
ANCOVAには平行性以外にもいくつかの重要な仮定があります。それぞれの仮定、違反の影響、診断法を整理します。
1. 共変量の測定誤差なし: 共変量に測定誤差があると、調整が不完全になります。真の回帰係数 $\beta$ に対して、測定誤差があると推定値は $\beta \cdot \rho_{XX’}$ に縮小されます（$\rho_{XX’}$ は信頼性係数）。これを回帰希釈バイアス（regression dilution bias）と呼びます。信頼性が0.8の共変量を使うと、調整の20%が失われます。対処法として、測定誤差の大きさが既知の場合はSimex法やDeming回帰を利用できます。
2. 残差の正規性と等分散性: ANOVAと同じ仮定です。正規性からの逸脱に対してはロバストですが、等分散性の違反は深刻な影響を与えうるため、レヴェン検定や残差プロットで事前に確認すべきです。等分散性が満たされない場合は、ウェルチ型の修正やブートストラップ法を検討します。
3. 共変量と目的変数の線形関係: 非線形の場合は残差プロットで確認できます。多項式項（$X^2, X^3$）を追加するか、一般化加法モデル（GAM）を検討します。
4. 共変量と群の独立性（ランダム化実験の場合）: ランダム化実験では期待的に成り立ちます。観察研究では共変量の群間差が存在しうるため、ANCOVAの調整が主目的となります。ただし共変量の群間差が極端に大きい場合は、外挿の問題が生じ調整の信頼性が低下します。
5. 処理による共変量への影響なし: 共変量は処理（群の割り当て）の前に測定されたものに限ります。処理後に測定された変数を共変量にすると、処理効果の一部を調整してしまいます。これはLord’s paradoxとして知られる問題と関連しています。
これらの仮定を検証する手順として、以下の診断法が推奨されます。
- 散布図で共変量と目的変数の関係を群別に確認（平行性・線形性）
- 残差プロットで等分散性を確認
- Q-Qプロットで正規性を確認
- 交互作用の検定で平行性を統計的に検定
ANCOVAの検出力の理論
ANOVAとの検出力比較
ANCOVAの検出力がANOVAより高くなるのは、共変量による説明で誤差分散が減少するためです。具体的には、共変量と目的変数の群内相関を $r_w$ とすると、ANCOVAの誤差分散はANOVAの誤差分散の $(1 – r_w^2)$ 倍になります。
$$ \sigma^2_{\text{ANCOVA}} = (1 – r_w^2) \sigma^2_{\text{ANOVA}} $$
$r_w = 0.5$ ならば誤差分散は75%に減少し、$r_w = 0.7$ ならば51%に減少します。相関が高いほどANCOVAの利点は大きくなります。
検出力の観点から、これがどの程度の影響を持つか具体的に計算しましょう。F統計量の非心度パラメータ（non-centrality parameter）は、ANOVAの場合、
$$ \lambda_{\text{ANOVA}} = \frac{n \sum_j \tau_j^2}{\sigma^2} $$
ANCOVAの場合、
$$ \lambda_{\text{ANCOVA}} = \frac{n \sum_j \tau_j^2}{(1 – r_w^2)\sigma^2} = \frac{\lambda_{\text{ANOVA}}}{1 – r_w^2} $$
ANCOVAの非心度パラメータは、ANOVAの $1/(1-r_w^2)$ 倍になります。たとえば $r_w = 0.6$ のとき、非心度パラメータは $1/0.64 \approx 1.56$ 倍になり、検出力が大幅に向上します。
ただし、自由度の損失（分母の自由度が $N-k$ から $N-k-1$ に減少）を考慮する必要があります。この損失は小さなサンプルでは無視できませんが、通常は分散の減少効果が上回ります。厳密には、ANCOVAがANOVAより有利になる条件は $r_w^2 > 1/(N-k-1)$ であり、サンプルサイズが十分大きければ（$N > k+2$）、わずかな相関でもANCOVAが有利です。
サンプルサイズの節約
ANCOVAの検出力向上は、実質的にサンプルサイズの節約と等価です。ANOVAで必要なサンプルサイズの $(1 – r_w^2)$ 倍で同じ検出力が得られます。
| $r_w$ | 分散減少率 $1 – r_w^2$ | 必要サンプル比 | 節約率 |
|---|---|---|---|
| 0.3 | 0.91 | 91% | 9% |
| 0.5 | 0.75 | 75% | 25% |
| 0.6 | 0.64 | 64% | 36% |
| 0.7 | 0.51 | 51% | 49% |
| 0.8 | 0.36 | 36% | 64% |
$r_w = 0.6$ の場合、ANCOVAに必要なサンプルサイズはANOVAの64%で済みます。これは臨床試験の設計において非常に重要な実用的利点です。たとえば、ANOVAで1群100人必要な試験でも、ベースライン測定値との相関が $r_w = 0.7$ であれば、ANCOVAを使うことで1群51人で同等の検出力が得られます。これは被験者数の削減、ひいてはコストと時間の大幅な削減につながります。
ランダム化実験でのANCOVAの利点
ランダム化実験では共変量の群間差は期待的にゼロですが、ANCOVAは依然として有用です。その理由は以下の2点です。
1. 誤差分散の減少: ランダム化により群間差がなくても、共変量が目的変数の変動を説明することで残差が小さくなり、検出力が向上します。臨床試験でベースライン値を共変量として含めることで、個人差による変動を除去し、処理効果をより精密に推定できます。
2. 偶然の偏りの補正: ランダム化は期待的に均等な群を作りますが、有限サンプルでは偶然の偏りが生じます。ANCOVAはこの偶然の偏りも補正します。たとえばランダム化で一方の群にやや高齢の被験者が偏った場合、年齢を共変量とするANCOVAがその影響を調整します。
このため、FDAのガイダンスでも臨床試験におけるベースライン調整としてのANCOVAの使用が推奨されています。EMAも同様のガイダンスを出しており、ランダム化実験でもANCOVAを主要な分析手法として使用することが国際的な標準となっています。
ランダム化が不完全な場合の注意
共変量の群間差が大きい場合（ランダム化が不十分な場合）、ANCOVAは共変量の調整に多くの「統計的パワー」を使うため、かえって検出力が低下する可能性があります。調整平均の標準誤差の式で見たように、$(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot})^2 / SS_{X, \text{within}}$ の項が大きくなり、推定の精度が下がります。
さらに、観察研究でANCOVAを適用する場合、共変量による調整は観測された共変量のみに基づくため、未観測の交絡因子は調整されません。因果推論の文脈では、これは重要な限界です。ANCOVAは「観測された共変量を同じ値に揃えたときの比較」であり、「もしランダム化していたら得られたであろう結果」とは必ずしも一致しません。
複数共変量への拡張
モデル
複数の共変量 $X^{(1)}, X^{(2)}, \dots, X^{(p)}$ を持つANCOVAモデルは、
$$ Y_{ij} = \mu + \tau_j + \sum_{m=1}^{p} \beta_m (X_{ij}^{(m)} – \bar{X}_{\cdot\cdot}^{(m)}) + \varepsilon_{ij} $$
行列表現では、デザイン行列に $p$ 個の共変量列が追加されます。
$$ \bm{X}_{\text{design}} = \begin{pmatrix} \bm{1}_N & \bm{D} & \bm{X}_c^{(1)} & \cdots & \bm{X}_c^{(p)} \end{pmatrix} $$
ここで $\bm{X}_c^{(m)} = \bm{X}^{(m)} – \bar{X}_{\cdot\cdot}^{(m)} \bm{1}_N$ は中心化された共変量ベクトルです。推定と検定の枠組みは1つの共変量の場合と同じですが、回帰係数ベクトル $\bm{\beta} = (\beta_1, \dots, \beta_p)^{\top}$ の推定が多変量回帰になります。
F統計量の一般化
複数共変量の場合のF統計量は、
$$ F = \frac{(SS_{\text{res, reduced}} – SS_{\text{res, full}}) / (k – 1)}{SS_{\text{res, full}} / (N – k – p)} $$
分母の自由度は $N – k – p$ です（$k$ 個の群パラメータと $p$ 個の共変量パラメータ）。
自由度の調整と共変量選択のトレードオフ
共変量が $p$ 個の場合、F統計量の分母の自由度は $N – k – p$ になります。共変量の数が増えると自由度が失われるため、不必要に多くの共変量を投入すると検出力がかえって低下する可能性があります。
一般に、共変量を追加することで得られる分散の減少（検出力の向上）と、自由度の損失（検出力の低下）のトレードオフがあります。共変量の追加が有利になる条件は、その共変量と目的変数の偏相関の2乗が $1/(N – k – p)$ を超える場合です。実用的には、偏相関が $0.1$ 未満の共変量は投入しない方がよいことが多いです。
共変量の選択基準
ANCOVAに投入する共変量の選択は、以下の基準に基づきます。
-
目的変数との相関: 共変量と目的変数の間に強い相関がある場合にのみ投入します。$r < 0.2$ 程度の共変量は自由度の損失のデメリットが利点を上回ることが多いです
-
事前の計画: 共変量はデータを見る前に（プロトコルで）決定すべきです。データを見てから共変量を選ぶとp-hackingになり、第一種過誤率が名目水準を超えます
-
処理前の測定: 共変量は処理（群の割り当て）の前に測定されたものに限ります。処理後に測定された変数を共変量にすると、処理効果の一部を調整してしまいます
-
共変量間の多重共線性の回避: 高い相関を持つ複数の共変量を投入すると、推定が不安定になります。VIF（分散膨張係数）が10を超える共変量は除外するか、主成分を使った次元削減を検討します
ここまでの理論を踏まえて、Pythonで実装し、具体的なデータで結果を確認しましょう。
Pythonによる実装
ANCOVAのスクラッチ実装
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
def ancova(y, x, groups, group_labels=None):
"""
一元配置ANCOVA（共変量1つ）
Parameters
----------
y : array 目的変数
x : array 共変量
groups : array 群のラベル（整数）
Returns
-------
result : dict
"""
y = np.asarray(y, dtype=float)
x = np.asarray(x, dtype=float)
groups = np.asarray(groups)
unique_groups = np.unique(groups)
k = len(unique_groups)
N = len(y)
x_grand = np.mean(x)
# --- フルモデル（群効果 + 共変量）---
# セル平均モデル: ダミー変数 + 共変量
X_full = np.zeros((N, k))
for j, g in enumerate(unique_groups):
X_full[:, j] = (groups == g).astype(float)
X_full = np.column_stack([X_full, x - x_grand])
beta_full = np.linalg.lstsq(X_full, y, rcond=None)[0]
y_pred_full = X_full @ beta_full
SS_res_full = np.sum((y - y_pred_full)**2)
beta_common = beta_full[-1] # 共通回帰係数
# --- 制約モデル（共変量のみ）---
X_red = np.column_stack([np.ones(N), x - x_grand])
beta_red = np.linalg.lstsq(X_red, y, rcond=None)[0]
y_pred_red = X_red @ beta_red
SS_res_red = np.sum((y - y_pred_red)**2)
# --- F統計量 ---
df1 = k - 1
df2 = N - k - 1
MS_res = SS_res_full / df2
F_val = ((SS_res_red - SS_res_full) / df1) / MS_res
p_val = 1 - stats.f.cdf(F_val, df1, df2)
# --- 調整平均と標準誤差 ---
# 群内平方和 SS_X_within
SS_X_within = sum(
np.sum((x[groups == g] - np.mean(x[groups == g]))**2)
for g in unique_groups
)
adj_means = {}
raw_means = {}
adj_se = {}
for j, g in enumerate(unique_groups):
mask = groups == g
n_j = np.sum(mask)
raw_means[g] = np.mean(y[mask])
x_group_mean = np.mean(x[mask])
adj_means[g] = raw_means[g] - beta_common * (x_group_mean - x_grand)
adj_se[g] = np.sqrt(
MS_res * (1/n_j + (x_group_mean - x_grand)**2 / SS_X_within)
)
# --- 通常のANOVA（比較用）---
group_data = [y[groups == g] for g in unique_groups]
F_anova, p_anova = stats.f_oneway(*group_data)
return {
"F": F_val, "p": p_val, "df": (df1, df2),
"beta": beta_common,
"adj_means": adj_means, "raw_means": raw_means,
"adj_se": adj_se,
"F_anova": F_anova, "p_anova": p_anova,
"SS_res_full": SS_res_full, "SS_res_red": SS_res_red,
"MS_res": MS_res
}
# --- データ生成 ---
np.random.seed(42)
n_per_group = 25
# 群A: 事前スコアが高い → 事後スコアも高い傾向
x_A = np.random.normal(70, 10, n_per_group)
y_A = 0.6 * x_A + 5 + np.random.normal(0, 5, n_per_group)
# 群B: 事前スコアが低い → しかし教授法Bの効果あり
x_B = np.random.normal(60, 10, n_per_group)
y_B = 0.6 * x_B + 10 + np.random.normal(0, 5, n_per_group)
y = np.concatenate([y_A, y_B])
x = np.concatenate([x_A, x_B])
groups = np.array([0] * n_per_group + [1] * n_per_group)
result = ancova(y, x, groups)
print("=" * 60)
print("ANCOVA vs ANOVA の比較")
print("=" * 60)
print(f"\n共変量の回帰係数 beta = {result['beta']:.4f}")
print(f"\n{'':20} {'群A':>10} {'群B':>10}")
print(f"{'観測平均':20} {result['raw_means'][0]:>10.2f} {result['raw_means'][1]:>10.2f}")
print(f"{'調整平均':20} {result['adj_means'][0]:>10.2f} {result['adj_means'][1]:>10.2f}")
print(f"{'調整平均のSE':20} {result['adj_se'][0]:>10.4f} {result['adj_se'][1]:>10.4f}")
print(f"\nANOVA: F({result['df'][0]}, {n_per_group*2-2}) = "
f"{result['F_anova']:.3f}, p = {result['p_anova']:.6f}")
print(f"ANCOVA: F({result['df'][0]}, {result['df'][1]}) = "
f"{result['F']:.3f}, p = {result['p']:.6f}")
print(f"\n残差平方和: フルモデル = {result['SS_res_full']:.2f}, "
f"制約モデル = {result['SS_res_red']:.2f}")
print(f"残差平均平方 (MS_res) = {result['MS_res']:.4f}")
この結果から、ANCOVAとANOVAの決定的な違いが確認できます。
-
観測平均と調整平均が異なる: 群Aは事前スコア（共変量）の平均が高いため、観測された事後スコアも高くなっています。しかし、共変量の影響を調整すると、群Bの方が実際には事後スコアが高い可能性があります。これは教授法Bの真の効果を反映しています
-
ANCOVAはANOVAと異なる結論を導く可能性がある: 共変量の調整により、ANOVAでは検出できなかった差がANCOVAで有意になったり、逆にANOVAで見えた差がANCOVAでは消失したりします
-
共変量の回帰係数: $\hat{\beta} \approx 0.6$ は事前テストの1点増加が事後テストの0.6点増加に対応することを意味し、データ生成時に設定した真の値と一致しています
-
残差平均平方が大幅に小さくなる: ANCOVAのMS_resはANOVAの残差平均平方と比較して小さくなっており、共変量が目的変数の変動を有効に説明していることを示しています
可視化
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
fig, axes = plt.subplots(1, 2, figsize=(13, 5.5))
# 散布図 + 回帰直線
ax = axes[0]
mask_A = groups == 0
mask_B = groups == 1
ax.scatter(x[mask_A], y[mask_A], color="blue", alpha=0.6, s=40,
label="Group A")
ax.scatter(x[mask_B], y[mask_B], color="red", alpha=0.6, s=40,
label="Group B")
# 共通傾斜の回帰直線
x_range = np.linspace(35, 95, 100)
x_grand = np.mean(x)
beta = result["beta"]
y_line_A = result["adj_means"][0] + beta * (x_range - x_grand)
y_line_B = result["adj_means"][1] + beta * (x_range - x_grand)
ax.plot(x_range, y_line_A, "b-", linewidth=2)
ax.plot(x_range, y_line_B, "r-", linewidth=2)
# 調整平均の表示
ax.axvline(x_grand, color="gray", linestyle=":", linewidth=1,
label=f"Grand mean of X = {x_grand:.1f}")
ax.plot(x_grand, result["adj_means"][0], "bs", markersize=12, zorder=5)
ax.plot(x_grand, result["adj_means"][1], "rs", markersize=12, zorder=5)
ax.set_xlabel("Pre-test score (covariate)", fontsize=12)
ax.set_ylabel("Post-test score (response)", fontsize=12)
ax.set_title("ANCOVA: Parallel regression lines", fontsize=13)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)
# 平均の比較
ax = axes[1]
x_pos = np.array([0, 1])
width = 0.3
ax.bar(x_pos - width/2, [result["raw_means"][0], result["raw_means"][1]],
width, color=["blue", "red"], alpha=0.4, edgecolor="white",
label="Unadjusted mean")
ax.bar(x_pos + width/2, [result["adj_means"][0], result["adj_means"][1]],
width, color=["blue", "red"], alpha=0.8, edgecolor="white",
label="Adjusted mean")
ax.set_xticks(x_pos)
ax.set_xticklabels(["Group A", "Group B"])
ax.set_ylabel("Post-test score", fontsize=12)
ax.set_title("Unadjusted vs Adjusted means", fontsize=13)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3, axis="y")
plt.tight_layout()
plt.show()
この可視化から、以下のことが確認できます。
-
左図の平行な回帰直線がANCOVAの本質を示している。2本の直線は同じ傾き（共通回帰係数 $\beta$）を持ち、切片のみが異なります。共変量の全体平均の位置（灰色の点線）での各直線上の値が調整平均です。四角いマーカーが各群の調整平均を表しています
-
右図は調整の効果を明確に示している。共変量の調整前（薄い色）と調整後（濃い色）で、群間差の方向や大きさが変わる可能性があります。これが交絡変数の影響を除去した真の処理効果です
平行性の検定と検出力のシミュレーション
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
np.random.seed(42)
# --- 平行性の検定 ---
def test_parallelism(y, x, groups):
"""群 x 共変量の交互作用検定"""
unique_groups = np.unique(groups)
k = len(unique_groups)
N = len(y)
x_grand = np.mean(x)
# 平行モデル（共通傾き）
X_par = np.zeros((N, k))
for j, g in enumerate(unique_groups):
X_par[:, j] = (groups == g).astype(float)
X_par = np.column_stack([X_par, x - x_grand])
b_par = np.linalg.lstsq(X_par, y, rcond=None)[0]
SS_par = np.sum((y - X_par @ b_par)**2)
# 交互作用モデル（群別傾き）
X_int = np.zeros((N, k))
for j, g in enumerate(unique_groups):
X_int[:, j] = (groups == g).astype(float)
for j, g in enumerate(unique_groups):
X_int = np.column_stack([
X_int,
(groups == g).astype(float) * (x - x_grand)
])
b_int = np.linalg.lstsq(X_int, y, rcond=None)[0]
SS_int = np.sum((y - X_int @ b_int)**2)
df1 = k - 1
df2 = N - 2 * k
F_int = ((SS_par - SS_int) / df1) / (SS_int / df2)
p_int = 1 - stats.f.cdf(F_int, df1, df2)
return F_int, p_int
F_par, p_par = test_parallelism(y, x, groups)
print(f"平行性の検定: F = {F_par:.3f}, p = {p_par:.4f}")
if p_par > 0.05:
print("-> 平行性の仮定は棄却されない（ANCOVAを適用可能）")
else:
print("-> 平行性の仮定が棄却された（ANCOVA注意）")
# --- 検出力のシミュレーション ---
n_sim = 2000
correlations = [0.0, 0.3, 0.5, 0.7, 0.9]
n_per = 20
true_effect = 3.0
power_anova = []
power_ancova = []
for r in correlations:
sig_anova = 0
sig_ancova = 0
for _ in range(n_sim):
x_all = np.random.normal(0, 10, 2 * n_per)
noise = np.random.normal(0, 10 * np.sqrt(1 - r**2), 2 * n_per)
y_all = r * x_all + noise
y_all[n_per:] += true_effect
g_all = np.array([0]*n_per + [1]*n_per)
# ANOVA
_, p_a = stats.f_oneway(y_all[:n_per], y_all[n_per:])
if p_a < 0.05:
sig_anova += 1
# ANCOVA
res = ancova(y_all, x_all, g_all)
if res['p'] < 0.05:
sig_ancova += 1
power_anova.append(sig_anova / n_sim)
power_ancova.append(sig_ancova / n_sim)
fig, ax = plt.subplots(figsize=(8, 5))
ax.plot(correlations, power_anova, 'bo-', markersize=8, linewidth=2,
label='ANOVA')
ax.plot(correlations, power_ancova, 'rs-', markersize=8, linewidth=2,
label='ANCOVA')
ax.set_xlabel('Correlation between covariate and outcome (r)', fontsize=12)
ax.set_ylabel('Power', fontsize=12)
ax.set_title('Power Comparison: ANOVA vs ANCOVA', fontsize=13)
ax.axhline(y=0.05, color='gray', linestyle=':', linewidth=1,
label='alpha = 0.05')
ax.legend(fontsize=11)
ax.grid(True, alpha=0.3)
ax.set_ylim(0, 1)
plt.tight_layout()
plt.show()
print("\n検出力の比較:")
for r, pa, pc in zip(correlations, power_anova, power_ancova):
ratio = (1 - r**2)
print(f" r = {r:.1f}: ANOVA = {pa:.3f}, ANCOVA = {pc:.3f}, "
f"理論的分散比 = {ratio:.3f}")
検出力シミュレーションのグラフから、共変量と目的変数の相関が高いほどANCOVAの優位性が顕著になることが確認できます。相関が0（$r = 0$）のときはANOVAとANCOVAの検出力はほぼ同じですが（わずかにANCOVAが自由度1の損失分だけ劣る可能性がある）、$r = 0.7$ 以上ではANCOVAの検出力がANOVAを大きく上回ります。これは理論的な予測 — 誤差分散が $(1-r^2)$ 倍に減少する — と一致しています。
また、各相関値での理論的分散比と検出力の実際の向上が概ね一致していることも確認できます。理論が予測する効率向上がシミュレーションで裏付けられているのです。
ANCOVAの実践的なガイドライン
いつANCOVAを使うべきか
以下のすべてが当てはまる場合に、ANCOVAが推奨されます。
- 目的変数と有意に相関する連続的な共変量が存在する
- 共変量は処理の前（ベースライン）に測定されている
- 回帰の平行性の仮定が（おおよそ）成り立つ
- 共変量は事前に決定されている（データを見てから選ばない）
ANCOVAとANOVAの選択フローチャート
実務での判断手順は以下の通りです。
- 連続的な共変量がありますか？ → No → ANOVA を使う
- 共変量と目的変数に有意な相関がありますか？ → No → ANOVA で十分（自由度の損失のみで利点なし）
- 平行性の仮定は成り立ちますか？ → No → Johnson-Neyman法 or 別手法
- すべてYes → ANCOVA を使う
よくある誤りとその対策
ANCOVAを適用する際に陥りやすい誤りを整理します。
1. 処理後の変数を共変量にする: 処理の結果として変わりうる変数を共変量に含めると、処理効果の一部を調整してしまいます。たとえば、薬の効果を調べる試験で服薬後の血圧を共変量にすることは不適切です。これはメディエーターの調整（mediation bias）と呼ばれ、処理効果の一部がブロックされてしまいます。
2. 共変量の多重選択（p-hacking）: データを見てから「有意な共変量」を選ぶと、多重比較問題が生じ、第一種の過誤率が膨らみます。共変量は分析プロトコルで事前に決定し、データ解析後に変更しないことが鉄則です。
3. 平行性の仮定の無視: 平行性の検定を行わずにANCOVAを適用すると、誤った調整平均を報告する可能性があります。特に交互作用が大きい場合、単一の調整平均は無意味であり、誤解を招きます。
4. 共変量の群間差が極端に大きい場合の過信: 共変量の群間差が非常に大きい場合（たとえば共変量の分布がほとんど重ならない場合）、調整は実質的に外挿となり、モデルの仮定への依存度が高くなります。このような状況では、ANCOVAの結果を慎重に解釈する必要があります。
5. 回帰希釈バイアスの無視: 共変量に測定誤差がある場合、調整が不完全になり、交絡の影響が残ります。特に観察研究では、この問題が深刻になりうるため、共変量の信頼性を報告し、感度分析を行うことが推奨されます。
ANCOVAと関連手法の比較
ANCOVAはいくつかの関連手法と比較されることがあります。
変化量（差分スコア）のANOVA: 目的変数として $Y_{\text{post}} – Y_{\text{pre}}$ を使うANOVAと、$Y_{\text{pre}}$ を共変量としたANCOVAは一般に異なる結果を与えます。変化量のANOVAは「事前テストから事後テストへの変化量」の群間差を見ますが、回帰効果（regression to the mean）の影響を受けやすいです。ANCOVAの方が統計的効率が良く、回帰効果への頑健性も高いことが多くの研究で示されています。
傾向スコアによる調整: 観察研究で多数の共変量を扱う場合、ANCOVAに直接すべての共変量を投入するのではなく、傾向スコア（propensity score）を1つの共変量として使う方法があります。傾向スコアは群への割り当て確率の推定値であり、多数の共変量を1次元に集約します。
重回帰分析との関係: ANCOVAは、ダミー変数と共変量を説明変数とする重回帰分析と数学的に等価です。したがって、ANCOVAは重回帰分析の特殊な場合として位置づけられます。ただし、ANCOVAでは「群間差の検定」が主目的であるのに対し、重回帰では「予測」や「変数の効果の推定」が主目的であることが多いです。
まとめ
本記事では、共分散分析（ANCOVA）の理論を直感的な説明から数学的な定式化、実装まで解説しました。
- ANCOVAはANOVAと回帰分析を組み合わせた手法であり、連続的な共変量の影響を調整した上で群間差を検定します。交絡変数の影響を統計的に除去できるため、処理効果の正確な評価が可能です
- パラメータ推定は最小二乗法で行われ、共変量の回帰係数 $\hat{\beta}$ は群内変動のみから推定されます。これはフリッシュ・ウォーの定理に基づいています
- 調整平均は共変量の群間差を統計的に除去した後の推定値であり、$\bar{Y}_{j, \text{adj}} = \bar{Y}_{\cdot j} – \hat{\beta}(\bar{X}_{\cdot j} – \bar{X}_{\cdot\cdot})$ で計算されます。その標準誤差は共変量の群間差が大きいほど増加します
- 回帰の平行性（共通傾斜）はANCOVAの最も重要な仮定であり、群ごとに回帰係数が異なる場合はJohnson-Neyman法などの代替手法を検討します
- ANCOVAはANOVAよりも高い検出力を持つことが多く、誤差分散が $(1 – r_w^2)$ 倍に減少します。$r_w = 0.7$ なら必要サンプルサイズが約半分で済みます
- 平方和の分解にはType I/II/III があり、ANCOVAではType III が標準的です。不均衡データでは特にType IIIの使用が重要です
- 複数共変量への拡張は自然に行えますが、共変量の追加による検出力向上と自由度損失のトレードオフを考慮する必要があります
次のステップとして、以下の記事も参考にしてください。