[Returnanalytics-commits] r2243 - in pkg/PerformanceAnalytics/sandbox/Meucci: R demo
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun Aug 19 03:54:20 CEST 2012
Author: mkshah
Date: 2012-08-19 03:54:17 +0200 (Sun, 19 Aug 2012)
New Revision: 2243
Added:
pkg/PerformanceAnalytics/sandbox/Meucci/demo/S_CheckDiagonalization.R
pkg/PerformanceAnalytics/sandbox/Meucci/demo/S_CovarianceEvolution.R
Modified:
pkg/PerformanceAnalytics/sandbox/Meucci/R/TheoryMultivariateOUnCointegration.R
Log:
Correcting functions and adding demos for TheoryMultivariateOUnCointegration.R
Modified: pkg/PerformanceAnalytics/sandbox/Meucci/R/TheoryMultivariateOUnCointegration.R
===================================================================
--- pkg/PerformanceAnalytics/sandbox/Meucci/R/TheoryMultivariateOUnCointegration.R 2012-08-19 00:36:14 UTC (rev 2242)
+++ pkg/PerformanceAnalytics/sandbox/Meucci/R/TheoryMultivariateOUnCointegration.R 2012-08-19 01:54:17 UTC (rev 2243)
@@ -4,20 +4,20 @@
N = ncol( X_0 )
# location
- ExpM = expm( Matrix( -Th %*% t ) )
+ ExpM = expm( -Th * t )
# repmat = function(X,m,n) - R equivalent of repmat (matlab)
X = t( Mu - ExpM %*% Mu )
mx = dim( X )[1]
nx = dim( X )[2]
- Mu_t = matrix( t ( matrix( X , mx , nx*1 ) ), mx * NumSimul, nx * 1, byrow = T ) + t( X_0 %*% ExpM )
+ Mu_t = matrix( t ( matrix( X , mx , nx*1 ) ), mx * NumSimul, nx * 1, byrow = T ) + X_0 %*% ExpM
# scatter
TsT = kronecker( Th , diag( N ) ) + kronecker( diag( N ) , Th )
VecSig = Sig
dim( VecSig ) = c( N^2 , 1 )
- VecSig_t = solve( TsT ) %*% ( diag( N^2 ) - expm( Matrix( -TsT %*% t ) ) ) %*% VecSig
+ VecSig_t = solve( TsT ) %*% ( diag( N^2 ) - expm( -TsT * t ) ) %*% VecSig
Sig_t = VecSig_t
dim( Sig_t ) = c( N , N )
Sig_t = ( Sig_t + t( Sig_t ) ) / 2
@@ -36,13 +36,13 @@
N = ncol( X_0 )
# location
- ExpM = expm( Matrix( -Th %*% Dt ) )
+ ExpM = expm( as.matrix( -Th %*% Dt ) )
# repmat = function(X,m,n) - R equivalent of repmat (matlab)
X = t( Mu - ExpM %*% Mu )
mx = dim( X )[1]
nx = dim( X )[2]
- Mu_t = matrix( t ( matrix( X , mx , nx*1 ) ), mx * NumSimul, nx * 1, byrow = T ) + t( X_0 %*% ExpM )
+ Mu_t = matrix( t ( matrix( X , mx , nx*1 ) ), mx * NumSimul, nx * 1, byrow = T ) + X_0 %*% ExpM
# scatter
Sig_t = Sig %*% Dt
Added: pkg/PerformanceAnalytics/sandbox/Meucci/demo/S_CheckDiagonalization.R
===================================================================
--- pkg/PerformanceAnalytics/sandbox/Meucci/demo/S_CheckDiagonalization.R (rev 0)
+++ pkg/PerformanceAnalytics/sandbox/Meucci/demo/S_CheckDiagonalization.R 2012-08-19 01:54:17 UTC (rev 2243)
@@ -0,0 +1,27 @@
+# This script verifies the correctness of the eigenvalue-eigenvector representation in terms of real matrices
+# for the transition matrix of an OU process discussed in A. Meucci (2009)
+# "Review of Statistical Arbitrage, Cointegration, and Multivariate Ornstein-Uhlenbeck"
+# available at ssrn.com
+
+# Code by A. Meucci, April 2009
+# Most recent version available at www.symmys.com > Teaching > MATLAB
+
+N = 5
+Theta = matrix( runif( N^2 ), 5, byrow = T )
+
+tmp = eigen( Theta )
+B = tmp$vectors
+L = tmp$values
+
+A = Re( B ) - Im( B )
+
+Index_j = as.matrix( which( Im ( L ) != 0 ) )
+L = diag( L )
+
+G = L
+for ( s in c( seq( from = 1, to = length( Index_j ), by = 2 ) ) )
+{
+ G[ Index_j[ s : (s+1) ] , Index_j[ s : (s+1) ] ] = rbind( c( 1 , 0 ), c( 0 , 1 ) ) * Re( L[ Index_j[s],Index_j[s] ] ) +
+ rbind( c( 0 , 1 ), c( -1 , 0 ) ) * Im( L [ Index_j[s] , Index_j[s] ] )
+}
+Theta_ = A %*% G %*% solve( A )
Added: pkg/PerformanceAnalytics/sandbox/Meucci/demo/S_CovarianceEvolution.R
===================================================================
--- pkg/PerformanceAnalytics/sandbox/Meucci/demo/S_CovarianceEvolution.R (rev 0)
+++ pkg/PerformanceAnalytics/sandbox/Meucci/demo/S_CovarianceEvolution.R 2012-08-19 01:54:17 UTC (rev 2243)
@@ -0,0 +1,75 @@
+# This script represents the evolution of the covariance of an OU process in terms of the dispersion ellipsoid
+# see A. Meucci (2009)
+# "Review of Statistical Arbitrage, Cointegration, and Multivariate Ornstein-Uhlenbeck"
+# available at ssrn.com
+
+# Code by A. Meucci, April 2009
+# Most recent version available at www.symmys.com > Teaching > MATLAB
+
+# input parameters of multivariate OU process
+K = 1
+J = 1
+
+x0 = matrix( runif(K + 2*J), ncol = 1 )
+
+Mu = matrix( runif(K + 2*J), ncol = 1 )
+
+A = matrix( runif( ( K + 2*J )^2 ), nrow = K + 2*J, byrow = T ) - .5
+ls = matrix( runif( K ), ncol = 1 ) - .5
+gs = matrix( runif( J ), ncol = 1 ) - .5
+os = matrix( runif( J ), ncol = 1 ) - .5
+
+S = matrix( runif( ( K + 2*J )^2 ), nrow = K + 2*J, byrow = T ) - .5
+
+ts_0 = .01 * as.matrix( seq( from = 0, to = 100, by = 10 ) )
+NumSimul = 10000
+
+# process inputs
+Gamma = diag( ls )
+for ( j in 1 : J )
+{
+ G = rbind( cbind( gs[j], os[j] ), cbind( -os[j], gs[j] ) )
+ Gamma = as.matrix( bdiag( Gamma , G ) )
+}
+
+Theta = A %*% Gamma %*% solve( A )
+
+# one-step exact simulation
+Sigma = S %*% t( S )
+X = t( x0 )
+mx = dim( X )[ 1 ]
+nx = dim( X )[ 2 ]
+X_0 = matrix( t ( matrix( X , mx , nx ) ), mx * NumSimul , nx , byrow=T )
+OUstepResult = OUstep( X_0 , ts_0[ nrow( ts_0 ) ] , Mu , Theta , Sigma )
+
+# multi-step simulation: exact and Euler approximation
+X_t = X_0
+X_tE = X_t
+for ( s in 1:length( ts_0 ) )
+{
+ Dt = ts_0[ 1 ]
+ if ( s > 1 )
+ {
+ Dt = ts_0[s] - ts_0[ s - 1 ]
+ }
+ PostOUStepResult = OUstep( X_t , Dt , Mu , Theta , Sigma )
+ #[X_tE,MuHat_tE,SigmaHat_tE]=OUstepEuler(X_tE,Dt,Mu,Theta,Sigma);
+}
+
+# plots
+Pick = cbind( K + 2*J - 1, K + 2*J )
+
+# horizon simulations
+plot( OUstepResult$X_t[ , Pick[ 1 ] ] , OUstepResult$X_t[ , Pick[ 2 ] ] )
+
+# horizon location
+plot( OUstepResult$Mu_t[ Pick[ 1 ] ] , OUstepResult$Mu_t[ Pick[ 2 ] ] )
+
+# horizon dispersion ellipsoid
+# TwoDimEllipsoid(MuHat_t1(Pick),SigmaHat_t1(Pick,Pick),2,0,0);
+
+# starting point
+plot( x0[ Pick[ 1 ] ] , x0[ Pick[ 2 ] ] )
+
+# starting generating dispersion ellipsoid
+# TwoDimEllipsoid(x0(Pick),Sigma(Pick,Pick),2,0,0);
\ No newline at end of file
More information about the Returnanalytics-commits
mailing list