Last active
October 31, 2016 01:15
-
-
Save BioSciEconomist/c767496a959d623650860dece4cf5222 to your computer and use it in GitHub Desktop.
Modeling Dependence with Copulas and Quantmod in R
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
# ------------------------------------------------------------------ | |
# |PROGRAM NAME: copual and quantmod ex | |
# |DATE: 10/30/16 | |
# |CREATED BY: MATT BOGARD | |
# |PROJECT FILE: | |
# |---------------------------------------------------------------- | |
# | PURPOSE: R code supporting: http://econometricsense.blogspot.com/2016/10/modeling-dependence-with-copulas-and.html | |
# | | |
# | REFERENCE: Modelling Dependence with Copulas in R | |
# | http://datascienceplus.com/modelling-dependence-with-copulas/ | |
# |------------------------------------------------------------------ | |
library(copula) | |
# use R quantmod to get returns | |
library(quantmod) # load quantmod | |
getSymbols("MON") | |
#get last 3 years of data | |
mon_dat <- last(MON, '36 months') | |
head(mon_dat) # view data | |
# calculate returns | |
mon_ret <- as.data.frame(dailyReturn(mon_dat)) # returns by day as an R data frame | |
head(mon_ret) # view data | |
dim(mon_ret) # check dimensions | |
# get returns for adm | |
getSymbols("ADM") # get stock quote data | |
#get last 3 years of data | |
adm_dat <- last(ADM, '36 months') | |
head(adm_dat) # view data | |
# calculate returns | |
adm_ret <- as.data.frame(dailyReturn(adm_dat)) # returns by day as an R data frame | |
print(adm_ret$daily.returns) | |
head(adm_ret) # view file | |
dim(adm_ret) # check dimensions | |
# Plot and correlation checks | |
plot(mon_ret$daily.returns,adm_ret$daily.returns,pch='.') | |
abline(lm(adm_ret$daily.returns~mon_ret$daily.returns),col='red',lwd=1) | |
cor(mon_ret$daily.returns,adm_ret$daily.returns,method='spearman') | |
# Fit a t-copula | |
t.cop <- tCopula(dim=2) | |
set.seed(500) | |
m <- pobs(as.matrix(cbind(mon_ret$daily.returns,adm_ret$daily.returns))) | |
fit <- fitCopula(t.cop,m,method='ml') | |
coef(fit) | |
# Set the parameters | |
rho <- coef(fit)[1] | |
df <- coef(fit)[2] | |
# Plot the density | |
persp(tCopula(dim=2,rho,df=df),dCopula) | |
# Sample random observation from the copula | |
u <- rCopula(3965,tCopula(dim=2,rho,df=df)) | |
plot(u[,1],u[,2],pch='.',col='blue') | |
cor(u,method='spearman') | |
# Estimate marginals parameters | |
mon_mu <- mean(mon_ret$daily.returns) | |
mon_sd <- sd(mon_ret$daily.returns) | |
adm_mu <- mean(adm_ret$daily.returns) | |
adm_sd <- sd(adm_ret$daily.returns) | |
# Build the distribution from the copula | |
copula_dist <- mvdc(copula=tCopula(rho,dim=2,df=df), margins=c("norm","norm"), | |
paramMargins=list(list(mean=mon_mu, sd=mon_sd), | |
list(mean=adm_mu, sd=adm_sd))) | |
# Sample from the distribution | |
sim <- rmvdc(copula_dist, 3965) | |
# Compare observed and simulated values | |
plot(mon_ret$daily.returns,adm_ret$daily.returns,main='Returns') | |
points(sim[,1],sim[,2],col='red') | |
legend('bottomright',c('Observed','Simulated'),col=c('black','red'),pch=21) | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment