SIM_WASC_2D Simulates a 2d WASC model for the given time points in t. [y, vol] = sim_wasc_2d(y0, mu, Sigma_0, M, Q, rho, beta, t) Performs an Euler simulation of the WASC model specified by the given parameters for the time points in t. Assuming that n is the dimension of the problem, the parameters expected are INPUT y_0: n-dimensional vector with initial log prices mu: n-dimensional mean vector mu Sigma_0: n x n covariance matrix M: n x n matrix Q: orthogonal n x n matrix rho: n-dimensional correlation vector beta: a natural number t: T-dimensional vector of the time points for which a price is to be simulated OUTPUT y: T x n matrix of the simulated log prices vol: T x 3 matrix of the simulated variances and correlation See also SIM_PCSV. created by Benedikt Rudolph DATE: 12-Aug-2012
0001 function [y, vol] = sim_wasc_2d(y0, mu, Sigma_0, M, Q, rho, beta, t) 0002 %SIM_WASC_2D Simulates a 2d WASC model for the given time points in t. 0003 % 0004 % [y, vol] = sim_wasc_2d(y0, mu, Sigma_0, M, Q, rho, beta, t) 0005 % Performs an Euler simulation of the WASC model specified by the given 0006 % parameters for the time points in t. 0007 % Assuming that n is the dimension of the problem, the parameters expected are 0008 % 0009 % INPUT y_0: n-dimensional vector with initial log prices 0010 % mu: n-dimensional mean vector mu 0011 % Sigma_0: n x n covariance matrix 0012 % M: n x n matrix 0013 % Q: orthogonal n x n matrix 0014 % rho: n-dimensional correlation vector 0015 % beta: a natural number 0016 % t: T-dimensional vector of the time points for which 0017 % a price is to be simulated 0018 % 0019 % OUTPUT y: T x n matrix of the simulated log prices 0020 % vol: T x 3 matrix of the simulated variances and correlation 0021 % 0022 % See also SIM_PCSV. 0023 % 0024 % created by Benedikt Rudolph 0025 % DATE: 12-Aug-2012 0026 0027 mu = reshape(mu, 1, 2); 0028 T = length(t); % number of time grid points 0029 dt = reshape(diff(t),T-1,1); % time increments 0030 0031 y = zeros(T, 2); % pre-allocation 0032 y(1,:) = y0; 0033 0034 vol = zeros(T, 3); 0035 vol(1,:) = [ sqrt(Sigma_0(1,1)), sqrt(Sigma_0(2,2)), ... 0036 Sigma_0(1,2)/sqrt(Sigma_0(1,1)*Sigma_0(2,2)) ]; 0037 0038 % choose values for X_0 using eigenvalues and 0039 % eigenvectors of Sigma_0 0040 [V, Lambda] = eig(Sigma_0); 0041 X = repmat(sqrt(diag(Lambda)),1,2).*V'; 0042 X = [X; zeros(beta-2, 2)]; % add zeros in the remaining part 0043 0044 % init matrices which are updated in every step 0045 Sigma = Sigma_0; 0046 sqrtm_Sigma = sqrtm(Sigma); 0047 sigma_tilde = inv(sqrtm_Sigma); 0048 0049 for k=2:T 0050 0051 % simulate noise 0052 dN = sqrt(dt(k-1))*normrnd(0,1,beta,2); 0053 dB = sqrt(dt(k-1))*normrnd(0,1,2,1); 0054 dW = sigma_tilde*X'*dN; 0055 0056 % Euler approximation of y 0057 y(k,:) = y(k-1,:) + (mu - 0.5*diag(Sigma)')*dt(k-1) ... 0058 + (sqrtm_Sigma*(dW*rho + sqrt(1-rho'*rho)*dB))'; 0059 0060 X = X + X*M*dt(k-1) + dN*Q; % Euler approximation of X 0061 % updates of Sigma 0062 Sigma = X'*X; 0063 sqrtm_Sigma = sqrtm(Sigma); 0064 sigma_tilde = inv(sqrtm_Sigma); 0065 0066 % store volatilites and correlation 0067 v1 = sqrt(Sigma(1,1)); 0068 v2 = sqrt(Sigma(2,2)); 0069 vol(k,:) = [ v1, v2, Sigma(1,2)/(v1*v2) ]; 0070 end 0071 end