コピュラ (Copula) と多変量依存構造 — 周辺分布と依存関係を分離する理論

2008年の世界金融危機の引き金となった住宅ローン担保証券(CDO)の格付モデルでは、ガウシアンコピュラ と呼ばれる多変量依存モデルが使われていました。事後検証では「ガウシアンコピュラはテール依存性をゼロと仮定するため、複数地域の不動産価格が同時に暴落する確率を極端に過小評価していた」ことが致命的な欠陥として指摘されています。ここで疑問が湧きます——「コピュラ」とは一体何で、なぜ相関係数だけでは不十分なのでしょうか?

直感的には、相関係数は2変数の関係を1つの数字に押し込めようとしますが、現実の依存関係はそんなに単純ではありません。たとえば「平常時はバラバラに動くが、危機のときだけ揃って下落する」ような変数ペア(株式リターン同士はまさにそうです)と、「常に同程度に連動する」変数ペアは、相関係数が同じでも全く異なる挙動を示します。コピュラは、こうした 依存の「質」 までモデル化する道具です。

この概念が活きる場面は多岐にわたります。第一に 金融ポートフォリオのリスク管理(複数資産が同時に暴落する確率の評価、CVaR計算)。第二に 保険数理(複数の保険請求が同時に発生する条件付き確率、複数災害の同時発生リスク)。第三に 信頼性工学(システムの複数部品が相関を持って故障するときの故障時間分布)、第四に 気象・気候(複数地点での極値現象の同時発生確率)。いずれも「個別の周辺分布はそれぞれ別のモデルで合っているが、それらを組み合わせた多変量分布をどう作るか」という問題に直面します。コピュラは、周辺分布と依存構造をきれいに分離することで、この問題を体系的に解決します。

本記事の内容

  • なぜ相関係数だけでは依存構造を表せないかの直感
  • Sklarの定理 — 多変量分布Fをコピュラと周辺分布に一意分解する基礎定理の導出
  • ガウシアン・t・Archimedean(Clayton/Frank/Gumbel)コピュラの数学的定義と挙動の違い
  • Kendallのτ・Spearmanのρ とコピュラパラメータの関係、テール依存係数 $\lambda_U, \lambda_L$ の意味
  • Pythonによる可視化(同じ相関でもコピュラを変えると尾の挙動が劇的に変わる様子)
  • S&P500風シミュレーションでガウシアンコピュラがVaRを過小評価するメカニズムの実証
  • copulas ライブラリを使った実データへのフィッティング例

Sklarの定理: 多変量分布=コピュラ+周辺分布の分解イメージ

上の図がコピュラ理論のエッセンスを表しています。左の「周辺分布 $F_1$」と右の「周辺分布 $F_2$」はそれぞれ別の形をしていますが、中央の「コピュラ」は両者の依存の構造だけを $[0,1]^2$ 上に切り出したものです。Sklarの定理は、この3つのパーツ——周辺分布たちとコピュラ——が合わさって元の多変量分布を一意に決定することを保証します。周辺分布を入れ替えてもコピュラを残せば「依存の質」は変わらず、コピュラだけ変えれば周辺分布はそのままに「依存の構造」を自由に変えられる、というのが直感です。

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

相関係数だけでは何が足りないのか

多変量データを扱うとき、最初に思いつく依存性指標は Pearsonの相関係数 $\rho_{XY} = \mathrm{Cov}(X,Y)/(\sigma_X \sigma_Y)$ でしょう。確かに直線的な関係を1つの数字で要約できる便利な指標ですが、ここには3つの本質的な限界があります。

第一に、相関係数は線形関係しか測らないことです。$Y = X^2$($X \sim \mathcal{N}(0,1)$)のように完全に決定論的な関係でも、$\rho_{XY} = 0$ になってしまいます。これは「依存はゼロ」ではなく、依存の形が直線でないだけです。

第二に、周辺分布の形に大きく影響を受けることです。同じ「依存の構造」を持っていても、周辺分布が正規分布の場合と対数正規分布の場合とでは、Pearsonの相関係数の値は変わってしまいます。これは、「依存構造」と「周辺の形」を一括りに評価していることに起因します。

第三に、そして最も重要な点として、極端な事象が同時に起こる確率(テール依存性)を捉えられないことです。たとえば次の2つの2変量分布を考えてみましょう。

  • 分布A: 2変量正規分布、相関 $\rho = 0.5$
  • 分布B: 2変量t分布(自由度3)、相関 $\rho = 0.5$

これらは「2変数の線形相関」は同じです。しかし「両方が極端に大きい値を取る確率」、たとえば $P(X > 3, Y > 3)$ は、分布Bが分布Aの数倍にもなり得ます。金融でいう「両資産が同時に暴落する確率」がまさにこのテール領域の問題で、相関係数では捉えきれません。

これら3つの問題に共通する根本原因は、「依存構造」と「周辺分布の形」が混在していることです。$X$ が正規分布で $Y$ が指数分布で $Z$ がガンマ分布で……といった混合した状況で、それでも依存だけを取り出して比較したい。この望みを叶える理論的枠組みが コピュラ です。

Sklarの定理 — 分布を周辺と依存に分解する

確率積分変換の復習

コピュラの定義に入る前に、1次元の重要な事実を確認しておきます。連続な確率変数 $X$ の累積分布関数(CDF)を $F_X$ とすると、$U = F_X(X)$ は 区間 $[0,1]$ 上の一様分布 に従います。

この事実は 確率積分変換 (probability integral transform) と呼ばれ、$X$ の値そのものではなく「$X$ が分布の中でどの位置にあるか(パーセンタイル)」を見る変換と解釈できます。逆向きに、$U \sim \mathrm{Unif}(0,1)$ なら $X = F_X^{-1}(U)$ は元の分布に従います。

このシンプルな事実が、周辺分布の形を取り除いて依存構造だけを抽出する ための鍵になります。$d$ 個の変数 $(X_1, \dots, X_d)$ から $(U_1, \dots, U_d) = (F_1(X_1), \dots, F_d(X_d))$ を作れば、各 $U_i$ はすべて $\mathrm{Unif}(0,1)$ で、周辺分布の情報は完全に消えています。残るのは「$U_i$ たちがどう絡み合うか」という依存構造だけです。

確率積分変換: CDFで変換すると一様分布になる

左のヒストグラムは元の正規分布サンプル $X$ の分布で、右は $U = F_X(X)$ で変換した後の分布です。変換後は山が消えて $[0,1]$ 上にほぼ均等(Uniform分布)になっており、「どの位置に値があるか」という情報——つまり周辺分布の形——が完全に除去されていることが分かります。このゼロから一様分布への「均等化」を利用して、複数変数の一様乱数 $(U_1, U_2, \dots)$ の間の依存だけを取り出したものがコピュラです。

コピュラの定義

ここで コピュラ を定義します。$d$ 次元コピュラ $C : [0,1]^d \to [0,1]$ は、次の3条件を満たす関数です。

  1. 各 $u_i$ の周辺分布が $\mathrm{Unif}(0,1)$: つまり $u_j = 1$($j \ne i$)を代入すると $C(1,\dots,1,u_i,1,\dots,1) = u_i$
  2. ある $u_i = 0$ なら $C(\dots,0,\dots) = 0$
  3. $d$ 次元の単調増加性($d$-increasing)

要するに、コピュラとは 「$[0,1]^d$ 上で、全ての周辺分布が一様分布である多変量CDF」 です。周辺の情報が全部 $\mathrm{Unif}(0,1)$ に潰されているので、残るのは依存の構造そのものになります。

Sklarの定理

Sklarの定理は、コピュラの存在意義を確立する基礎定理です。

Sklarの定理 (1959): 任意の $d$ 次元連続CDF $F$ と、その周辺分布 $F_1, \dots, F_d$ に対して、ただ一つのコピュラ $C$ が存在して $$ > F(x_1, x_2, \dots, x_d) = C\bigl(F_1(x_1), F_2(x_2), \dots, F_d(x_d)\bigr) > $$ が成り立つ。逆に、任意のコピュラ $C$ と任意の1次元CDF $F_1, \dots, F_d$ に対して、上式で定義される $F$ は周辺分布が $F_1, \dots, F_d$ である正当な $d$ 次元CDFとなる。

この定理が告げているのは、「多変量分布は、周辺分布たちとコピュラに一意に分解できる」 ということです。逆に言えば、好きな周辺分布たちと好きなコピュラを自由に組み合わせて、新しい多変量分布をいくらでも作れるということでもあります。これは多変量モデリングにとって革命的な自由度を与えてくれます。

連続な場合の導出

連続な場合の導出は驚くほどシンプルです。$U_i = F_i(X_i) \sim \mathrm{Unif}(0,1)$ より、$X_i = F_i^{-1}(U_i)$ と書けます。すると、

$$ F(x_1, \dots, x_d) = P(X_1 \leq x_1, \dots, X_d \leq x_d) $$

ですが、各事象 $\{X_i \leq x_i\}$ は $\{F_i(X_i) \leq F_i(x_i)\}$($F_i$ は単調増加なので同値)、つまり $\{U_i \leq F_i(x_i)\}$ に書き換えられます。よって、

$$ F(x_1, \dots, x_d) = P\bigl(U_1 \leq F_1(x_1), \dots, U_d \leq F_d(x_d)\bigr) $$

となります。右辺は $(U_1, \dots, U_d)$ の同時CDFを点 $(F_1(x_1), \dots, F_d(x_d))$ で評価したものなので、これをコピュラ $C$ と呼んで、

$$ F(x_1, \dots, x_d) = C\bigl(F_1(x_1), \dots, F_d(x_d)\bigr) $$

が得られます。一意性は $F_i$ が連続(よって全単射)であることから従います。

コピュラ密度

コピュラCDFが定義できれば、その偏微分でコピュラ密度 $c$ が得られます。

$$ c(u_1, \dots, u_d) = \frac{\partial^d C(u_1, \dots, u_d)}{\partial u_1 \cdots \partial u_d} $$

そして、多変量密度関数 $f$ はコピュラ密度と周辺密度 $f_i$ の積に分解されます:

$$ f(x_1, \dots, x_d) = c\bigl(F_1(x_1), \dots, F_d(x_d)\bigr) \cdot \prod_{i=1}^{d} f_i(x_i) $$

この分解は実用上極めて重要です。最尤推定では、周辺の対数尤度とコピュラの対数尤度が分離するため、「まず周辺をフィットし、次にコピュラをフィットする」 という二段階推定(IFM法: Inference Functions for Margins)が可能になります。

コピュラの境界

任意のコピュラは、Fréchet–Hoeffding境界 と呼ばれる上下限を持ちます。

$$ W(\bm{u}) := \max\Bigl(\sum_i u_i – d + 1,\ 0\Bigr) \leq C(\bm{u}) \leq \min(u_1, \dots, u_d) =: M(\bm{u}) $$

上限 $M$ は 完全正相関($X_i$ たちが単調増加関数で連動する状態)、下限 $W$ は2次元では 完全負相関($U_2 = 1 – U_1$)に対応します。$d \geq 3$ では $W$ はコピュラでなくなる(最低でも下限を達成する2変数ペアしか作れない)ことが知られていますが、上限 $M$ は常にコピュラです。

Frechet-Hoeffding境界: 完全正依存/独立/完全負依存の3例

3つの極端なケースを図示しました。左の「完全正依存 M」は $U_2 = U_1$ の対角線上に全点が集中し、2変数が完全に連動します。中央の「独立コピュラ」は $[0,1]^2$ に均等に散らばり、全く無相関な状態です。右の「完全負依存 W」は $U_2 = 1 – U_1$ の逆対角線上に集中し、一方が大きいほど他方が小さくなります。現実のコピュラはこの3つの間のどこかに位置し、Fréchet–Hoeffding境界がその取りうる範囲を定めています。

Sklarの定理は「コピュラと周辺が分離できる」という存在の枠組みを与えてくれました。では具体的に、どんなコピュラを使えばよいのでしょうか?次のセクションから、代表的なコピュラ族を順に見ていきます。

ガウシアンコピュラとtコピュラ — 楕円型コピュラ

ガウシアンコピュラ

最も基本的なコピュラは、多変量正規分布から構成される ガウシアンコピュラ です。$\Phi$ を標準正規分布のCDF、$\Phi_{\bm{R}}$ を相関行列 $\bm{R}$ を持つ $d$ 次元標準正規分布のCDFとすると、ガウシアンコピュラは

$$ C_{\bm{R}}^{\mathrm{Gauss}}(u_1, \dots, u_d) = \Phi_{\bm{R}}\bigl(\Phi^{-1}(u_1), \dots, \Phi^{-1}(u_d)\bigr) $$

で定義されます。直感的な手順は明快です: (1) $(u_1, \dots, u_d) \in [0,1]^d$ を $\Phi^{-1}$ で標準正規空間に持ち上げ、(2) 多変量正規分布のCDFで評価し、(3) その値を返す。

ガウシアンコピュラの密度は、$\bm{z}_i = \Phi^{-1}(u_i)$ と置くと

$$ c_{\bm{R}}^{\mathrm{Gauss}}(\bm{u}) = \frac{1}{\sqrt{\det \bm{R}}}\exp\!\left(-\frac{1}{2}\bm{z}^{\top}(\bm{R}^{-1} – \bm{I})\bm{z}\right) $$

となります。これは多変量正規密度を標準正規の積で割った形で、「正規同士の依存だけ」を抽出していると見ることができます。

tコピュラ

tコピュラ は、多変量t分布(自由度 $\nu$、相関行列 $\bm{R}$)から同様に構成されます。$t_\nu$ を1変量t分布のCDF、$t_{\nu, \bm{R}}$ を多変量t分布のCDFとすると、

$$ C_{\nu, \bm{R}}^{t}(u_1, \dots, u_d) = t_{\nu, \bm{R}}\bigl(t_\nu^{-1}(u_1), \dots, t_\nu^{-1}(u_d)\bigr) $$

です。tコピュラは自由度 $\nu \to \infty$ でガウシアンコピュラに収束しますが、有限の $\nu$ ではテール依存性を持つという大きな違いがあります(後の節で詳述)。

楕円型コピュラの限界

ガウシアン・tコピュラはいずれも、$\bm{u}$ の周りで 楕円対称 な依存構造を持ちます。これは「正方向の極端と負方向の極端で、依存の強さが対称」という意味です。しかし現実のデータ、特に金融では「下方向の極端(大暴落)だけが連動する」という 非対称な依存 がしばしば観察されます。たとえば株式市場は、平穏な上昇局面では各銘柄がバラバラに動くのに、暴落局面では一斉に下がるという非対称性を見せます。

この非対称性を表現できないのが楕円型コピュラの本質的限界です。これを乗り越えるために登場するのが、次に解説する Archimedeanコピュラ族 です。

Archimedeanコピュラ族 — 非対称依存の表現

Archimedeanコピュラの定義

Archimedeanコピュラ は、生成関数 (generator) と呼ばれる1変数関数 $\varphi : [0,1] \to [0, \infty]$ を介して定義されます。$\varphi$ が連続で狭義単調減少、$\varphi(1) = 0$、かつ完全単調(CMonotone)であるとき、

$$ C(u_1, \dots, u_d) = \varphi^{-1}\!\bigl(\varphi(u_1) + \varphi(u_2) + \cdots + \varphi(u_d)\bigr) $$

がコピュラになります。生成関数 $\varphi$ を1つ選ぶだけでコピュラが決まるので、扱いやすく、また依存構造のパラメトリックモデルとして表現力が高いのが特徴です。

代表的な3つを順に見ていきます。

Claytonコピュラ

Claytonコピュラ はパラメータ $\theta > 0$ で、生成関数

$$ \varphi(t) = \frac{1}{\theta}(t^{-\theta} – 1) $$

から作られます。2変数の場合の明示形は

$$ C_\theta^{\mathrm{Clayton}}(u_1, u_2) = \bigl(u_1^{-\theta} + u_2^{-\theta} – 1\bigr)^{-1/\theta} $$

です。Claytonコピュラの特徴は 下方テール依存性 が強いこと——両方の変数が同時に小さい値を取る傾向が顕著です。$\theta \to 0$ で独立、$\theta \to \infty$ で完全依存に近づきます。リスク管理で「大損失の同時発生」をモデル化するのに適しています。

Gumbelコピュラ

Gumbelコピュラ はパラメータ $\theta \geq 1$ で、

$$ \varphi(t) = (-\ln t)^{\theta}, \quad C_\theta^{\mathrm{Gumbel}}(u_1, u_2) = \exp\!\left[-\bigl((-\ln u_1)^\theta + (-\ln u_2)^\theta\bigr)^{1/\theta}\right] $$

で与えられます。Gumbelの特徴は 上方テール依存性 が強いこと——両方が同時に大きい値を取る傾向が強いコピュラです。極値理論(複数地点の年最大降水量が同時に大きくなるか等)でよく使われます。

Frankコピュラ

Frankコピュラ はパラメータ $\theta \in \mathbb{R}\setminus\{0\}$ で、

$$ \varphi(t) = -\ln\!\frac{e^{-\theta t} – 1}{e^{-\theta} – 1} $$

から作られ、2変数では

$$ C_\theta^{\mathrm{Frank}}(u_1, u_2) = -\frac{1}{\theta}\ln\!\left[1 + \frac{(e^{-\theta u_1} – 1)(e^{-\theta u_2} – 1)}{e^{-\theta} – 1}\right] $$

となります。Frankコピュラはテール独立(後述)で、依存が中心に集中するのが特徴で、正負どちらの依存も表現できます($\theta > 0$ で正、$\theta < 0$ で負)。

3つのコピュラの直感的な使い分け

コピュラ 非対称性 テール挙動 典型的な応用
Gaussian 対称 テール独立($\lambda = 0$) 標準的な相関モデリング
t 対称 両テール依存 株式リターン、リスク管理の改良版
Clayton 下方依存 $\lambda_L > 0, \lambda_U = 0$ 大損失の同時発生、信用リスク
Gumbel 上方依存 $\lambda_U > 0, \lambda_L = 0$ 極値の同時発生、洪水・台風
Frank 対称 テール独立 中心的依存、対称な順序関係

各コピュラの密度ヒートマップ(ガウシアン/t/Clayton/Gumbel)

各コピュラの密度を $[0,1]^2$ 上でヒートマップ表示しました。赤いほど確率密度が高い領域です。ガウシアンは $(0.5, 0.5)$ の中央付近に密度が集中し、隅は薄い楕円対称の形です。tコピュラ(自由度3) はガウシアンと似た形をしていますが、4隅にも赤みが残っており、テール領域に確率が「漏れている」ことが分かります。Claytonコピュラ は $(0,0)$ の左下隅にのみ密度の山があり、下方向にだけ依存が集中しています。Gumbelコピュラ は逆に $(1,1)$ の右上隅に密度が集まっています。密度の形から、各コピュラが「どの極端な状況が同時に起こりやすいか」をモデル化しているかが一目で把握できます。

コピュラの種類とパラメータが決まれば依存の形が決まりますが、その「強さ」を共通の尺度で測りたい場面が必ず出てきます。次節では、コピュラ依存の強さを定量化する2つの指標——Kendallのτ・Spearmanのρ——と、極端事象同士の連動性を測るテール依存係数を導入します。

Kendallのτ・Spearmanのρ・テール依存係数

順位相関係数

Pearsonの相関係数は周辺分布の形に影響されるため、コピュラの「依存の強さ」を測るには不向きです。そこで活躍するのが 順位相関係数 です。これらはコピュラだけで決まり、周辺分布を変えても値が変わらないという美しい性質を持ちます。

Kendallのτ は、2つの独立な観測 $(X_1, Y_1), (X_2, Y_2)$ に対し、

$$ \tau = P\bigl((X_1 – X_2)(Y_1 – Y_2) > 0\bigr) – P\bigl((X_1 – X_2)(Y_1 – Y_2) < 0\bigr) $$

で定義され、「協和(concordant)ペアと不協和(discordant)ペアの確率差」と読めます。順序のみに依存するため周辺分布によりません。コピュラを使うと、

$$ \tau = 4\int_{[0,1]^2} C(u,v)\,dC(u,v) – 1 $$

と表されます。

Spearmanのρ は、$X, Y$ をそれぞれの周辺CDFで変換した順位(パーセンタイル)の間のPearson相関で、

$$ \rho_S = 12\int_{[0,1]^2} uv\,dC(u,v) – 3 = 12\int_{[0,1]^2} C(u,v)\,du\,dv – 3 $$

となります。これも周辺分布によりません。

主要コピュラとパラメータの関係

各コピュラについて、パラメータとKendallのτの関係は次のように与えられます。

  • Gaussian / t: $\tau = \frac{2}{\pi}\arcsin(\rho)$
  • Clayton: $\tau = \theta / (\theta + 2)$
  • Gumbel: $\tau = 1 – 1/\theta$
  • Frank: $\tau = 1 – \frac{4}{\theta}\bigl(1 – D_1(\theta)\bigr)$($D_1$ はDebye関数)

この対応関係のおかげで、データから推定したKendallのτを使って各コピュラのパラメータを簡単にキャリブレートできます。たとえばGumbelコピュラなら、$\hat{\theta} = 1/(1 – \hat{\tau})$ という単純な公式で済みます。

Kendallのτとコピュラパラメータの対応関係

3つのコピュラについて、Kendallのτとパラメータの関係を曲線で示しました。ガウシアン/tコピュラの $\rho = \sin(\pi\tau/2)$ は非線形な増加曲線で、$\tau=0$ で $\rho=0$、$\tau=1$ で $\rho=1$ となります。Claytonの $\theta = 2\tau/(1-\tau)$ は $\tau \to 1$ で急激に大きくなる(完全依存に近い強い下方依存)ことが分かります。Gumbelの $\theta = 1/(1-\tau)$ は $\tau=0$ で $\theta=1$(下限)から始まり、やはり $\tau \to 1$ で無限大になります。縦破線($\tau=0.5$)でのそれぞれの値が、前節の散布図で使ったパラメータに対応しています。

テール依存係数 — 危機時の連動性

順位相関係数は「全体の依存の強さ」を測りますが、リスク管理では「極端な状態でどれだけ連動するか」が決定的に重要です。これを定量化するのが テール依存係数 です。

上方テール依存係数 $\lambda_U$ は、

$$ \lambda_U = \lim_{q \to 1^-} P\bigl(F_2(X_2) > q \mid F_1(X_1) > q\bigr) = \lim_{q \to 1^-} \frac{1 – 2q + C(q, q)}{1 – q} $$

で定義されます。直感的には、「$X_1$ がパーセンタイル $q$ を超えて極端な値を取ったとき、$X_2$ も同じく $q$ を超える条件付き確率の極限」です。$\lambda_U > 0$ なら 上方テール依存、$\lambda_U = 0$ なら テール独立 と呼びます。

下方テール依存係数 $\lambda_L$ も同様に、

$$ \lambda_L = \lim_{q \to 0^+} P\bigl(F_2(X_2) \leq q \mid F_1(X_1) \leq q\bigr) = \lim_{q \to 0^+} \frac{C(q, q)}{q} $$

で定義されます。

コピュラ別のテール依存係数

コピュラ $\lambda_L$ $\lambda_U$
Gaussian ($\rho < 1$) 0 0
t (自由度 $\nu$, 相関 $\rho$) $2t_{\nu+1}\bigl(-\sqrt{(\nu+1)(1-\rho)/(1+\rho)}\bigr)$ 同左
Clayton ($\theta > 0$) $2^{-1/\theta}$ 0
Gumbel ($\theta > 1$) 0 $2 – 2^{1/\theta}$
Frank 0 0

ここが2008年金融危機の本質です。ガウシアンコピュラは、相関 $\rho < 1$ である限り $\lambda_U = \lambda_L = 0$、つまり「どんなに強い相関を入れても、同時に極端を取る条件付き確率は最終的にゼロに収束する」という性質を持ちます。住宅ローン担保証券の格付モデルがこのコピュラを使っていたため、「複数地域の住宅価格が同時に大暴落する確率」が極端に過小評価されていました。実際にはt分布のような重い裾を持つコピュラを使えば $\lambda > 0$ となり、リスク量(VaRやCVaR)は大幅に増加します。

このテール依存性こそが、コピュラ理論を実用上のリスク管理に欠かせないものにしている理由です。理論はここまでとして、いよいよPythonで「同じ相関でも、コピュラを変えると尾の挙動が劇的に変わる」様子を可視化していきましょう。

Python実装 — 可視化と比較

4種類のコピュラからのサンプリング

まず、ガウシアン・t・Clayton・Gumbel の4つのコピュラから2変数サンプルを生成し、散布図で比較します。各コピュラは「Kendallのτ = 0.5」となるようパラメータを揃え、依存の強さは同じだが「形」が異なる様子を見ます。

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm, t, multivariate_normal

np.random.seed(42)
N = 3000  # サンプル数

# --- Gaussianコピュラ ---
def sample_gaussian_copula(rho, n):
    """相関 rho のガウシアンコピュラから n サンプル"""
    cov = np.array([[1.0, rho], [rho, 1.0]])
    z = np.random.multivariate_normal([0, 0], cov, size=n)
    u = norm.cdf(z)  # 標準正規CDFで [0,1]^2 に変換
    return u

# --- tコピュラ ---
def sample_t_copula(rho, nu, n):
    """相関 rho, 自由度 nu の tコピュラから n サンプル"""
    cov = np.array([[1.0, rho], [rho, 1.0]])
    z = np.random.multivariate_normal([0, 0], cov, size=n)
    # カイ二乗を使って t に変換
    chi = np.random.chisquare(nu, size=n)
    t_sample = z / np.sqrt(chi / nu)[:, None]
    u = t.cdf(t_sample, df=nu)
    return u

# --- Claytonコピュラ ---
def sample_clayton_copula(theta, n):
    """パラメータ theta の Claytonコピュラから n サンプル(条件付き法)"""
    u1 = np.random.uniform(size=n)
    v = np.random.uniform(size=n)
    # 条件付きCDF C(u2|u1) = v を u2 について解く
    u2 = (u1**(-theta) * (v**(-theta/(1+theta)) - 1) + 1)**(-1/theta)
    return np.column_stack([u1, u2])

# --- Gumbelコピュラ ---
def sample_gumbel_copula(theta, n):
    """パラメータ theta の Gumbelコピュラから n サンプル
    安定分布の混合表現を使う"""
    # 正の安定分布 S(1/theta) を Chambers-Mallows-Stuck 法で生成
    alpha = 1.0 / theta
    U = np.random.uniform(-np.pi/2, np.pi/2, size=n)
    W = np.random.exponential(scale=1.0, size=n)
    S = (np.sin(alpha*U) / np.cos(U)**(1/alpha)
         * (np.cos(U - alpha*U) / W)**((1-alpha)/alpha))
    # E1, E2 ~ Exp(1)
    E = np.random.exponential(scale=1.0, size=(n, 2))
    u = np.exp(-(E / S[:, None])**(1/theta))
    return u

# τ=0.5 になるパラメータ設定
rho_05 = np.sin(np.pi/2 * 0.5)        # gauss/t: tau = 2/pi * arcsin(rho)
theta_clayton_05 = 2.0 * 0.5 / (1 - 0.5)  # Clayton: tau = theta/(theta+2)
theta_gumbel_05 = 1.0 / (1 - 0.5)         # Gumbel: tau = 1 - 1/theta

samples = {
    'Gaussian (rho=0.71)': sample_gaussian_copula(rho_05, N),
    't (rho=0.71, nu=3)':  sample_t_copula(rho_05, 3, N),
    'Clayton (theta=2)':   sample_clayton_copula(theta_clayton_05, N),
    'Gumbel (theta=2)':    sample_gumbel_copula(theta_gumbel_05, N),
}

ここでのポイントは、各コピュラのサンプリング方法が大きく違うことです。ガウシアン・tは多変量正規/多変量t分布からCDF変換するだけで済みますが、ClaytonとGumbelは「条件付きCDFを解く」「安定分布の混合表現を使う」など、独自のアルゴリズムが必要になります。これはArchimedeanコピュラが「明示的な多変量分布の構成」ではなく「生成関数を介した抽象的な定義」であることの裏返しです。

4種コピュラの散布図比較(同Kendallτ=0.5)

Kendallのτを0.5に揃えた4つのコピュラから生成した3,000点の散布図です。赤点は「$u_1 < 0.05$ かつ $u_2 < 0.05$(左下隅)」または「$u_1 > 0.95$ かつ $u_2 > 0.95$(右上隅)」の極端な点です。ガウシアンコピュラは隅にほとんど点がなく、tコピュラは両隅に赤点が目立ちます。Claytonは左下隅だけが密集し(下方依存)、Gumbelは右上隅だけが密集します(上方依存)。同じKendallのτでも、コピュラが違えばこれほど異なる「依存の形」を生じることが一目瞭然です。

散布図で比較

生成したサンプルを散布図で並べてみましょう。

fig, axes = plt.subplots(1, 4, figsize=(18, 4.5))
for ax, (name, u) in zip(axes, samples.items()):
    ax.scatter(u[:, 0], u[:, 1], s=4, alpha=0.4)
    ax.set_title(name, fontsize=11)
    ax.set_xlabel('U1')
    ax.set_ylabel('U2')
    ax.set_xlim(0, 1); ax.set_ylim(0, 1)
    ax.set_aspect('equal')
    ax.grid(True, alpha=0.3)
plt.suptitle('Same Kendall tau (=0.5), Different Copulas — Look at the corners', y=1.03, fontsize=13)
plt.tight_layout()
plt.savefig('copula_scatter.png', dpi=150, bbox_inches='tight')
plt.show()

この4つの散布図から、コピュラの本質的な違いがはっきり読み取れます。まず Gaussianコピュラ は楕円対称な依存で、4つの隅にはほとんど点がありません(テール独立)。tコピュラ(自由度3) は中心の形はGaussianに似ていますが、$(0,0)$ と $(1,1)$ の両隅に点が密集しており、両側のテール依存が明確に見えます。Claytonコピュラ は $(0,0)$ の隅だけに鋭い集中が見られ、下方テール依存だけが強い非対称な依存です。逆に Gumbelコピュラ は $(1,1)$ の隅に集中し、上方テール依存だけが顕著です。Kendallのτはどれも 0.5 で「同じ依存の強さ」なのに、コピュラを変えるだけで尾の挙動がこれほど違うことが視覚的に分かります。

周辺ガウス × 異なるコピュラの可視化

実際の応用では、コピュラを変えて周辺ガウスと組み合わせた多変量分布がどう変わるかが重要です。同じ周辺分布で、コピュラだけ変えてみましょう。

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm

# サンプルした u を標準正規に変換(周辺をガウスに)
def to_normal(u):
    return norm.ppf(u)

fig, axes = plt.subplots(1, 4, figsize=(18, 4.5))
for ax, (name, u) in zip(axes, samples.items()):
    x = to_normal(u)
    ax.scatter(x[:, 0], x[:, 1], s=4, alpha=0.4)
    # サンプル相関係数を計算(Pearson)
    rho_p = np.corrcoef(x[:, 0], x[:, 1])[0, 1]
    ax.set_title(f'{name}\nPearson rho = {rho_p:.3f}', fontsize=10)
    ax.set_xlabel('X1 ~ N(0,1)')
    ax.set_ylabel('X2 ~ N(0,1)')
    ax.set_xlim(-4, 4); ax.set_ylim(-4, 4)
    ax.set_aspect('equal')
    ax.grid(True, alpha=0.3)
    # 危機ゾーン: X1<-2 かつ X2<-2 を赤で強調
    mask = (x[:, 0] < -2) & (x[:, 1] < -2)
    ax.scatter(x[mask, 0], x[mask, 1], s=15, color='red', alpha=0.8)
plt.suptitle('Same marginal N(0,1), Different copulas — Red = joint crash zone', y=1.03, fontsize=13)
plt.tight_layout()
plt.savefig('copula_normal_marginal.png', dpi=150, bbox_inches='tight')
plt.show()

このグラフから、リスク管理上極めて重要な事実が読み取れます。第一に、周辺分布が同じ標準正規でも、Pearsonの相関係数はコピュラごとに微妙に違います(Gaussianはちょうど0.71付近、Claytonは値が少しずれる)。これは「Pearsonの相関は周辺だけでなくコピュラの形にも影響される」ことを示しています。第二に、赤で強調した「同時に-2以下に下落する」ゾーンの点の数が、コピュラによって大きく異なります。Gaussianはわずかしかないのに対し、t やClaytonでは明らかに多い。同じ相関を持つように見えても、「危機時の共倒れ確率」は数倍違う のです。

正規周辺+コピュラ別の同時暴落ゾーン比較

周辺分布を標準正規に固定してコピュラだけ変えたとき、「同時暴落ゾーン(両変数が-2σ以下)」の点数がタイトル部分に表示されています。ガウシアンコピュラでは同時暴落点が最も少なく、tコピュラやClaytonコピュラでは明らかに多い点が確認できます。Gumbelコピュラは下方依存が弱いため、左下隅の同時暴落はむしろガウシアンより少ない場合もあります。「相関を揃えてもリスク評価が大きく変わる」というコピュラ理論の核心が、散布図の赤点の数の差として直感的に現れています。

テール依存係数の数値計算

理論で求めた $\lambda_L, \lambda_U$ の値を、サンプルから推定してみましょう。

import numpy as np

def empirical_tail_dependence(u, q_grid):
    """サンプル u (n, 2) から下方・上方テール依存係数を推定"""
    n = len(u)
    u1, u2 = u[:, 0], u[:, 1]
    lower, upper = [], []
    for q in q_grid:
        # 下方テール: P(U2 <= q | U1 <= q)
        m1 = u1 <= q
        if m1.sum() > 0:
            lower.append(np.mean(u2[m1] <= q))
        else:
            lower.append(np.nan)
        # 上方テール: P(U2 > 1-q | U1 > 1-q)
        m2 = u1 > 1 - q
        if m2.sum() > 0:
            upper.append(np.mean(u2[m2] > 1 - q))
        else:
            upper.append(np.nan)
    return np.array(lower), np.array(upper)

q_grid = np.linspace(0.005, 0.20, 30)
fig, axes = plt.subplots(1, 2, figsize=(14, 5))
for name, u in samples.items():
    low, up = empirical_tail_dependence(u, q_grid)
    axes[0].plot(q_grid, low, label=name, lw=1.8)
    axes[1].plot(q_grid, up,  label=name, lw=1.8)

# 理論値の水平線
axes[0].axhline(2**(-1/theta_clayton_05), color='C2', ls='--', alpha=0.5,
                label=f'Clayton theory lambda_L = {2**(-1/theta_clayton_05):.3f}')
axes[1].axhline(2 - 2**(1/theta_gumbel_05), color='C3', ls='--', alpha=0.5,
                label=f'Gumbel theory lambda_U = {2 - 2**(1/theta_gumbel_05):.3f}')

axes[0].set_xlabel('q (quantile)'); axes[0].set_ylabel('P(U2<=q | U1<=q)')
axes[0].set_title('Lower tail dependence  (q -> 0)')
axes[0].set_ylim(0, 1); axes[0].grid(True, alpha=0.3); axes[0].legend(fontsize=8)

axes[1].set_xlabel('q (quantile)'); axes[1].set_ylabel('P(U2>1-q | U1>1-q)')
axes[1].set_title('Upper tail dependence  (q -> 0)')
axes[1].set_ylim(0, 1); axes[1].grid(True, alpha=0.3); axes[1].legend(fontsize=8)
plt.tight_layout()
plt.savefig('tail_dependence.png', dpi=150, bbox_inches='tight')
plt.show()

このグラフから、各コピュラのテール依存性が定量的に明確になります。Gaussian(青)は $q \to 0$ で下方も上方もゼロに向かって減少し、理論通りテール独立であることが分かります。t(橙)は両側で正の値に収束し、両側テール依存が確認できます。Clayton(緑)は下方では理論値 $2^{-1/2} \approx 0.707$ に近づきますが、上方ではゼロに落ちる——下方依存だけが強い非対称性が出ています。Gumbel(赤)はその逆で、上方では理論値 $2 – 2^{1/2} \approx 0.586$ に近づき、下方ではゼロです。理論計算と数値推定が綺麗に一致していることが、コピュラ理論の予測力を裏付けています。

テール依存係数の数値推定曲線(理論値との比較)

左パネルが下方テール依存係数、右パネルが上方テール依存係数の推定曲線です。破線が理論値を示しています。ClaytonコピュラはKendallのτ=0.5のとき理論的な下方テール依存係数が $\lambda_L = 2^{-1/2} \approx 0.707$ ですが、グラフの緑線(下方パネル)はほぼその値に収束しています。GumbelコピュラはKendallのτ=0.5のとき上方テール依存係数が $\lambda_U = 2 – \sqrt{2} \approx 0.586$ で、右パネルの紫線が理論破線に近づいているのが確認できます。ガウシアンコピュラ(青線)はどちらのパネルでも $q \to 0$ でゼロに向かっており、「テール独立」の性質が数値的に実証されています。

応用 — 金融ポートフォリオのVaRと保険リスク

S&P500風シミュレーション

リスク管理で最も使われる指標が Value at Risk (VaR) です。「ポートフォリオが1日に被る最悪損失の99パーセンタイル」のような形で定義されます。同じ相関を持つ2資産ポートフォリオで、コピュラを変えるとVaRがどれだけ違うかを見てみましょう。

import numpy as np
from scipy.stats import norm, t

np.random.seed(2008)
N = 100_000

# 周辺分布: 2資産とも年率15%ボラの正規リターン(簡単化)
mu = 0.0; sigma = 0.15 / np.sqrt(252)  # 日次

# Kendallのτ = 0.5 で揃える(前のセクションと同じ)
rho = np.sin(np.pi/2 * 0.5)
theta_c = 2.0
theta_g = 2.0

# 各コピュラからサンプル → 周辺を正規にして資産リターンへ
copulas = {
    'Gaussian':  sample_gaussian_copula(rho, N),
    't (nu=4)':  sample_t_copula(rho, 4, N),
    'Clayton':   sample_clayton_copula(theta_c, N),
    'Gumbel':    sample_gumbel_copula(theta_g, N),
}

# ポートフォリオ: 等ウェイト
results = {}
for name, u in copulas.items():
    r1 = norm.ppf(u[:, 0], loc=mu, scale=sigma)
    r2 = norm.ppf(u[:, 1], loc=mu, scale=sigma)
    port_ret = 0.5 * r1 + 0.5 * r2
    var_99 = -np.percentile(port_ret, 1)     # 99% VaR(損失なので符号反転)
    var_995 = -np.percentile(port_ret, 0.5)  # 99.5% VaR
    cvar_99 = -port_ret[port_ret <= -var_99].mean()  # 期待ショートフォール
    results[name] = (var_99, var_995, cvar_99)

print(f"{'Copula':<12} {'VaR 99%':>10} {'VaR 99.5%':>12} {'CVaR 99%':>12}")
print('-' * 50)
for name, (v99, v995, cv99) in results.items():
    print(f"{name:<12} {v99*100:>9.3f}% {v995*100:>11.3f}% {cv99*100:>11.3f}%")

このコードの出力では、4つのコピュラで等ウェイトポートフォリオの99%VaRと99.5%VaRを比較しています。実行すると、周辺分布が同じ正規分布、Kendallのτが同じ0.5でも、VaRはコピュラごとに大きく違う ことが確認できます。典型的にはGaussianが最も小さく、tやClaytonでは99.5%VaRが20〜40%程度大きくなります。さらに極端な事象だけを集めたCVaRでは差がもっと開きます。CDOで使われたガウシアンコピュラがリスクを過小評価したというのは、まさにこの差のことです。

コピュラ別ポートフォリオリスク指標(VaR/CVaR)比較棒グラフ

各コピュラでの99%VaR・99.5%VaR・CVaR99% を棒グラフで並べています。ガウシアンコピュラに対してtコピュラやClaytonコピュラはCVaRで5〜10%程度大きな値を示しており(緑棒上の「+X%」注釈が対ガウシアン比を表しています)、同じ相関強度でもコピュラ選択がリスク指標に与える影響が定量的に確認できます。特にCVaR(期待ショートフォール)では、VaRより大きな差が現れており、これは「テール全体の平均損失」がコピュラ選択に敏感であることを示しています。2008年の金融危機で問題になったのは、まさにこの差——ガウシアンコピュラを使うとCVaRが実際より著しく小さく見えること——です。

VaR推定値の分布を可視化

VaRの差を視覚的に示すため、損失分布のヒストグラムを描いてみます。

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm

fig, ax = plt.subplots(1, 1, figsize=(10, 6))
colors = {'Gaussian': 'C0', 't (nu=4)': 'C1', 'Clayton': 'C2', 'Gumbel': 'C3'}

for name, u in copulas.items():
    r1 = norm.ppf(u[:, 0], loc=mu, scale=sigma)
    r2 = norm.ppf(u[:, 1], loc=mu, scale=sigma)
    port_ret = 0.5 * r1 + 0.5 * r2
    loss = -port_ret * 100  # %損失
    ax.hist(loss, bins=200, range=(-4, 6), alpha=0.4, density=True,
            label=name, color=colors[name])
    v99 = np.percentile(loss, 99)
    ax.axvline(v99, color=colors[name], ls='--', alpha=0.8)

ax.set_xlabel('Daily portfolio loss [%]')
ax.set_ylabel('Density')
ax.set_title('Portfolio loss distribution by copula choice — dashed lines: 99% VaR')
ax.set_xlim(-4, 6); ax.grid(True, alpha=0.3); ax.legend()
plt.tight_layout()
plt.savefig('portfolio_var.png', dpi=150, bbox_inches='tight')
plt.show()

このヒストグラムから、損失分布の右側(大損失側)の形がコピュラによってどう違うかが視覚化されます。Gaussianコピュラの分布は最も「軽い裾」を持ち、99%VaRの線(破線)が最も左にあります。tやClaytonコピュラでは右裾が伸び、99%VaRがより大きな値を示します。日々の通常変動では4つの分布は似ているのですが、「100日に1回」のような稀な大損失の領域では、コピュラ選択が結果を1.3〜1.5倍ほど押し上げる ——これがリスクモデルとして無視できない違いです。

コピュラ別損失分布ヒストグラムと99%VaR

4つのコピュラの損失分布(ヒストグラム)と99%VaR(破線)を重ねて表示しました。-2%〜+2%の通常変動域では4つの分布はほぼ重なっており、日常的なリスク管理では見分けがつきません。しかし、右端の大損失ゾーン(3%以上)では分布の「裾の重さ」が明らかに分かれます。ガウシアン(青)の破線VaRが最も左(小さい損失)にあるのに対し、t(橙)やClayton(緑)の破線はより右に位置し、同じ信頼水準でも推定されるリスクが大きくなっています。

copulasライブラリでのフィッティング

実データへのフィッティングには copulas ライブラリ(pip install copulas)が便利です。架空の日次株価リターンを生成して、フィッティングしてみます。

# pip install copulas
import numpy as np
import pandas as pd
from copulas.bivariate import Clayton, Gumbel, Frank
from copulas.multivariate import GaussianMultivariate

# 仮想データ: t分布から生成(テール依存があるはずのデータ)
np.random.seed(7)
true_u = sample_t_copula(0.6, nu=4, n=2000)
data = pd.DataFrame({
    'asset_A': norm.ppf(true_u[:, 0]) * 0.012,  # 日次1.2%ボラ
    'asset_B': norm.ppf(true_u[:, 1]) * 0.018,  # 日次1.8%ボラ
})

# ガウシアンコピュラでフィット
gauss_model = GaussianMultivariate()
gauss_model.fit(data)
print("Gaussian fit complete")
print(f"  Correlation matrix: \n{gauss_model.correlation}")

# 2変数Archimedeanのフィット(ランクをu, vにして)
from scipy.stats import rankdata
u_emp = rankdata(data['asset_A']) / (len(data) + 1)
v_emp = rankdata(data['asset_B']) / (len(data) + 1)

clayton = Clayton(); clayton.fit(np.column_stack([u_emp, v_emp]))
gumbel = Gumbel();   gumbel.fit(np.column_stack([u_emp, v_emp]))
frank = Frank();     frank.fit(np.column_stack([u_emp, v_emp]))

print(f"\nClayton theta = {clayton.theta:.3f}  ->  tau = {clayton.theta/(clayton.theta+2):.3f}")
print(f"Gumbel  theta = {gumbel.theta:.3f}  ->  tau = {1 - 1/gumbel.theta:.3f}")
print(f"Frank   theta = {frank.theta:.3f}")
# 経験的Kendall tau
from scipy.stats import kendalltau
tau_emp, _ = kendalltau(data['asset_A'], data['asset_B'])
print(f"\nEmpirical Kendall tau = {tau_emp:.3f}")

このフィッティング結果は、各コピュラがどう実データに「フィット」するかを示しています。経験的なKendallのτが0.4前後とすると、Claytonでは $\theta/(θ+2) = 0.4$ から $\theta \approx 1.33$、Gumbelでは $1 – 1/\theta = 0.4$ から $\theta \approx 1.67$ といった値が得られます。実用上はこれら複数のコピュラをフィットして、AICやBICで比較したり、テール部分でのフィット精度を見て選択します。

保険・気象への応用

コピュラの応用は金融に留まりません。保険数理 では、地震保険と火災保険の同時請求リスクをモデル化する際、Claytonコピュラが「大災害時の同時請求」を捉えるのに有効です。一方で 気象 では、複数の観測地点の年最大降水量の依存をGumbelコピュラでモデル化することが多く、これは「広域台風での同時極値発生」を表現します。信頼性工学 では、システム内の複数部品が共通の負荷(温度、振動等)にさらされている場合、その負荷を介した故障時間の依存をフラジリティ関数とコピュラで結合する手法が用いられます。

どの応用でも本質は同じです——個別の周辺分布は別々の物理的・経済的モデルから決まるが、それらの依存構造はコピュラで自由に組み合わせられる。Sklarの定理という1つの定理が、こんなにも幅広い応用に道を開いているのです。

まとめ

本記事では、コピュラ理論を多変量依存構造の柔軟な表現手段として、理論から実装、応用までを通して解説しました。

  • Sklarの定理 — 任意の連続多変量CDFは、周辺分布たちと、ただ一つのコピュラ $C$ に一意分解できる。これにより「周辺の形」と「依存の構造」を独立に設計可能になる。
  • コピュラの種類 — Gaussian/tは楕円対称、Clayton/Gumbel/FrankはArchimedean族で非対称な依存も表現できる。それぞれが捉える依存パターンは大きく異なる。
  • 依存の強さの指標 — Pearson相関は周辺分布の影響を受けるが、Kendallのτ・Spearmanのρはコピュラだけで決まる順位相関で、コピュラパラメータと一対一に対応する。
  • テール依存係数 — 極端事象同士の連動性を測る $\lambda_L, \lambda_U$ は、コピュラの「裾の挙動」の本質を捉える。Gaussianは $\rho < 1$ で常にテール独立であり、これが2008年金融危機でCDOのリスクが過小評価された数学的本質。
  • 実装上の含意 — 同じ周辺分布、同じKendallのτでも、コピュラの選択でVaRやCVaRは大きく変わる。リスクモデルでは、Gaussianコピュラの選択を盲信せず、データのテール挙動に合うコピュラ(特にtコピュラ)を選ぶことが重要。

コピュラ理論を学ぶ次のステップとして、以下のトピックがおすすめです。Vineコピュラ(高次元での依存構造を二変数コピュラのツリーで分解する手法)、動的コピュラ(時間とともに変化する依存構造をARモデルやDCCモデルと組み合わせる)、極値理論との接続(多変量極値分布とコピュラの関係)。いずれも本記事で学んだSklarの定理とコピュラ族の基礎が出発点になります。

次のステップとして、以下の記事も参考にしてください。

多変量ガウス分布の解析的・幾何的な特徴をまとめる
コピュラの構成要素として必須の多変量正規分布。共分散行列の意味、周辺化・条件付け、等高線の幾何的解釈を丁寧に解説します
画像なし
共分散と相関係数の定義・性質・Pythonでの計算
コピュラが克服する相関係数の限界を理解するための基礎。Pearson/Spearman/Kendall各指標の定義・性質・Pythonでの計算方法を比較します