r_rapport.Rmd

---
title: "SIMULATION DE CRISE ET IMPACTS SUR UN PORTEFEUILLE D'ACTIFS"
output:
  slidy_presentation: default
  fontsize: 5pt
  geometry: margin=1in
  ioslides_presentation: default
  beamer_presentation: default
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = FALSE)
```

## Données

Nous avons utlisé un portefeuille de 10 actifs pour un total de 100 actions avec les pondérations suivantes:

1. Air Liquide : 9,97% 
2. Total : 18,40%
3. Axa: 12,05%
4. Bnp Paribas : 22,63% 
5. Société Générale : 1,85% 
6. Vinci : 9,79%
7. Compagnie Saint Gobin : 1,92% 
8. LMVH : 2,87%
9. Orange : 3,48%
10. Crédit Agricole : 17,03%

Données : Cours moyen journalier de l'action entre 01/01/2010 et 31/12/2015. (yahoo finance)

```{r cars1, cache=FALSE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center", message=FALSE}
#Chargement des librairies
library('ggplot2')
library('PerformanceAnalytics')
library('zoo')
library('xts')
library('tseries')
library('corrplot')
library('matrixStats')
library('dygraphs')
library('reshape2')

setwd("F:/FEDORA_WINDOWS10/data_stress_test/journaliers")
# Chargement des données

df_airliquide=read.csv('airliquide_jour_2010_2015.csv')

df_axa=read.csv('axa_jour_2010_2015.csv')

df_bnp=read.csv('bnp_jour_2010_2015.csv')

df_compagniesaintgobin=read.csv('compagniesaintgobin_jour_2010_2015.csv')

df_lmvh=read.csv('lmvh_jour_2010_2015.csv')

df_orange=read.csv('orange_jour_2010_2015.csv')

df_sg=read.csv('societegenerale_jour_2010_2015.csv')

df_total=read.csv('total_jour_2010_2015.csv')

df_vinci=read.csv('vinci_jour_2010_2015.csv')

df_creditagricole=read.csv('creditagricole_jour_2010_2015.csv')

date2016 = read.csv('date2016.csv')

# Nombre total d'actifs
nb_actifs= 100

# le portefeuille des actifs; 10 actifs différents.
# Axa,bnp,sg,credit agricole,orange,total,lmvh,csg,vinci,air liquide
portefeuille_parts = c(0.20,0.18,0.18,0.11,0.05,0.099,0.01,0.02,0.12,0.04)

# Portefeuille en nombre d'actifs

portefeuille_nb_actifs = portefeuille_parts*nb_actifs

# Valeurs historiques des actions dans le pportefeuille
#Après avoir chargé les données on ne retient que deux colonnes, la date et la colonne valeur qui est la moyenne des High et Low

df_vairliquide=data.frame(row.names =as.Date(df_airliquide$Date),Valeur=(df_airliquide$High+df_airliquide$Low)/2)

df_vaxa=data.frame(row.names=as.Date(df_axa$Date),Valeur=(df_axa$High+df_axa$Low)/2)

df_vbnp=data.frame(row.names=as.Date(df_bnp$Date),Valeur=(df_bnp$High+df_bnp$Low)/2)

df_vcompagniesaintgobin=data.frame(row.names=as.Date(df_compagniesaintgobin$Date),Valeur=(df_compagniesaintgobin$High+df_compagniesaintgobin$Low)/2)

df_vlmvh=data.frame(row.names=as.Date(df_lmvh$Date),Valeur=(df_lmvh$High+df_lmvh$Low)/2)

df_vorange=data.frame(row.names=as.Date(df_orange$Date),Valeur=(df_orange$High+df_orange$Low)/2)

df_vsg=data.frame(row.names=as.Date(df_sg$Date),Valeur=(df_sg$High+df_sg$Low)/2)

df_vtotal=data.frame(row.names=as.Date(df_total$Date),Valeur=(df_total$High+df_total$Low)/2)

df_vvinci=data.frame(row.names=as.Date(df_vinci$Date),Valeur=(df_vinci$High+df_vinci$Low)/2)

df_vcreditagricole=data.frame(row.names=as.Date(df_creditagricole$Date),Valeur=(df_creditagricole$High+df_creditagricole$Low)/2)

# historique du Portefeuille en terme de valeur total des actifs...
portefeuille_historique=data.frame(row.names =as.Date(df_airliquide$Date),
                                   Valeur=df_vairliquide$Valeur*portefeuille_nb_actifs[10]+
                                     df_vaxa$Valeur*portefeuille_nb_actifs[1]+
                                     df_vbnp$Valeur*portefeuille_nb_actifs[2]+
                                     df_vcompagniesaintgobin$Valeur*portefeuille_nb_actifs[8]+
                                     df_vlmvh$Valeur*portefeuille_nb_actifs[7]+
                                     df_vorange$Valeur*portefeuille_nb_actifs[5]+
                                     df_vsg$Valeur*portefeuille_nb_actifs[3]+
                                     df_vtotal$Valeur*portefeuille_nb_actifs[6]+
                                     df_vvinci$Valeur*portefeuille_nb_actifs[9]+
                                     df_vcreditagricole$Valeur*portefeuille_nb_actifs[4])

# Graphique historique du portefeuille ....
```

## Evolution du cours historique du portefeuille

```{r cars2, cache=FALSE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
dygraph(as.xts(portefeuille_historique), main = "Evolution historique du cours du PF", 
        ylab = "Cours du PF")

```


```{r cars3, cache=FALSE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}

#Poids 
vecteur_poids=data.frame(airliquide=df_vairliquide$Valeur[1]*portefeuille_nb_actifs[10]/portefeuille_historique$Valeur[1],
                         axa=df_vaxa$Valeur[1]*portefeuille_nb_actifs[1]/portefeuille_historique$Valeur[1],
                         bnp=df_vbnp$Valeur[1]*portefeuille_nb_actifs[2]/portefeuille_historique$Valeur[1],
                         csg=df_vcompagniesaintgobin$Valeur[1]*portefeuille_nb_actifs[8]/portefeuille_historique$Valeur[1],
                         crediagricole=df_vcreditagricole$Valeur[1]*portefeuille_nb_actifs[4]/portefeuille_historique$Valeur[1],
                         lmvh=df_vlmvh$Valeur[1]*portefeuille_nb_actifs[7]/portefeuille_historique$Valeur[1],
                         orange=df_vorange$Valeur[1]*portefeuille_nb_actifs[5]/portefeuille_historique$Valeur[1],
                         sg=df_vsg$Valeur[1]*portefeuille_nb_actifs[3]/portefeuille_historique$Valeur[1],
                         total=df_vtotal$Valeur[1]*portefeuille_nb_actifs[6]/portefeuille_historique$Valeur[1],
                         vinci= df_vvinci$Valeur[1]*portefeuille_nb_actifs[9]/portefeuille_historique$Valeur[1])

vecteur_poids = c(vecteur_poids$airliquide[1],
                  vecteur_poids$axa[1],
                  vecteur_poids$bnp[1],
                  vecteur_poids$csg[1],
                  vecteur_poids$crediagricole[1],
                  vecteur_poids$lmvh[1],
                  vecteur_poids$orange[1],
                  vecteur_poids$sg[1],
                  vecteur_poids$total[1],
                  vecteur_poids$vinci[1])


# Rendement des actifs
r_airliquide = na.omit(Return.calculate(df_vairliquide))

r_axa = na.omit(Return.calculate(df_vaxa))

r_bnp = na.omit(Return.calculate(df_vbnp))

r_csg = na.omit(Return.calculate(df_vcompagniesaintgobin))

r_crediagricole = na.omit(Return.calculate(df_vcreditagricole))

r_lmvh = na.omit(Return.calculate(df_vlmvh))

r_orange = na.omit(Return.calculate(df_vorange))

r_sg = na.omit(Return.calculate(df_vsg))

r_total = na.omit(Return.calculate(df_vsg))

r_vinci =na.omit(Return.calculate(df_vvinci))

# Matrice des rendements

matrice_rendement = data.frame(row.names = as.Date(row.names(r_airliquide)),airliquide=r_airliquide$Valeur,
                               axa=r_axa$Valeur, bnp=r_bnp$Valeur, csg=r_csg$Valeur, crediagricole = r_crediagricole$Valeur
                               ,lmvh=r_lmvh$Valeur,orange=r_orange$Valeur,sg=r_sg$Valeur,total=r_total$Valeur,
                               vinci=r_vinci$Valeur)

# Rendement historique du portefeuille

rendement_histo_port = as.data.frame(Return.portfolio(matrice_rendement,vecteur_poids))
```

## Evolution du Rendement historique du portefeuille d'actif

```{r cars4, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}

dygraph(as.xts(rendement_histo_port), main = "Evolution historique du rendement du PF", 
        ylab = "Rendement du PF")
```

## Ratio de sharpe annualisé et MaxDrowDown

```{r cars5, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}

rendement_moy_port = (mean(r_airliquide$Valeur)*vecteur_poids[1]+mean(x = r_axa$Valeur)*vecteur_poids[2]+
  mean(r_bnp$Valeur)*vecteur_poids[3]+
  mean(r_csg$Valeur)*vecteur_poids[4]+
  mean(r_crediagricole$Valeur)*vecteur_poids[5]+
  mean(r_lmvh$Valeur)*vecteur_poids[6]+
  mean(r_orange$Valeur)*vecteur_poids[7]+
  mean(r_sg$Valeur)*vecteur_poids[8]+
  mean(r_total$Valeur)*vecteur_poids[9]+
  mean(r_vinci$Valeur)*vecteur_poids[10])

maxdd=maxDrawdown(rendement_histo_port)
cat("Le MaxDrowDown: ", maxdd)


cat("Le Ratio de sharpe annualisé: ", SharpeRatio.annualized(as.xts(rendement_histo_port),Rf=0.007, scale = 255))
sigma_portefeuille = vecteur_poids%*%cov(matrice_rendement)%*%vecteur_poids
cat("Volatilité du portefeuille", sigma_portefeuille)
ratio_sharpe = mean(rendement_moy_port*255-0.007)/(sigma_portefeuille*sqrt(255))
#cat("Ratio de sharpe", ratio_sharpe)

```

## Matrice des corrélations

- Il existe une forte dépendance entre les rendements des actifs du portefeuille.

```{r cars7, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}



mat_correl = cor(matrice_rendement)

mat_covar = cov(matrice_rendement)

# Plot matrice des correlation et des rendements avec ggplot et melt de reshape2
mat_correl[lower.tri(mat_correl)]=NA
mat_covar[lower.tri(mat_covar)]=NA
melted_cor = melt(mat_correl)
melted_cov = melt(mat_covar)
ggplot(data = melted_cor, aes(x=Var1, y=Var2, fill=value)) + 
  geom_tile(color="white") +scale_fill_gradient2(low = "blue", high = "red", mid = "white", 
   midpoint = 0, limit = c(-1,1), space = "Lab")+ labs(title="Matrice des correlations", x="variable 1", y="variable 2")+theme(axis.text.x = element_text(angle = 45, vjust = 1, 
    size = 13, hjust = 1))
```

## Matrice des covariances
```{r cars8, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
ggplot(data = melted_cov, aes(x=Var1, y=Var2, fill=value)) + 
  geom_tile() + scale_fill_gradient2(low = "blue", high = "red", mid = "white", 
   midpoint = 0, limit = c(0.0001,0.0005), space = "Lab")+labs(title="Matrice des covariances", x="variable 1", y="variable 2")+theme(axis.text.x = element_text(angle = 45, vjust = 1, 
    size = 13, hjust = 1))

```

## Test de normalité des rendements

- Les rendements historiques journaliers ne suivent pas une loi normale.
- Pour but de simplification nous supposons le contraire.
```{r cars10, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
qqplot(quantile(rnorm(1561,mean((as.data.frame(rendement_histo_port))$portfolio.returns),sd=sd((as.data.frame(rendement_histo_port))$portfolio.returns)), probs = seq(0, 0.99975, 1/2000)),
       quantile((as.data.frame(rendement_histo_port))$portfolio.returns, probs = seq(0, 0.99975, 1/2000)),xlab = "quantile théorique", ylab = "quantile empirique", main="Q-rendements et Q-loi normale")

abline(a=0,b=1,col='red')
```

```{r cars11, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
shapiro.test((as.data.frame(rendement_histo_port))$portfolio.returns)

jarque.bera.test((as.data.frame(rendement_histo_port))$portfolio.returns)
```

## VaR à 99% et CTE

```{r cars12, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
varhistoriq = VaR(rendement_histo_port,p=0.99,method = "historical")

cat("VaR historique à 99% : ", varhistoriq)

vargauss = VaR(rendement_histo_port,p=0.99,method = "gaussian")

cat("VaR théorique à 99% de la loi normale: ", vargauss)
```

```{r cars13, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
t_varhistoriq=mean(rendement_histo_port$portfolio.returns[rendement_histo_port$portfolio.returns<rep(varhistoriq,1561)])

cat("CTE à partir de la VaR historique à 99%: ", t_varhistoriq)
```

## Simulation : Projection du portefeuille sur un an.

Qu'est ce que l'avenir réserve à ce portefeuille ?

- Projection du portefeuille du 01/01/2016 au 31/12/2016

- Modèle de projection utilisé : Black and Scholes. (wikipedia)


$$S_{t}=S_{0}exp{(\sigma W_{t}+\mu t)}$$
$$S_{t}:\: prix\: du\: portefeuille\: à\: t$$ 
$$S_{0} : prix\: au\: 31/12/2015$$
$$\sigma : volatilité$$ 
$$W_{t} : Mouvement\: brownien\: standard$$
$$\mu : rendement$$
$$t : temps\: t=0 \leq t \leq T=1$$
- Nous supposons que le portefeuille ne paie pas de dividende.

- Estimation des paramètres par les moments d'ordre 1 et 2.


## Simulation : Mouvement Brownien standard

```{r cars14, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
# La fonction black and scholes
# Paramètres
# le taux sans risque r=0.7%

##############################################
# MOUVEMENT BROWNIEN 
#############################################

nsimu=100
ar=0
mat_simu_brownien=matrix(0,261,nsimu)

for(ar in seq(1,nsimu,1)){
  
temps=seq(0,1,length=261)

pas=1/260

Bacc = rnorm(260,sd=sqrt(pas))

# Trajectoire
Bsim = c(0,cumsum(Bacc))

mat_simu_brownien[,ar]=Bsim
  
# plot(temps,Bsim, type="l")

}

test_data_long = mat_simu_brownien[1:115,1:3]
test_data_long <- melt(test_data_long)  # convert to long format

ggplot(data=test_data_long,
       aes(x=Var1, y=value, color=Var2)) + 
  geom_line()+labs(title="Mouvements browniens standard", x="Simulation au temp t", y="Valeur du MB")
```

## Simulation : Rendements simulées

```{r cars15, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
mat_simu_brownien_df = data.frame(mat_simu_brownien)

########### Projection du prix du portefeuille par ,St de black and scholes #######################"
# S0 = 10; 

i = 0
mat_predict=matrix(0,261,nsimu)

for (i in seq(1,nsimu,1)){
  
  mat_predict[,i]=portefeuille_historique$Valeur[1562]*exp(sd((as.data.frame(rendement_histo_port))$portfolio.returns)*mat_simu_brownien[,i]+mean((as.data.frame(rendement_histo_port))$portfolio.returns)*(1/255))

  }  

# Transormation de la matrice de prediction qui contien les prix projet pour chaque mouvement brownien simulé 

matrice_prix_predict_df = data.frame(row.names = as.Date(date2016$Date),mat_predict)


############### calcul des rendements projetés à partir de la matrice des prix  

Rendement_simu = na.omit(Return.calculate(matrice_prix_predict_df))

var99_predict = data.frame(v=colQuantiles(Rendement_simu, probs=(1-0.99)))

var99_predict_moy = mean(var99_predict$v)

############ Grphique des rendement projetés

dygraph(as.xts(Rendement_simu), main = "Evolution du rendement du PF avant Stress", ylab = "Rendement du PF")%>%
  dyLegend(show = "follow")

```

## Simulation : Evolution des VaR

```{r cars30, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
ggplot(var99_predict) +
  geom_line(aes(seq(1,nsimu,1), v)) +
  labs(title="Evolution des VaR à 99% simulées", x="simulation n", y="VaR à 99% par simulation")


varminim = min(var99_predict)

cat("VaR minimale à 99%, de toute les simulations: ", varminim)

cat("VaR moyenne à 99%, de toute les simulations: ", var99_predict_moy)

cat("VaR historique à 99% : ", varhistoriq)



```

## Simulation : Volatilités observées

```{r cars16, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
ev_volatilite_proj=data.frame(v=(apply(as.matrix(Rendement_simu),2,sd)))

ggplot(ev_volatilite_proj) +
  geom_line(aes(seq(1,nsimu,1), v)) +
  labs(title="Evolution des volatilités simulées", x="simulation n", y="Volatilité par simulation")

vola_proj_max = max(ev_volatilite_proj$v)

cat("Volatilité maximale observée , de toute les simulations (avant crise): ", vola_proj_max)

cat("Volatilité moyenne observée , de toute les simulations (avant crise): ", mean(ev_volatilite_proj$v))

```

## Simulation de crise : Scénario extrême

Pour la simulation de crise nous partons du scénario imaginatif suivant : 

- Baisse du rendement moyen de 80%
- Augmentation de la volatilité de 100%

Nous ne tenons pas compte de la probabilité de survenance de crise.

## Simulation de crise : VaR stressée

```{r cars17, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}
i = 0
mat_predict2=matrix(0,261,nsimu)

for (i in seq(1,nsimu,1)){
  
  mat_predict2[,i]=portefeuille_historique$Valeur[1562]*exp(sd((as.data.frame(rendement_histo_port)*2)$portfolio.returns)*mat_simu_brownien[,i]+mean((as.data.frame(rendement_histo_port))$portfolio.returns)*0.2*(1/255))

  }  

# Transormation de la matrice de prediction qui contien les prix projet pour chaque mouvement brownien simulé 

matrice_prix_predict_df2 = data.frame(row.names = as.Date(date2016$Date),mat_predict2)

############### calcul des rendements projetés à partir de la matrice des prix  

Rendement_simu2 = na.omit(Return.calculate(matrice_prix_predict_df2))

var99_predict2 = data.frame(v=colQuantiles(Rendement_simu2, probs=(1-0.99)))

ggplot(var99_predict2) +
  geom_line(aes(seq(1,nsimu,1), v)) +
  labs(title="Evolution des VaR à 99% simulées et stressées", x="simulation n", y="VaR à 99% par simulation puis stress")


###########" Calcul des cte par colonnes de rendement simulées ##########################
j=0
mat_var=matrix(0,260,nsimu)
for(j in seq(1,nsimu,1)){
  
  mat_var[,j] = Rendement_simu2[,j]<var99_predict2$v[j]
  
}
mat_var=mat_var*Rendement_simu2

t_var_mat_var = data.frame(v=apply(mat_var,2,mean))

min_tvar_stressed = min(t_var_mat_var)

# Plot de l'évolution des cte

# ggplot(t_var_mat_var) +
#   geom_line(aes(seq(1,100,1), v))+
#   labs(title="Evolution des CTE simulées", x="simulation n", y="CTE simulation")

#Moyenne des Tvar simulées

moy_tvar_simu = mean(t_var_mat_var$v)

# Var minimal de la simulation de VaR

varminim2 = min(var99_predict2)

var99_predict2=as.matrix(round(var99_predict2$v,digits = 2))

cat("VaR minimale à 99% de toutes simulations (stressée): ", min(var99_predict2))

cat("VaR minimale à 99%, de toute les simulations: ", varminim)

cat("VaR moyenne de toutes simulations (stressée): ", mean(var99_predict2))

cat("VaR moyenne à 99%, de toute les simulations: ", var99_predict_moy)

cat("VaR historique à 99% : ", varhistoriq)

# cat("CTE minimale de toutes simulations: ",min_tvar_stressed)
# 
# cat("CTE moyenne de toutes simulations: ",moy_tvar_simu)

```


## Volatilités en simulation de crise

```{r cars19, cache=TRUE, eval=TRUE, dpi=100, echo=FALSE, cache.comments=FALSE, warning=FALSE, fig.align="center"}

############ Graphique des volatilites projetés au lieu de supersossé les 100 colonnes de rendement simuléées ####

ev_volatilite_proj2=data.frame(v=(apply(as.matrix(Rendement_simu2),2,sd)))

ggplot(ev_volatilite_proj2) +
  geom_line(aes(seq(1,nsimu,1), v)) +
  labs(title="Evolution des volatilités en simulation de crise", x="simulation n", y="Volatilité par simulation")

# Volatilité projetée maximale

vola_proj_max2 = max(ev_volatilite_proj2$v)


cat("Volatilité maximale observée après simulation de crise: ",vola_proj_max2)

cat("Volatilité maximale observée , de toute les simulations (avant crise): ", vola_proj_max)

cat("Volatilité moyenne observée après simulation de crise: ",mean(vola_proj_max2))

cat("Volatilité moyenne observée , de toute les simulations (avant crise): ", mean(ev_volatilite_proj$v))

```

## Conclusion

1. Nous avons projeté le portefeuille sur un an puis en avons calculé la VaR.

2. Nous avons imaginé puis impacté les paramètres de projection du portefeuille et avons calculé à nouveau la VaR dite stressée.

3. On en déduit la VaR en situation sans crise et la VaR si toutefois une crise se produit sur les marchés financiers.

4. Nous observons ensuite des indicateurs de risque de marché avant et après crise : VaR et Volatilité.

5. Nous avons posé un ensemble d'hypothèse à des buts de simplification.

6. l'estimation des paramètres et les modèles de VaR peuvent être très complexes.

7. La simulation de crise pour le risque de marché permet d'estimer la perte qui pourrait survenir si un scénario extrême se produisait sur les marchés financiers.

8. Les entreprises peuvent donc estimer leur résistance et mettre en place les processus de gestion de risque adéquats.