[Returnanalytics-commits] r2236 - pkg/PerformanceAnalytics/sandbox/Meucci/R

[noreply at r-forge.r-project.org]

Author: mkshah
Date: 2012-08-16 23:32:08 +0200 (Thu, 16 Aug 2012)
New Revision: 2236

Modified:
   pkg/PerformanceAnalytics/sandbox/Meucci/R/EmpiricalMultivariateOUnCointegration.R
Log:
Remaining functions added

Modified: pkg/PerformanceAnalytics/sandbox/Meucci/R/EmpiricalMultivariateOUnCointegration.R
===================================================================
--- pkg/PerformanceAnalytics/sandbox/Meucci/R/EmpiricalMultivariateOUnCointegration.R	2012-08-16 21:13:09 UTC (rev 2235)
+++ pkg/PerformanceAnalytics/sandbox/Meucci/R/EmpiricalMultivariateOUnCointegration.R	2012-08-16 21:32:08 UTC (rev 2236)
@@ -20,8 +20,54 @@
 
   VecSig_tau = Sig_tau
   dim( VecSig_tau ) = c( N^2 , 1 )
-  VecSig = solve( diag( N^2 ) - expm(- TsT * tau ) ) %*% TsT %*% VecSig_tau
+  VecSig = solve( diag( N^2 ) - expm( Matrix( -TsT * tau ) ) ) %*% TsT %*% VecSig_tau
   Sig = VecSig
   dim( Sig ) = c( N , N )
 }
 
+OUstep = function( X_0 , t , Mu , Th , Sig )
+{
+  NumSimul = nrow( X_0 )
+  N = ncol( X_0 )
+
+  # location
+  ExpM = expm( Matrix( -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 )
+              
+  # 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
+  Sig_t = VecSig_t
+  dim( Sig_t ) = c( N , N )
+  Sig_t = ( Sig_t + t( Sig_t ) ) / 2
+
+  Eps = mvrnorm( NumSimul, rep( 0 , N ), Sig_t )
+                     
+  X_t = Mu_t + Eps
+  Mu_t = t( colMeans( Mu_t ) )
+}
+
+ProjectOU = function( x_0 , t , Mu , Th , Sig )
+{
+  N = length( x_0 )
+
+  # location
+  Mu_t = Mu + expm( Matrix( -Th %*% t ) ) %*% ( x_0 - Mu )
+
+  # 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
+  Sig_t = VecSig_t
+  dim( Sig_t ) = c( N , N )
+}
\ No newline at end of file