このページは Cloud Translation API によって翻訳されました。
Switch to English

コピュラプライマー

TensorFlow.orgで見る Google Colabで実行 GitHubでソースを表示する ノートブックをダウンロード
import numpy as np
import matplotlib.pyplot as plt
import tensorflow.compat.v2 as tf
tf.enable_v2_behavior()

import tensorflow_probability as tfp
tfd = tfp.distributions
tfb = tfp.bijectors

[コピュラ](https://en.wikipedia.org/wiki/Copula_(probability_theory%29)は、確率変数間の依存関係を取得するための古典的なアプローチです。より正式には、コピュラは多変量分布$ C(U_1、U_2 、....、U_n)$マージナライズすると$ U_i \ sim \ text {Uniform}(0、1)$が得られます。

コピュラは、それらを使用して任意の周辺をもつ多変量分布を作成できるので興味深いです。これがレシピです:

  • 確率積分変換を使用すると、任意の連続RV $ X $が均一な1つの$ F_X(X)$に変わります。ここで、$ F_X $は$ X $のCDFです。
  • コピュラ(たとえば、2変量)$ C(U、V)$が与えられると、$ U $と$ V $の周辺分布は均一になります。
  • 次に、関心のあるRVが$ X、Y $であるとして、新しい分布$ C '(X、Y)= C(F_X(X)、F_Y(Y))$を作成します。 $ X $と$ Y $のマージナルは私たちが望んだものです。

マージナルは単変量であるため、測定やモデル化が簡単になる場合があります。コピュラを使用すると、境界から開始することができますが、次元間の任意の相関も実現できます。

ガウスコピュラ

コピュラがどのように構築されるかを説明するために、多変量ガウス相関に従って依存関係を取得する場合を考えます。ガウスコピュラは、$ C(u_1、u_2、... u_n)= \ Phi_ \ Sigma(\ Phi ^ {-1}(u_1)、\ Phi ^ {-1}(u_2)、... \ Phi ^ {-1}(u_n))$ここで、$ \ Phi_ \ Sigma $は、共分散$ \ Sigma $と平均0のMultivariateNormalのCDFを表し、$ \ Phi ^ {-1} $は逆CDFです標準の通常。

法線の逆CDFを適用すると、均一に分散する次元が正規分布します。次に、多変量正規のCDFを適用すると、分布が潰されて、わずかに均一になり、ガウス相関が適用されます。

したがって、ガウスコピュラは単位ハイパーキューブ$ [0、1] ^ n $上の辺縁が均一な分布です。

このように定義すると、ガウスコピュラはtfd.TransformedDistributionと適切なBijector実装できます。つまり、 tfb.NormalCDFバイジェクターによって実装された正規分布の逆CDFを使用して、MultivariateNormalを変換しています。

以下では、1つの単純化の仮定を使用してガウスコピュラを実装します。共分散はコレスキー因子によってパラメーター化される(したがってMultivariateNormalTriL共分散)。 (他のtf.linalg.LinearOperatorsを使用して、行列のないさまざまな仮定をエンコードできます。)

class GaussianCopulaTriL(tfd.TransformedDistribution):
  """Takes a location, and lower triangular matrix for the Cholesky factor."""
  def __init__(self, loc, scale_tril):
    super(GaussianCopulaTriL, self).__init__(
        distribution=tfd.MultivariateNormalTriL(
            loc=loc,
            scale_tril=scale_tril),
        bijector=tfb.NormalCDF(),
        validate_args=False,
        name="GaussianCopulaTriLUniform")


# Plot an example of this.
unit_interval = np.linspace(0.01, 0.99, num=200, dtype=np.float32)
x_grid, y_grid = np.meshgrid(unit_interval, unit_interval)
coordinates = np.concatenate(
    [x_grid[..., np.newaxis],
     y_grid[..., np.newaxis]], axis=-1)

pdf = GaussianCopulaTriL(
    loc=[0., 0.],
    scale_tril=[[1., 0.8], [0., 0.6]],
).prob(coordinates)

# Plot its density.

plt.contour(x_grid, y_grid, pdf, 100, cmap=plt.cm.jet);

png

ただし、このようなモデルの能力は、確率積分変換を使用して任意のRVでコピュラを使用することです。このようにして、任意の周辺を指定し、コピュラを使用してそれらをつなぎ合わせることができます。

モデルから始めます:

$$\begin{align*} X &\sim \text{Kumaraswamy}(a, b) \\ Y &\sim \text{Gumbel}(\mu, \beta) \end{align*}$$

そして、コピュラを使用して、限界のクマラスワミーガンベルを持つ2変量RV $ Z $を取得します。

まず、これら2つのRVによって生成された製品分布をプロットすることから始めます。これは、コピュラを適用するときの比較ポイントとして機能するためです。

a = 2.0
b = 2.0
gloc = 0.
gscale = 1.

x = tfd.Kumaraswamy(a, b)
y = tfd.Gumbel(loc=gloc, scale=gscale)

# Plot the distributions, assuming independence
x_axis_interval = np.linspace(0.01, 0.99, num=200, dtype=np.float32)
y_axis_interval = np.linspace(-2., 3., num=200, dtype=np.float32)
x_grid, y_grid = np.meshgrid(x_axis_interval, y_axis_interval)

pdf = x.prob(x_grid) * y.prob(y_grid)

# Plot its density

plt.contour(x_grid, y_grid, pdf, 100, cmap=plt.cm.jet);

png

限界が異なる共同分布

次に、ガウスコピュラを使用して分布を結合し、それをプロットします。ここでも、選択したツールは、適切なBijectorを適用して選択された限界を取得するTransformedDistributionです。

具体的には、ベクトルの異なる部分に異なるBlockwiseイジェクターを適用するBlockwiseワイズBlockwiseイジェクターを使用します(これはまだ全単射変換です)。

これで、必要なコピュラを定義できます。ターゲットバイマージナル(バイジェクターとしてエンコードされた)のリストを指定すると、コピュラを使用し、指定されたマージナルを持つ新しい分布を簡単に構築できます。

class WarpedGaussianCopula(tfd.TransformedDistribution):
  """Application of a Gaussian Copula on a list of target marginals.

  This implements an application of a Gaussian Copula. Given [x_0, ... x_n]
  which are distributed marginally (with CDF) [F_0, ... F_n],
  `GaussianCopula` represents an application of the Copula, such that the
  resulting multivariate distribution has the above specified marginals.

  The marginals are specified by `marginal_bijectors`: These are
  bijectors whose `inverse` encodes the CDF and `forward` the inverse CDF.

  block_sizes is a 1-D Tensor to determine splits for `marginal_bijectors`
  length should be same as length of `marginal_bijectors`.
  See tfb.Blockwise for details
  """
  def __init__(self, loc, scale_tril, marginal_bijectors, block_sizes=None):
    super(WarpedGaussianCopula, self).__init__(
        distribution=GaussianCopulaTriL(loc=loc, scale_tril=scale_tril),
        bijector=tfb.Blockwise(bijectors=marginal_bijectors,
                               block_sizes=block_sizes),
        validate_args=False,
        name="GaussianCopula")

最後に、実際にこのガウスコピュラを使用してみましょう。 $ \ begin {bmatrix} 1&0 \\ rho&\ sqrt {(1- \ rho ^ 2)} \ end {bmatrix} $のコレスキーを使用します。これは、分散1と相関$ \に対応します多変量正規のrho $。

いくつかのケースを見てみましょう。

# Create our coordinates:
coordinates = np.concatenate(
    [x_grid[..., np.newaxis], y_grid[..., np.newaxis]], -1)


def create_gaussian_copula(correlation):
  # Use Gaussian Copula to add dependence.
  return WarpedGaussianCopula(
      loc=[0.,  0.],
      scale_tril=[[1., 0.], [correlation, tf.sqrt(1. - correlation ** 2)]],
      # These encode the marginals we want. In this case we want X_0 has
      # Kumaraswamy marginal, and X_1 has Gumbel marginal.

      marginal_bijectors=[
          tfb.Invert(tfb.KumaraswamyCDF(a, b)),
          tfb.Invert(tfb.GumbelCDF(loc=0., scale=1.))])


# Note that the zero case will correspond to independent marginals!
correlations = [0., -0.8, 0.8]
copulas = []
probs = []
for correlation in correlations:
  copula = create_gaussian_copula(correlation)
  copulas.append(copula)
  probs.append(copula.prob(coordinates))


# Plot it's density

for correlation, copula_prob in zip(correlations, probs):
  plt.figure()
  plt.contour(x_grid, y_grid, copula_prob, 100, cmap=plt.cm.jet)
  plt.title('Correlation {}'.format(correlation))

png

png

png

最後に、必要な限界を実際に取得していることを確認しましょう。

def kumaraswamy_pdf(x):
    return tfd.Kumaraswamy(a, b).prob(np.float32(x))

def gumbel_pdf(x):
    return tfd.Gumbel(gloc, gscale).prob(np.float32(x))


copula_samples = []
for copula in copulas:
  copula_samples.append(copula.sample(10000))

plot_rows = len(correlations)
plot_cols = 2  # for 2  densities [kumarswamy, gumbel]
fig, axes = plt.subplots(plot_rows, plot_cols, sharex='col', figsize=(18,12))

# Let's marginalize out on each, and plot the samples.

for i, (correlation, copula_sample) in enumerate(zip(correlations, copula_samples)):
  k = copula_sample[..., 0].numpy()
  g = copula_sample[..., 1].numpy()


  _, bins, _ = axes[i, 0].hist(k, bins=100, normed=True)
  axes[i, 0].plot(bins, kumaraswamy_pdf(bins), 'r--')
  axes[i, 0].set_title('Kumaraswamy from Copula with correlation {}'.format(correlation))

  _, bins, _ = axes[i, 1].hist(g, bins=100, normed=True)
  axes[i, 1].plot(bins, gumbel_pdf(bins), 'r--')
  axes[i, 1].set_title('Gumbel from Copula with correlation {}'.format(correlation))
  

png

結論

さあ、行きましょう! Bijector APIを使用してガウスコピュラを構築できることを示しました。

より一般的には、 Bijector APIを使用してBijectorイジェクターを記述し、それらをディストリビューションで構成すると、柔軟なモデリングのための豊富なディストリビューションファミリーを作成できます。