1、Monte Carlo SimulationFE5107 Risk Analysis and ManagementInstructor: Xianhua PengDepartment of MathematicsHKUSTmaxhpengust.hkFall 2014Xianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 1 / 59Outline1 Generating Multivariate Normal and t Random Vectors2 Simulating Stock Price Under Black-Scholes M
2、odel3 The Monte Carlo Framework and Output Analysis4 Calculating VaR via Simulation5 Importance Sampling6 Relative ErrorsXianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 2 / 59Generating Multivariate Normal and t Random VectorsOutline1 Generating Multivariate Normal and t Random Vectors2 Simula
3、ting Stock Price Under Black-Scholes Model3 The Monte Carlo Framework and Output Analysis4 Calculating VaR via Simulation5 Importance Sampling6 Relative ErrorsXianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 3 / 59Generating Multivariate Normal and t Random VectorsGenerate Multivariate Normal R
4、andom VectorsWe want to simulate X = (X1,.,Xn)T d Nn(,)Let Z = (Z1,.,Zn)T d Nn(0,In). We know how to simulate Z.Let A Rnn. We know + AZ d Nn(,AAT)If we find A such that AAT = , then we obtain + AZ d Nn(,).How to find such a A?Xianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 4 / 59Generating Mul
5、tivariate Normal and t Random VectorsSpectral Decomposition of Spectral decomposition of real-valued symmetric matrix : = T, where = diag1,.,nA can be chosen to beA = diagradicalbig1,.,radicalbignXianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 5 / 59Generating Multivariate Normal and t Random
6、VectorsCholesky Decomposition of Cholesky decomposition: any positive semi-definite symmetric matrix can be decomposed into = AAT,where A is a lower triangular matrix with all elements above diagonalbeing zero.The 3-dimensional case:11 12 1321 22 2331 32 33=a11 0 0a21 a22 0a31 a32 a33a11 a21 a310 a2
7、2 a320 0 a33Matlab function for Cholesky decomposition: cholBe careful! If C = chol(), then = CT C!Matlab function for simulating Nn(,): mvnrndXianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 6 / 59Generating Multivariate Normal and t Random VectorsGenerate Multivariate t Random VectorsThe rand
8、om vector X has a multivariate t distribution with dof (denoted by tn(,), ifX = + HradicalbigV/, whereH d Nn(0,)V d 2 and V is independent of ZX can be simulated by simulating H and V respectively.Matlab function: mvtrndXianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 7 / 59Simulating Stock Pri
9、ce Under Black-Scholes ModelOutline1 Generating Multivariate Normal and t Random Vectors2 Simulating Stock Price Under Black-Scholes Model3 The Monte Carlo Framework and Output Analysis4 Calculating VaR via Simulation5 Importance Sampling6 Relative ErrorsXianhua Peng (HKUST) Monte Carlo Simulation F
10、all 2014 8 / 59Simulating Stock Price Under Black-Scholes ModelBrownian MotionA stochastic process Wt,t 0 is a called a Brownian motion withdrift and volatility , if it has the following propertiesW0 = 0Independent increments: for any 0 0 and t 0, the distribution ofWt+s Ws does not depend on s.Norm
11、al distribution: for any s 0 and t 0,Wt+s Ws d N(t,t2)Xianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 9 / 59Simulating Stock Price Under Black-Scholes ModelBlack-Scholes Model of Stock PricesBlack-Scholes model: stock price at time t isSt = S0 exp( 22 )t + Wt), for all t,where Wt is a standard
12、 Brownian motion.Overall log return over 0,t is a Brownian motion:log StS0= ( 22 )t + Wtd N(t( 22 ),t2)Logarithm of stock price is a constant plus a Brownian motion:log St = log S0 + ( 22 )t + Wt, for all tLog return of the stock from time s to time s + t:log Ss+tSs= ( 22 )t + (Ws+t Ws)d N( 22 )t,2t
13、)Xianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 10 / 59Simulating Stock Price Under Black-Scholes ModelSimulate Stock Prices Under Black-Scholes ModelGiven S0, how to simulate the future stock price at timet1 0) = E1(0,)(L) = E1(0,)(V(y) V(y + X) 99% VaR of the loss: VaR99%(L) = F1L (99%)Xian
14、hua Peng (HKUST) Monte Carlo Simulation Fall 2014 14 / 59The Monte Carlo Framework and Output AnalysisThe Monte Carlo ApproachThe Monte Carlo approach for estimating = Eh(X):1 Simulate i.i.d. X 1,X 2,.,X m d X2 Calculate Hi := h(X i), i = 1,2,.,m.3 Estimate bym := 1mmsummationdisplayi=1HiXianhua Pen
15、g (HKUST) Monte Carlo Simulation Fall 2014 15 / 59The Monte Carlo Framework and Output AnalysisThe Monte Carlo Approach (continued)Is m a good estimator for ?m is unbiased, i.e.,Em = E 1mmsummationdisplayi=1Hi = 1mmsummationdisplayi=1EHi = .m is consistent, i.e., by strong law of large numbers,m = 1
16、mmsummationdisplayi=1Hi EH1 = , as m .Xianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 16 / 59The Monte Carlo Framework and Output AnalysisOutput Analysis of Monte Carlo ApproachHow far away is m different from the true value ? The answer lies incentral limit theorem:TheoremSuppose H1,H2,.,Hm a
17、re i.i.d. random variables with mean andvariance 2, thenm /md N(0,1), as m ,where m := 1m summationtextmi=1 Hi is the sample average.More precisely,PparenleftBiggm /m xparenrightBigg (x), as m , for any x R,where (x) is the CDF of N(0,1) distribution.Xianhua Peng (HKUST) Monte Carlo Simulation Fall
18、2014 17 / 59The Monte Carlo Framework and Output AnalysisOutput Analysis of Monte Carlo Approach (continued)Let z = 1(), i.e., the quantile of N(0,1) distribution, e.g.,z0.025 = 1.96 and z0.975 = 1.96.When m is large,P(z2 0, then define =p(1p)fL(xp) .Then L(mp) xp/md N(0,1), as m .More precisely,Ppa
19、renleftbiggL(mp) xp/m yparenrightbigg (y), as m , for any y R,where (y) is the CDF of N(0,1) distribution.Xianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 28 / 59Calculating VaR via SimulationCentral Limit Theorem for F1L (p) (continued)By the CLT for sample quantile, when m is large, the distr
20、ibution ofL(mp) is approximately normal with mean xp and variance 1m p(1p)fL(xp)2,i.e.,L(mp) d Nparenleftbiggxp, 2mparenrightbigg, where 2 = p(1 p)fL(xp)2.Xianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 29 / 59Calculating VaR via SimulationCalculating VaR via Simulation (continued)How to const
21、ruct a confidence interval for VaRp = F1L (p) = xp?By Bahadur representation, the 1 “confidence interval“ for xpis given bybracketleftbiggL(mp) z12 m,L(mp) z2 mbracketrightbigg, =radicalbigp(1 p)fL(xp) (1)If we can find a good estimate fL(xp) for fL(xp), then we canconstruct a 1 confidence interval for xp by replacing fL(xp) byfL(xp) in (1).However, estimating probability density is a notoriously difficultstatistical problem.Alternative method to construct confidence interval: SectioningXianhua Peng (HKUST) Monte Carlo Simulation Fall 2014 30 / 59