Thierry Moudiki — Nov 7, 2013, 10:03 PM
## Chargement du package
library(ESG)
## Courbe des prix zéro-coupons
data(ZC)
## 1 - Actions et taux court
simulStock <- rStock(horizon=9, nScenarios=7, ZC=ZC, vol=.1, k=2, volStock=.2, stock0=100, rho=.5)
par(mfrow=c(2,1))
matplot(t(simulStock$stockPaths), type='l')
matplot(t(simulStock$shortRatePaths), type='l')
## 2 - Taux courts, immo, défaut, liquidité
rt <- rShortRate(horizon=15, nScenarios=500, ZC=ZC, vol=.1, k=2)
re <- rRealEstate(horizon=25, nScenarios=1000, ZC=ZC, vol=.1, k=2, volRealEstate=.15, realEstate0=50)
ds <- rDefaultSpread(horizon=30, nScenarios=200, defaultSpread0=.01, volDefault=.2, alpha=.1, beta=1)
ls <- rLiquiditySpread(horizon=20, nScenarios=150, eta=.05, liquiditySpread0=.01)
x11()
par(mfrow=c(2,2))
matplot(t(rt), type='l', xlab = "horizon", ylab = "taux court")
matplot(t(re$realEstatePaths), type='l', xlab = "horizon", ylab = "immobilier")
matplot(t(ds), type='l', xlab = "horizon", ylab = "spread de défaut")
matplot(t(ls), type='l', xlab ="horizon", ylab = "spread de liquidité")
## 3 - Tous les facteurs de risque, sur un même horizon
simulAllRiskFactors <- rAllRisksFactors(horizon=5, nScenarios=10, ZC, vol=.1, k=2, volStock=.2, stock0=100, rho=.5, volRealEstate=.15, realEstate0=50, eta=.05, liquiditySpread0=.01, defaultSpread0=.01, volDefault=.2, alpha=.1, beta=1)
par(mfrow=c(2,2))
matplot(t(simulAllRiskFactors$shortRate), type='l')
matplot(t(simulAllRiskFactors$s), type='l')
matplot(t(simulAllRiskFactors$realEstate), type='l')
matplot(t(simulAllRiskFactors$defaultSpread), type='l')
## 4 - Calcul de prix, temps d'exécution (dépend du processeur)
ptm1 <- proc.time()
put_stock_distr <- rAssetDistribution(type="EuroPut_Stock",5,25,Strike=98.5,vol=.1,k=2,ZC=ZC,volStock=.2, stock0=100, rho=.5,nScenarios=10000)
bond_distr <- rAssetDistribution(type="Bond",t=3,T=35,nCoupons=20, couponsRate=0.3,vol=.1, k=2, ZC=ZC, nScenarios=10000)
zcbond_distr <- rAssetDistribution(type="Zero-Coupon",t=2,T=30,vol=.1, k=2, ZC=ZC, nScenarios=10000)
call_zc_distr <- rAssetDistribution(type="EuroCall_ZC",4,4.5,s=5, Strike=.985,vol=.1, k=2, ZC=ZC,nScenarios=10000)
ptm2 <- proc.time()
duree <- ptm2 - ptm1
duree
user system elapsed
1.21 0.00 1.20
## 5 - Calcul de BE
source("2-fonctionsESG.R")
k <- 0.12 #Vitesse de retour à la moyenne du taux court
sTaux <- 0.05 #Volatilité du processus de taux court
sUC <- .16 #Volatilité de l'UC
H <- 10 #Horizon de projection
nSimulations <- 10000 #Nombre de simulations
tauxRachatS <- .02
tauxRachatC <- .05
# Trajectoire UC et taux. Fonction rStock du package ESG
traj <- rStock(horizon=H, nScenarios=nSimulations, ZC=ZC, vol=sTaux, k=k, volStock=sUC, stock0=1,
rho=.5)
# Simulation de taux courts
trajectoiresTaux <- traj$shortRatePaths
# Simulation du cours de actions
trajectoiresUC <- traj$stockPaths
# Flux futurs
Flux_futurs <- calculFlux(trajectoiresTaux,trajectoiresUC,tauxRachatS,tauxRachatC)
# Taux d'actualisation
ActuFlux_futurs <- Flux_futurs$flux*Flux_futurs$actu
# Calcul du Best Estimate
BEempirique <- sum(ActuFlux_futurs)/nSimulations
BEempirique
[1] 1.001
## 6 - Convergence du BE
# Vecteur des moyennes empiriques
moyemp <- rep(0, nSimulations)
# Calcul des moyennes empiriques successives
w <- rep(1, nSimulations)/(1:nSimulations)
temp_moyemp <- apply(ActuFlux_futurs,2,cumsum)*w
moyemp <- apply(temp_moyemp, 1, sum)
# Graphique de convergence
x11()
plot(x=(1:nSimulations), y=moyemp,type="l",xlab="Nombre de simulations",ylab="Best Estimate")
lines(moyemp,col="blue", lwd=2)
titre = paste("Convergence du Best Estimate", "\n", "en fonction du nombre de simulations")
title(titre)
lines(x=(1:nSimulations), y=rep(1, nSimulations), col='red', lwd=2, lty='dashed')
print(paste("Valeur de la moyenne des flux futurs actualisés par simulation : ", BEempirique, sep=""))
[1] "Valeur de la moyenne des flux futurs actualisés par simulation : 1.00066180689886"