/* ** Version 01 ** Date : 20/07/2001 ** -------------------------------------------------------------------------- ** Copyright K.Cuthbertson and D. Nitzsche ** "Financial Engineering:Derivatives and Risk Manangement" - J. Wiley 2001 ** ** VaR OF A TWO ASSET PORTFOLIO USING MCS ** ** Portfolio = long calls (on S1) and long stocks of S2 ** This is a similar (but slightly easier) problem to that found on p676-681 ** ** We use the procedure "rndnm" (see below) to generate multivariate normal ** series using Choleski. We use sqrt(T)-rule to obtain VaR over a T-day horizon. ** ** The program is written for the case of 2x2 underlying assets for pedagogic ** reasons, although it is not difficult to generalise it to the nxn case, using ** more matrix commands. ** */ new ; cls ; format /m1/rd 12,2 ; outwidth 225 ; output on ; /* output file = c:\docs\arun\out.out reset ; */ @ USER: Check path for output file @ ; @ ---------------------------------- USER INPUT ( n=2 asset case only) -------------------------------------- @ runs = 5000 ; @ Runs in MCS,this should be > 10,000 @ n = 2 ; @ Number of underlying assets (eg. stocks S1 and S2) @ mu = zeros(n,1) ; @ Annual mean growth of underlying assets, propn @ sigma = zeros(n,1) ; @ Annual stdv, propn, for Black-Scholes @ sigmaD = zeros(n,1) ; @ Daily stdv, propn, for VaR, MCS @ vcvD = zeros(n,n) ; @ Daily variance-covariance matrix @ sigma = { 0.30,0.30 } ; sigma[1,1] = 0.3 ; @ std deviation underlying asset-1, ANNUAL, propn @ sigma[2,1] = 0.3 ; @ std deviation underlying asset-2, ANNUAL, propn @ rho12 = 0.7 ; @ correlation coefficient( dimensionless) between returns on the 2 underlying assets S1,S2, @ mu[1,1] = 0.10 ; @ growth rate asset-1,annual, propn @ mu[2,1] = 0.10 ; @ growth rate asset-2,annual, propn @ T = 25 ; @ VaR horizon-DAYS @ NV0 = { 100, 10 } ; @ Initial holding of 100 long calls (on S1) plus 10 stocks of S2: 2x1 COLN vector @ /* ** Changing NV0 so you WRITE 100 calls and hold the stocks gives a v. different VaR - try it! ** NV0 = { -100, 10 } ; */ s0 = { 10, 10 } ; @ Initial $-price of underlying S1 and S2: 2x1 COLN vector @ @ --------------------------- Additional inputs for the option on underlying asset-1 ------------------------ @ K = 10 ; @ strike price: at-the-money @ tau = 0.5 ; @ time to maturity, years @ r = 0.05 ; @ risk free rate @ @ ---------------------------------------------- END OF USER INPUT ------------------------------------------- @ sT = zeros(runs,rows(mu)) ; @ vector for final values at T of the underlying assets @ s0 = ones(runs,rows(mu)) .* s0'; @ initial value of the underlying assets at t=0: s01,s02.. etc in colns @ @ ------------------------------- Calculate DAILY, mean and variance-covariance (vcvD) matrix ------------------------ @ mu = mu/365 ; @ daily growth, propn @ sigmaD[1,1] = sigma[1,1]./sqrt(252) ; @ daily std dev., propn @ sigmaD[2,1] = sigma[2,1]./sqrt(252) ; @ daily std dev., propn @ vcvD[1,2] = rho12*sigmaD[1,1]*sigmaD[2,1] ; @ formula for daily covariance, propn @ vcvD[2,1] = vcvD[1,2] ; @ symmetric 'SIGMA' for VCV matrix @ vcvD[1,1] = sigmaD[1,1]^2 ; @ Now make sigmaD a VARIANCE, propn @ vcvD[2,2] = sigmaD[2,1]^2 ; @ Now make sigmaD a VARIANCE, propn @ @ --------------------------- End of calculations of DAILY, variance-covariance (VCV) matrix ------------------------- @ /* ------------------------------------------------------------------------------------------------------------------- ** Procedure for the Choleski decomposition of the VCV matrix (Note: NOT Choleski of CORRELATION matrix) ** Output y is 2 series where y1=N(mu1,vcv11) and y2=N(mu2, vcvD22) and y1 and y2 have a correlation ** coefficient equal to rho12 ** Note: y1 and y1 do NOT have variance=1, as in the text, but already have variances, vcvDD11 and vcvD22 -------------------------------------------------------------------------------------------------------------------- */ {y} = rndmn(mu,vcvD,runs); @ procedure to give y= N(mu,vcvD), ( y is 'runs x 2') @ /* ------------------------------------------------------------------------------------------------------------------- Calculate the final values (at T) the 2 underlying assets : This uses root-T rule, saves a do-loop Note: We do not require 'sigmaD' in the formula below because y1 already incorporates the variance -------------------------------------------------------------------------------------------------------------------- */ j=1; do until j > rows(mu) ; sT[.,j] = s0[.,j].* exp( sqrt(T)*y[.,j] ) ; j = j+1 ; endo ; @ --------------------------------------- Black-Scholes uses ANNUAL standard deviations -------------------------------- @ {callp0} = bscall( K,s0[.,1],sigma[1,1],tau,r ) ; @ Uses ANNUAL sigma @ {callpT} = bscall( K,sT[.,1],sigma[1,1],tau-T/365,r ) ; dcall = callpT - callp0 ; ds = sT - s0 ; dP = dcall~ds[.,2] ; @ Chg price of call (on S1) and in value of S2: (runs x 2) @ Val0 = NV0'*( callp0[1,1] | s0[1,2] ) ; @ Initial $-Value, portfolio, calls and stocks S2, (1x2)x(2x1)=scalar @ dVal1_2 = NV0'.*dP ; @ Changes in $-value of all calls then all stocks S2: (runsx2) @ dVal = sumc( (dVal1_2)' ) ; @ Sums the ROWS of dVal1_2 and puts result into a coln: (runsx1) @ /* ---- Check the data if you wish ------ ** ?; " call premiums at 0 = " callp0' ; ** ?; " call premiums at T = " callpT' ; ------------------------------------------------ */ ?; " ----------------------------------------------------------------------------------------------------------------- " ; " Output: VaR of portfolio of long calls (on S1) and long position in stock S2 " ; " ----------------------------------------------------------------------------------------------------------------- " ; ?; " VaR horizon-days " T ; ?; " K, tau,r for option on asset-1 " K~tau~r ; ?; " No. of calls(on S1) and no. of stocks (S2) " NV0' ; ?; " Initial $-price of S1 and S2 " s0[1,1]~s0[1,2] ; ?; " Initial $-call premium " callp0[1,1] ; ?; " Initial $-value in all calls and S2 " NV0[1]*callp0[1,1]~NV0[2]*s0[1,2] ; ?; " Initial $-value of portfolio of long call(on S1)+stocks S2 " Val0 ; ?; " This is the $-VaR (over 25 days) of portfolio at 1% and 5% " quantile(dVal, 0.01|0.05) ; ?; " VaR at 1% and 5% as % of initial value of portfolio " 100*quantile(dVal, 0.01|0.05 )./Val0 ; ?; " ----------------------------------------------------------------------------------------------------------------- " ; " Check that the MCS gives similar variances and correlation coeff. as the inputs used " ; " ----------------------------------------------------------------------------------------------------------------- " ; ?; " Correlation coefficient as input data (rho12) " rho12 ; ?; " Corrln matrix of 2 stocks from: y= N(mu,sigma) " corrx(y) ; ?; " VCV matrix of 2 stocks (scaled to annual, propn) from MCS, y= N( mu,sigma) " 252.*vcx(y) ; ?; " VCV matrix of 2 stocks (scaled to annual propn) from input data (sigmaD) " 252.*vcvD ; ?; " Stdv of 2 stocks (scaled to annual, %) from MCS, y= N( mu,sigma) " 100.*sqrt( diag(252.*vcx(y)) ) ; ?; " Stdv of 2 stocks (scaled to annual, %) from input data (sigmaD) " 100.*sqrt( diag( 252.*vcvD) ) ; /* ----- Check more data if you wish ------- ** ?; " Chg in call premiums at T = " dcall' ; ** ?; " Chg in stock prices at T = " ds ; ** ?; " Chg in vals of elements of Portfolio (call,S2) " dP ; ** ?; " Chg in Vals of whole Portfolio (calls+S2) " dVal ; */ @ ------------- graphs ----------------- @ library gauss pgraph ; graphset ; screen off ; begwind; window(2,2,0); axmargin(1,1,1,1); xlabel("Change in value of a portfolio of calls on stock1 plus long stock2 "); ylabel("Frequency"); histp(dVal[.,1],20) ; nextwind; xlabel("Change in value of a single long call"); histp(dcall[.,1],20) ; nextwind; xlabel("Change in value of a single Stock-2"); histp(ds[.,2],20) ; endwind; @ -------------------------------------- END OF PROGRAM ----------------------------------------- @ @ ----------------------------------- PROCEDURES FOLLOW -------------------------------------- @ /* ** Procedure: Random generation of multivariate normal variables (See Appendix 24.2) ** ** INPUTS ** mu = mean, (mx1)-vector ** sigma = variance-cov matrix , (m x m) Note: input is VCV-martrix, NOT stdv or correlation matrix ** k = length of data series to generate ** ** OUTPUT ** y =(kxm) m-series of multivariate random numbers each of 'length' k ** y1~y2 have variances sigma11~sigma22 and correlation between them of rho */ proc rndmn(mu,sigma,k); local d,p,y,i; if not ( sigma == sigma'); errorlog "ERROR : Matrix sigma must be symmetric "; retp(error(0)); endif; d=rows(mu); p=chol(sigma); @ -- GAUSS code for Cholesky decompostion: p is upper triangular(mxm) -- @ y=zeros(k,d); i=1; do until i>k; y[i,.] = ( p'*rndn(d,1) )'; @ -- p' is therefore lower triangular -- @ @ -- the last transpose =' in above makes y=(kxm) -- @ i=i+1; endo; y=y+mu'; @ --- note the addition of the mean --- @ retp(y); endp; @ ------------------------------------- Black-Scholes Procedures Follow -------------------------- @ proc bscall(k,s,sigma,tau,r); local d1,d2,c ; d1=( ln(s./k) + ( r + sigma^2/2).*tau ) ./ ( sigma.*sqrt(tau) ) ; d2=d1 - sigma.*sqrt(tau); c=s.*cdfn(d1)-k.*exp(-r.*tau).*cdfn(d2); retp(c); endp; @ -------------------------------------------- Next procedure ------------------------------------ @ proc bsput(k,s,sigma,tau,r); local d1,d2,p ; d1=( ln(s./k) + ( r + sigma^2/2).*tau ) ./ ( sigma.*sqrt(tau) ) ; d2=d1 - sigma.*sqrt(tau); p= k.*exp(-r.*tau).*cdfn(-d2) - s.*cdfn(-d1); retp(p); endp; end ;