A Simple Copula-GARCH Example

In this example, we will load a dataset which contains returns from 3 ETF and attempt to simulate future returns. Instead of fitting a multivariate GARCH model, what we will do instead is to fit a univariate GARCH model to each returns stream and construct a dependency model among these returns streams with a copula.

The copulas package can be installed separately at

conda install -c conda-forge copulae  # for anaconda

pip install copulae  # pip

Model Overview.

We will assume an AR(1)-GARCH(1, 1)-Normal model for each returns stream. For their dependency structure, we will assume a Student (T) copula.

[1]:
from muarch import MUArch, UArch
from muarch.datasets import load_etf
from copulae import TCopula


returns = load_etf()  # load returns data
returns.head()
[1]:
VOO EEM VT
Date
2010-10-01 0.043390 0.030154 0.037955
2010-11-01 -0.001108 -0.029055 -0.026458
2010-12-01 0.064337 0.063868 0.056120
2011-01-01 0.026914 -0.030360 0.033688
2011-02-01 0.034664 -0.000436 0.029297
[2]:
num_assets = returns.shape[1]

# sets up a MUArch model collection where each model defaults to
# mean: AR(1)
# vol: GARCH(1, 1)
# dist: normal
models = MUArch(num_assets, mean='AR', lags=1)

We could overwrite each model in the MUArch instance. For example, if we believe that a skew-t distribution better describes the innovation of VOO (an ETF tracking S&P 500), we could overwrite it as follows.

[3]:
# set first model to AR(1)-GARCH(1, 1) with skewt innovations
models[0] = UArch('AR', lags=1, dist='skewt')

In fact, we could set all of the models separately. All we have to do is to call MUArch like this

models = MUArch(5)  # 5 models

for i in range(5):
    models[i] = make_uarch_model(...)

To fit the model, we just need to call the .fit() method. This applies to both the UArch and MUArch models

[4]:
models.fit(returns)

We can see the summary of the models using the .summary() method.

[5]:
models.summary()
[5]:

VOO

AR - GARCH Model Results
Dep. Variable: y R-squared: 0.003
Mean Model: AR Adj. R-squared: -0.007
Vol Model: GARCH Log-Likelihood: 209.083
Distribution: Standardized Skew Student's t AIC: -404.166
Method: Maximum Likelihood BIC: -385.860
No. Observations: 101
Date: Mon, Mar 18 2019 Df Residuals: 94
Time: 14:12:56 Df Model: 7
Mean Model
coef std err t P>|t| 95.0% Conf. Int.
Const 0.0158 3.120e-03 5.048 4.466e-07 [9.635e-03,2.187e-02]
y[1] -0.3131 7.794e-02 -4.017 5.884e-05 [ -0.466, -0.160]
Volatility Model
coef std err t P>|t| 95.0% Conf. Int.
omega 2.6889e-04 1.293e-04 2.079 3.759e-02 [1.543e-05,5.224e-04]
alpha[1] 0.2066 7.759e-02 2.663 7.735e-03 [5.458e-02, 0.359]
beta[1] 0.5258 0.131 4.025 5.702e-05 [ 0.270, 0.782]
Distribution
coef std err t P>|t| 95.0% Conf. Int.
nu 37.8296 132.062 0.286 0.775 [-2.210e+02,2.967e+02]
lambda -0.4572 0.110 -4.172 3.021e-05 [ -0.672, -0.242]


Covariance estimator: robust




EEM

AR - GARCH Model Results
Dep. Variable: y R-squared: 0.006
Mean Model: AR Adj. R-squared: -0.004
Vol Model: GARCH Log-Likelihood: 157.484
Distribution: Normal AIC: -304.968
Method: Maximum Likelihood BIC: -291.893
No. Observations: 101
Date: Mon, Mar 18 2019 Df Residuals: 96
Time: 14:12:56 Df Model: 5
Mean Model
coef std err t P>|t| 95.0% Conf. Int.
Const 5.0854e-03 5.935e-03 0.857 0.392 [-6.546e-03,1.672e-02]
y[1] -0.0628 0.109 -0.579 0.563 [ -0.275, 0.150]
Volatility Model
coef std err t P>|t| 95.0% Conf. Int.
omega 2.1956e-03 9.297e-04 2.362 1.820e-02 [3.734e-04,4.018e-03]
alpha[1] 0.1735 0.139 1.246 0.213 [-9.935e-02, 0.446]
beta[1] 0.0000 0.330 0.000 1.000 [ -0.647, 0.647]


Covariance estimator: robust




VT

AR - GARCH Model Results
Dep. Variable: y R-squared: 0.015
Mean Model: AR Adj. R-squared: 0.005
Vol Model: GARCH Log-Likelihood: 194.275
Distribution: Normal AIC: -378.551
Method: Maximum Likelihood BIC: -365.475
No. Observations: 101
Date: Mon, Mar 18 2019 Df Residuals: 96
Time: 14:12:56 Df Model: 5
Mean Model
coef std err t P>|t| 95.0% Conf. Int.
Const 0.0114 4.305e-03 2.650 8.039e-03 [2.972e-03,1.985e-02]
y[1] -0.1711 8.939e-02 -1.914 5.559e-02 [ -0.346,4.086e-03]
Volatility Model
coef std err t P>|t| 95.0% Conf. Int.
omega 3.2937e-04 4.396e-04 0.749 0.454 [-5.322e-04,1.191e-03]
alpha[1] 0.2420 0.189 1.280 0.201 [ -0.129, 0.613]
beta[1] 0.5257 0.471 1.117 0.264 [ -0.397, 1.448]


Covariance estimator: robust

Now that we have a model for each of the returns streams, the question is how can we create a dependency model amongst them? We can do so by fitting the residuals for each UArch model into a copula.

This fitting of residuals means for the copula to find a relationship among the different models. Subsequently, we can use the copula to randomly generate the residuals which can be used by the UArch model to simulate returns stream for the assets.

[6]:
residuals = models.residuals() # defaults to return the standardized residuals

cop = TCopula(dim=num_assets)
cop.fit(residuals)

print(cop.summary())
Student T Copula with 3 dimensions

Degrees of Freedom: 9.837148817580086

Correlation Matrix (P):
    [[1.         0.59016153 0.90300841]
 [0.59016153 1.         0.80087304]
 [0.90300841 0.80087304 1.        ]]

Log. Lik        : -149.63318356908047
Var. Est.       : Not Implemented Yet
Method          : Maximum pseudo-likelihood
Data Pts.       : 101

Optim Options
        bounds         : [(0.0, inf), (-1.0, 1.0), (-1.0, 1.0), (-1.0, 1.0)]
        options        : {'maxiter': 20000, 'ftol': 1e-06, 'iprint': 1, 'disp': False, 'eps': 1.5e-08}
        method         : SLSQP

Results
        x              : [9.83714882 0.59016153 0.90300841 0.80087304]
        fun            : -149.63318356908047
        jac            : [-1.70530257e-05 -1.11034145e-03  2.97859515e-03  1.59540529e-03]
        nit            : 13
        nfev           : 84
        njev           : 13
        status         : 0
        message        : Optimization terminated successfully.
        success        : True

We could of course overwrite the correlation matrix of the TCopula with the historical correlation. I’ll show an example below. For more details, check out the Copulae package documentation.

[7]:
cop[:] = returns.corr()

print(cop.summary())
Student T Copula with 3 dimensions

Degrees of Freedom: 9.837148817580086

Correlation Matrix (P):
    [[1.         0.72108633 0.94803336]
 [0.72108633 1.         0.85096406]
 [0.94803336 0.85096406 1.        ]]

Log. Lik        : -149.63318356908047
Var. Est.       : Not Implemented Yet
Method          : Maximum pseudo-likelihood
Data Pts.       : 101

Optim Options
        bounds         : [(0.0, inf), (-1.0, 1.0), (-1.0, 1.0), (-1.0, 1.0)]
        options        : {'maxiter': 20000, 'ftol': 1e-06, 'iprint': 1, 'disp': False, 'eps': 1.5e-08}
        method         : SLSQP

Results
        x              : [9.83714882 0.59016153 0.90300841 0.80087304]
        fun            : -149.63318356908047
        jac            : [-1.70530257e-05 -1.11034145e-03  2.97859515e-03  1.59540529e-03]
        nit            : 13
        nfev           : 84
        njev           : 13
        status         : 0
        message        : Optimization terminated successfully.
        success        : True

Notice the difference in the correlation matrix above. Note that correlation matrix is a parameter only for elliptical copulas. For Archimedean and others, change the parameters accordingly.

Let us now simulate 100 trials of returns, 10 steps into the future with the copula. All you have to do is to pass in the .random method to generate the innovations.

[8]:
horizon = 10
trials = 100

models.simulate_mc(horizon, trials, custom_dist=cop.random)
[8]:
array([[[ 2.36819180e-02,  5.91736570e-02,  4.76506242e-02],
        [ 1.44822161e-02,  6.27440217e-02,  1.37890054e-02],
        [-5.48242994e-02, -7.61327302e-02, -4.33654769e-02],
        ...,
        [ 2.16221301e-02, -3.66574706e-03,  5.80400414e-03],
        [ 1.24950196e-02, -6.65987110e-02, -1.91601512e-03],
        [ 4.35219595e-02,  1.42134542e-01,  1.18264754e-01]],

       [[ 1.14673295e-02, -7.24770972e-03,  1.12111664e-03],
        [ 1.15180971e-02,  3.51954043e-02,  1.45830445e-02],
        [ 5.44165821e-02,  2.19578347e-02,  3.13810703e-02],
        ...,
        [ 2.89393067e-02,  5.70089775e-02,  3.36656830e-02],
        [ 4.99865506e-04,  4.40031966e-02,  8.95130002e-05],
        [ 6.15984202e-02,  8.85007780e-02,  1.28838163e-01]],

       [[ 2.39743360e-02,  3.34112771e-02,  3.14031296e-02],
        [ 2.92497275e-02, -2.05959398e-03,  1.59909456e-02],
        [-8.70354146e-02, -8.66656129e-02, -6.46162277e-02],
        ...,
        [-6.16685073e-03, -2.02064819e-02, -7.07666634e-03],
        [ 2.26282761e-02, -3.40303169e-03,  1.63187322e-02],
        [ 1.32278891e-02,  2.03830411e-02,  1.30598931e-02]],

       ...,

       [[ 1.38241116e-02, -1.77706225e-02,  3.15010579e-03],
        [ 4.07693772e-02,  3.72229328e-02,  3.48073220e-02],
        [-5.15138745e-02, -3.47109312e-02, -3.92909818e-02],
        ...,
        [ 1.99715299e-02,  8.51965041e-02,  4.27344365e-02],
        [-1.78759597e-02,  6.98907216e-03, -1.19619336e-02],
        [ 2.82702866e-02,  2.89281519e-02,  3.18538287e-02]],

       [[ 1.41677467e-02, -1.76165401e-02,  1.20177353e-02],
        [-8.85086028e-02, -1.42317567e-01, -8.50680186e-02],
        [ 1.49899395e-02, -1.54490974e-02, -4.67326285e-03],
        ...,
        [ 2.33224844e-03, -5.60107187e-02, -2.00618802e-02],
        [ 1.67808821e-02,  3.39528074e-02,  2.07668467e-02],
        [ 7.08545343e-03,  3.61847390e-02,  1.77982130e-02]],

       [[ 1.70144077e-02, -2.83423287e-02,  9.56741286e-03],
        [ 2.73957468e-02, -1.38556690e-01, -3.99365353e-02],
        [ 4.69939240e-02,  2.94684632e-02,  4.31496651e-02],
        ...,
        [ 4.15675721e-02,  5.30936408e-02,  3.20754793e-02],
        [ 2.85645038e-02,  9.01270560e-02,  4.26091035e-02],
        [ 3.27343634e-02,  1.04311493e-02,  2.55611539e-02]]])