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
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 |
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] |
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] |
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
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 |
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] |
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
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 |
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] |
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]]])