If It Bleeds, We Can Kill It

Introduction

Welcome to the first post in my series on copulas in Stan. After StanCon 2024 I was inspired to start writing short blog posts about this both to help get other started and also because I often don’t really know how a thing works until I have to write about it or present it.

If you’ve ever felt intimidated by Sklar’s theorem or how the Frank Copula is defined

\[ C_\theta^F(\mathbf u) = -\frac1\theta\log\left(1 + \frac{(e^{-\theta u_1} - 1)(e^{-\theta u_2} - 1)}{e^{-\theta} - 1}\right), \]

just remember Arnold’s famous words from Predator

What Are Copulas?

At their core, copulas are functions that link univariate marginal distribution functions to form a multivariate distribution. According to Sklar’s Theorem, any multivariate joint distribution can be expressed in terms of its marginals and a copula that captures the dependence structure between variables.

If we let \(X = (X_1, \dots, X_D)\) be a multivariate random variable with marginal distribution functions \(F_i\), the joint distribution function of \(X\) can be written

\[ H(X) = C\left( F_1(X_1), \dots, F_D(X_D) \right) \]

where:

• \(H(X)\) is the joint cumulative distribution function (CDF) of the collection of random variables \(X\).
• \(F_i(X_i)\) are the marginal CDFs of each variate.
• \(C\) is the copula function.

Copulas as Densities

Copulas are multivariate distribution functions for random variables with uniform marginal distributions, i.e. they are functions that map the unit cube \([0,1]^D\) to \([0,1]\). They can also be described using copula density functions when the marginals are continuous. If \(H(X)\) is the CDF of \(X\), and the multivariate distribution has a PDF, \(h\), we write

\[ h(X) = c\left(F_1(X_1), \dots, F_D(X_D)\right) \prod_{i=1}^D f_i(x_i), \]

where \(c\) is the density of the copula. Most often, we’d model the log of the PDF

\[ \log h(X) = \log c\left(F_1(X_1), \dots, F_D(X_D)\right) + \sum_{i=1}^D \log f_i(X_i). \]

Notice that \(\sum_{i=1}^D \log f_i(X_i)\) is just the usual sum over marginal log-densities. Let’s rewrite the other term a little bit and explicitly write the parameters we’re conditioning on

\[ \begin{aligned} \log h(X) &= \log c\left(u_1, \dots, u_D \vert \theta_{c}\right) + \sum_{i=1}^D \log f_i(X_i \vert \theta_i) \\ u_i &= F_i(X_i \vert \theta_i) \end{aligned} \]

The main difference when modeling with a copula is

1. We need to use the CDFs \(F_i(X_i \vert \theta_i)\) as well as the pdfs.
2. We need to code up some function \(\log c\left(u_1, \dots, u_D\vert \theta_c\right)\) that takes as input the data \(X\) after it’s been transformed to \([0,1]^D\) by our CDFs and outputs a density.

The Copula We All Use

The simplest copula is the independence copula where we simply multiply together the uniform variates:

\[ \begin{aligned} C(\mathbf{u}) &= \prod_{i=1}^D F(X_i\vert \theta_i) = \prod_{i=1}^D u_i \\ c(\mathbf{u}) &= 1 \\ \log c(\mathbf{u}) &= 0 \end{aligned} \]

We see that if we use the independence copula, we just end up with the usual likelihood

\[ \begin{aligned} \log h(X) &= 0 + \sum_{i=1}^D \log f_i(X_i \vert \theta_i) \\ &= \log f_i(X_i \vert \theta_i) \end{aligned} \]

In this way, we all use copulas whether we want to or not!

An Imaginary Stan Model

The Stan code below is just to give an idea of what a barebones model that uses copulas might look like in Stan. In future posts I’ll write models that are based on this blueprint to implement different types of copulas.

The code is basically a simple implementation of this equation from above:

\[ \begin{aligned} \log h(X) &= \log c\left(u_1, \dots, u_D \vert \theta_{c}\right) + \sum_{i=1}^D \log f_i(X_i \vert \theta_i) \\ u_i &= F_i(X_i \vert \theta_i), \end{aligned} \]

except we’re allowing for \(N\) realisations of the \(D\)-dimensional random variable \(X\), so the input data becomes an \(N \times D\) matrix.

functions {
  real log_c(vector u, vector theta) {
    return copula_log_density
  }

  real marginal_lpdf(real x, vector theta) {
    return marginal_log_density
  }

  real marginal_cdf(real x, vector theta) {
    return uniform_variable
  }
}

data {
  int<lower = 0> N;
  int<lower = 0> D;
  matrix[n, d] X;
}

parameters {
  vector[marginal_params] theta_marginal;
  vector[copula_params] theta_copula;
}

model {
  matrix[N, D] U;
  for (i in 1:N) {
    for (j in 1:D) {
      target += marginal_lpdf(X[i, j] | theta_marginal)
      U[i, j] = marginal_cdf(X[i, j] | theta_marginal)
    }
    target += log_c(U[i, ] | theta_copula);
  }
}

Reading Material

For general purpose reading material on copulas, I recommend the following books:

• Elements of Copula Modeling with R (2018)
• An Introduction to Copulas (2006)
• Extremes in Nature: An Approach Using Copulas (2007)
• Dependence Modeling with Copulas (2014)

Coming up

The next post will be about the Gaussian copula with a couple of examples in Stan. I’ll aim to publish a new post every week.