[Returnanalytics-commits] r3984 - in pkg/Meucci: R demo man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu Aug 20 17:52:02 CEST 2015
Author: xavierv
Date: 2015-08-20 17:52:02 +0200 (Thu, 20 Aug 2015)
New Revision: 3984
Modified:
pkg/Meucci/R/InvariantProjection.R
pkg/Meucci/R/LognormalCopulaPdf.R
pkg/Meucci/R/NormalCopulaPdf.R
pkg/Meucci/R/PerformIidAnalysis.R
pkg/Meucci/R/StudentTCopulaPdf.R
pkg/Meucci/R/TwoDimEllipsoid.R
pkg/Meucci/demo/AnalyticalvsNumerical.R
pkg/Meucci/demo/ButterflyTrading.R
pkg/Meucci/demo/DetectOutliersviaMVE.R
pkg/Meucci/demo/FullFlexProbs.R
pkg/Meucci/demo/FullyFlexibleBayesNets.R
pkg/Meucci/demo/FullyIntegratedLiquidityAndMarketRisk.R
pkg/Meucci/demo/HermiteGrid_CVaR_Recursion.R
pkg/Meucci/demo/HermiteGrid_CaseStudy.R
pkg/Meucci/demo/HermiteGrid_demo.R
pkg/Meucci/demo/InvariantProjection.R
pkg/Meucci/demo/MeanDiversificationFrontier.R
pkg/Meucci/demo/Prior2Posterior.R
pkg/Meucci/demo/RobustBayesianAllocation.R
pkg/Meucci/demo/RobustBayesianCaseStudy.R
pkg/Meucci/demo/S_AutocorrelatedProcess.R
pkg/Meucci/demo/S_BondProjectionPricingNormal.R
pkg/Meucci/demo/S_BondProjectionPricingStudentT.R
pkg/Meucci/demo/S_BuyNHold.R
pkg/Meucci/demo/S_CallsProjectionPricing.R
pkg/Meucci/demo/S_CheckDiagonalization.R
pkg/Meucci/demo/S_CornishFisher.R
pkg/Meucci/demo/S_CorrelationPriorUniform.R
pkg/Meucci/demo/S_CovarianceEvolution.R
pkg/Meucci/demo/S_DerivativesInvariants.R
pkg/Meucci/demo/S_DeterministicEvolution.R
pkg/Meucci/demo/S_DisplayLognormalCopulaPdf.R
pkg/Meucci/demo/S_DisplayNormalCopulaCdf.R
pkg/Meucci/demo/S_DisplayNormalCopulaPdf.R
pkg/Meucci/demo/S_DisplayStudentTCopulaPdf.R
pkg/Meucci/demo/S_DynamicManagementCase1.R
pkg/Meucci/demo/S_DynamicManagementCase2.R
pkg/Meucci/demo/S_EigenvalueDispersion.R
pkg/Meucci/demo/logToArithmeticCovariance.R
pkg/Meucci/man/Central2Raw.Rd
pkg/Meucci/man/Cumul2Raw.Rd
pkg/Meucci/man/GenerateLogNormalDistribution.Rd
pkg/Meucci/man/LognormalCopulaPdf.Rd
pkg/Meucci/man/NormalCopulaPdf.Rd
pkg/Meucci/man/PerformIidAnalysis.Rd
pkg/Meucci/man/Raw2Central.Rd
pkg/Meucci/man/Raw2Cumul.Rd
pkg/Meucci/man/StudentTCopulaPdf.Rd
pkg/Meucci/man/SummStats.Rd
pkg/Meucci/man/TwoDimEllipsoid.Rd
pkg/Meucci/man/std.Rd
Log:
Fixed documentation and reformated demo scripts and its relating functions up to S_EigenvalueDispersion included
Modified: pkg/Meucci/R/InvariantProjection.R
===================================================================
--- pkg/Meucci/R/InvariantProjection.R 2015-08-20 09:49:54 UTC (rev 3983)
+++ pkg/Meucci/R/InvariantProjection.R 2015-08-20 15:52:02 UTC (rev 3984)
@@ -7,168 +7,176 @@
#' Note the first central moment defined as expectation.
#'
#' \deqn{\tilde{ \mu } ^ {\big(n\big)} _{X} \equiv E \big\{ X^{n} \big\},
-#' \\ \mu ^{ \big(n\big) }_{X} \equiv \sum_0^{n-1} \big(-1\big)^{n-k} \mu ^{n-k}_{X} \tilde{ \mu }^{k}_{X} + \tilde{ \mu }_{X}^{n} }
+#' \\ \mu ^{ \big(n\big) }_{X} \equiv \sum_0^{n-1} \big(-1\big)^{n-k}
+#' \mu ^{n-k}_{X} \tilde{ \mu }^{k}_{X} + \tilde{ \mu }_{X}^{n} }
#'
-#' @param mu_ : [vector] (length N corresponding to order N) corresponding raw moments
+#' @param mu_ : [vector] (length N corresponding to order N) corresponding
+#' raw moments
#'
#' @return mu : [vector] (length N corresponding to order N) central moments
#'
#' @author Ram Ahluwalia \email{rahluwalia@@gmail.com}
#'
#' @references
-#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management" \url{http://symmys.com/node/170},
-#' "E 16- Raw Moments to central moments".
+#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management"
+#' \url{http://symmys.com/node/170}, "E 16- Raw Moments to central moments".
#'
#' See Meucci's script for "Raw2Central.m"
#' @export
-Raw2Central = function( mu_ )
-{
- N = length( mu_ );
- mu = mu_;
- for( n in 1 : N )
- {
- mu[ n ] = ( (-1) ^ n ) * ( mu_[ 1 ] )^( n );
-
- for( k in 1 : (n-1) )
- {
- if( n != 1 ){ mu[ n ] = mu[ n ] + choose( n, k ) * ((-1)^(n-k)) * mu_[ k ] * (mu_[ 1 ])^(n-k) } ;
- }
+Raw2Central <- function(mu_) {
+ N <- length(mu_)
+ mu <- mu_
- mu[ n ] = mu[ n ] + mu_[ n ];
+ for (n in 1:N) {
+ mu[n] <- ((-1) ^ n) * (mu_[1]) ^ (n)
+
+ for(k in 1:(n - 1)) {
+ if (n != 1) {
+ mu[n] <- mu[n] + choose(n, k) * ((-1) ^ (n - k)) * mu_[k] *
+ (mu_[1]) ^ (n - k)
+ }
+ }
+ mu[n] <- mu[n] + mu_[n]
}
-
- return( mu = mu )
+ return(mu = mu)
}
-#' Map cumulative moments into raw moments.
+#' @title Map cumulative moments into raw moments.
#'
-#' Step 5 of the projection process:
-#'
-#' From the cumulants of Y we compute the raw non-central moments of Y
-#'
-#' We do so recursively by the identity in formula (24) which follows from applying (21) and re-arranging terms
+#' @description Step 5 of the projection process:\\
+#' From the cumulants of Y we compute the raw non-central moments of Y, and We
+#' do so recursively by the identity in formula (24) which follows from applying
+#' (21) and re-arranging terms
#'
-#' \deqn{ \tilde{ \mu } ^{ \big(n\big) }_{Y}
+#' @details \deqn{ \tilde{ \mu } ^{ \big(n\big) }_{Y}
#' \equiv \kappa^{ \big(n\big) }_{Y} + \sum_{k=1}^{n-1} (n-1)C_{k-1}
#' \kappa_{Y}^{ \big(k\big) } \tilde{ \mu } ^{n-k}_{Y} }
#'
-#' @param ka : [vector] (length N corresponding to order N) cumulative moments
+#' @param ka : [vector] (length N corresponding to order N) cumulative
+#' moments
#'
-#' @return mu_ : [vector] (length N corresponding to order N) corresponding raw moments
+#' @return mu_ : [vector] (length N corresponding to order N) corresponding
+#' raw moments
#'
#' @author Ram Ahluwalia \email{rahluwalia@@gmail.com}
#'
#' @references
-#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management" \url{http://symmys.com/node/170}.
-#' See Meucci's script for "Cumul2Raw.m".
+#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management"
+#' \url{http://symmys.com/node/170}. See Meucci's script for "Cumul2Raw.m".
+#' A. Meucci - "Annualization and General Projection of Skewness, Kurtosis and
+#' All Summary Statistics" - formula (24) \url{http://www.symmys.com/node/136}
#'
-#' A. Meucci - "Annualization and General Projection of Skewness, Kurtosis and All Summary Statistics" - formula (24) \url{http://www.symmys.com/node/136}
#' @export
-Cumul2Raw = function( ka )
-{
- N = length( ka );
- mu_ = ka;
+Cumul2Raw <- function(ka) {
+ N <- length(ka)
+ mu_ <- ka
- for( n in 1 : N )
- {
- ka[ n ] = mu_[ n ];
+ for (n in 1 : N) {
+ ka[n] <- mu_[n]
- for( k in 1 : (n-1) )
- {
- if( n != 1 ){ mu_[ n ] = mu_[ n ] + choose( n-1, k-1 ) * ka[ k ] * mu_[ n-k ] };
- }
+ for (k in 1:(n - 1)) {
+ if (n != 1) {
+ mu_[n] <- mu_[n] + choose(n - 1, k - 1) * ka[k] * mu_[n-k]
+ }
+ }
}
-
- return( mu_ );
+ return(mu_)
}
-#' Transforms raw moments into cumulants
+#' @title Transforms raw moments into cumulants
#'
-#' Step 3 of the projection process: From the non-central moments of X-t, we compute the cumulants.
-#'
-#'
-#' This process follows from the Taylor approximations for any small z and ln(1+x)~x for any small x,
-#' and from the definition of the first cumulant in (17). The we apply recursively the identity
-#' in formula (21). See Kendall and Stuart (1969)
+#' @description Step 3 of the projection process: From the non-central moments
+#' of X-t, we compute the cumulants.
+#' This process follows from the Taylor approximations for any small z and
+#' ln(1+x)~x for any small x, and from the definition of the first cumulant in
+#' (17). The we apply recursively the identity in formula (21). See Kendall and
+#' Stuart (1969)
#'
-#' \deqn{ \kappa^{ \big(n\big) }_{X} \equiv \tilde{ \mu } ^{ \big(n\big) }_{X} - \sum_{k=1}^{n-1} (n-1)C_{k-1} \kappa_{X}^{ \big(k\big) } \tilde{ \mu } ^{n-k}_{X} }
+#' @details \deqn{ \kappa^{ \big(n\big) }_{X} \equiv \tilde{ \mu } ^
+#' { \big(n\big) }_{X} - \sum_{k=1}^{n-1} (n-1)C_{k-1} \kappa_{X} ^
+#' { \big(k\big) } \tilde{ \mu } ^{n-k}_{X} }
#'
-#' @param mu_ : [vector] (length N corresponding to order N) corresponding raw moments
+#' @param mu_ : [vector] (length N corresponding to order N) corresponding
+#' raw moments
#'
-#' @return ka : [vector] (length N corresponding to order N) cumulative moments
+#' @return ka : [vector] (length N corresponding to order N) cumulative
+#' moments
#'
#' @author Ram Ahluwalia \email{rahluwalia@@gmail.com}
+#'
#' @references
-#' A. Meucci - "Annualization and General Projection of Skewness, Kurtosis and All Summary Statistics" - formula (21)
-#' Symmys site containing original MATLAB source code \url{http://www.symmys.com/node/136}
+#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management"
+#' \url{http://symmys.com/node/170}. See Meucci's script for "Cumul2Raw.m".
+#' A. Meucci - "Annualization and General Projection of Skewness, Kurtosis and
+#' All Summary Statistics" - formula (24) \url{http://www.symmys.com/node/136}
#'
-#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management" \url{http://symmys.com/node/170}.
-#' See Meucci's script for "Raw2Cumul.m"
#' @export
-Raw2Cumul = function( mu_ )
-{
- N = length( mu_ )
- ka = mu_
-
- for ( i in 1:N )
- {
- ka[i] = mu_[i];
+Raw2Cumul <- function(mu_) {
+ N <- length(mu_)
+ ka <- mu_
- for ( k in 1:(i-1) )
- {
- if ( i != 1 ) { ka[i] = ka[i] - choose(i-1,k-1) * ka[k] * mu_[i-k] }
+ for (i in 1:N) {
+ ka[i] <- mu_[i]
+
+ for (k in 1:(i - 1)) {
+ if (i != 1) {
+ ka[i] <- ka[i] - choose(i - 1,k - 1) * ka[k] * mu_[i - k]
+ }
}
}
-
- return( ka = ka )
+ return(ka = ka)
}
-#' Transforms first n central moments into first n raw moments (first central moment defined as expectation)
+#' @title Transforms first n central moments into first n raw moments (first
+#' central moment defined as expectation)
#'
-#' Step 2 of projection process: From the central moments of step 1, we compute the non-central moments. To do so we start
-#' with the first non-central moment and apply recursively an identity (formula 20)
+#' @description Step 2 of projection process: From the central moments of
+#' step 1, we compute the non-central moments. To do so we start with the first
+#' non-central moment and apply recursively an identity (formula 20)
#'
-#' \deqn{ \tilde{ \mu }^{ \big(1\big) }_{X} \equiv \mu ^{\big(1\big)}_{X}
-#' \\ \tilde{ \mu }^{ \big(n\big) }_{X} \equiv \mu ^{n}_{X} \sum_{k=0}^{n-1} \big(-1\big)^{n-k+1} \mu ^{n-k}_{X} \tilde{ \mu }^{\big(k\big)}_{X} }
+#' @details \deqn{ \tilde{ \mu }^{ \big(1\big) }_{X} \equiv
+#' \mu^{\big(1\big)}_{X} \\ \tilde{ \mu }^{ \big(n\big) }_{X} \equiv
+#' \mu^{n}_{X} \sum_{k=0}^{n-1} \big(-1\big)^{n-k+1} \mu^{n-k}_{X}
+#' \tilde{ \mu }^{\big(k\big)}_{X} }
#'
#' @param mu : [vector] (length N corresponding to order N) central moments
#'
-#' @return mu_ : [vector] (length N corresponding to order N) corresponding raw moments
+#' @return mu_ : [vector] (length N corresponding to order N) corresponding
+#' raw moments
#'
#' @author Ram Ahluwalia \email{rahluwalia@@gmail.com}
#'
#' @references
-#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management" \url{http://symmys.com/node/170},
-#' "E 16- Raw moments to central moments".
+#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management"
+#' \url{http://symmys.com/node/170}, "E 16- Raw moments to central moments".
#'
#' See Meucci's script for "Central2Raw.m"
#' @export
-Central2Raw = function( mu )
-{
- N = length( mu )
- mu_ = mu
-
- for ( i in 2:N ) # we use the notation 'i' instead of 'n' as in Meucci's code so we can browser
- {
- mu_[i] = ( ( -1 ) ^ ( i + 1 ) ) * ( mu[1] ) ^ (i)
-
- for ( k in 1:(i-1) )
- {
- mu_[i] = mu_[i] + choose(i,k) * ((-1)^(i-k+1)) * mu_[k] * (mu_[1])^(i-k)
- }
- mu_[i] = mu_[i] + mu[i]
+
+Central2Raw <- function(mu) {
+ N <- length(mu)
+ mu_ <- mu
+
+ for (i in 2:N) {
+ mu_[i] <- ((-1) ^ (i + 1)) * (mu[1]) ^ (i)
+
+ for (k in 1:(i - 1)) {
+ mu_[i] <- mu_[i] + choose(i,k) * ((-1) ^ (i - k + 1)) * mu_[k] *
+ (mu_[1]) ^ (i - k)
+ }
+ mu_[i] <- mu_[i] + mu[i]
}
-
- return ( mu_ = mu_ )
+ return (mu_ = mu_)
}
-#' Compute summary stats
+#' @title Compute summary stats
#'
-#' Step 0 in projection process: Compute summary stats (mean, skew, kurtosis, etc.) of the invariant X-t
-#' step 1 in the project process We collect the first 'n' central moments of the invariant X-t.
+#' @description Step 0 in projection process: Compute summary stats (mean, skew,
+#' kurtosis, etc.) of the invariant X-t step 1 in the project process We collect
+#' the first 'n' central moments of the invariant X-t.
#'
#' @param X an invariant
#' @param N the number of order statistics to collect
@@ -178,55 +186,71 @@
#' @export
#' @author Ram Ahluwalia \email{rahluwalia@@gmail.com}
#' @export
-SummStats = function( X , N )
-{
- suppressWarnings( library( matlab ) ) # otherwise, stops with "Error: (converted from warning) package 'matlab' was built under R version 2.13.2"
- library( moments )
-
+
+SummStats <- function(X, N) {
+ # otherwise, stops with "Error: (converted from warning) package 'matlab' was
+ # built under R version 2.13.2"
+ suppressWarnings(library(matlab))
+ library(moments)
+
# step 1: compute central moments based on formula 15
- mu = zeros( 1 , N )
- mu[ 1 ] = mean( X )
- for ( n in 2:N ) { mu[n] = moment( X , central = TRUE , order = n ) } # mu(2:n) contains the central moments. mu(1) is the mean
-
- # step 0: compute standardized statistics
- ga = mu
- ga[ 2 ] = sqrt( mu[2] )
- for ( n in 3:N ) # we focus on case n >= 3 because from the definition of central moments (15) and from (3) that i) the first central moment is the mean of the invariant X-t, and ii) the second central moment is standard deviaiton of the of the invariant X-t
- {
- ga[n] = mu[n] / ( ga[2]^n ) # based on formula 19. ga[1] = mean of invariant X-t, ga[2] = sd, ga[3] = skew, ga[4] = kurtosis...
- }
-
- return( list( ga = ga , mu = mu ) )
+ mu <- zeros(1, N)
+ mu[1] <- mean(X)
+ # mu(2:n) contains the central moments. mu(1) is the mean
+ for (n in 2:N) {
+ mu[n] <- moment(X, central = TRUE, order = n)
+ }
+
+ # step 0: compute standardized statistics
+ ga <- mu
+ ga[2] <- sqrt(mu[2])
+ # we focus on case n >= 3 because from the definition of central moments (15)
+ # and from (3) that i) the first central moment is the mean of the invariant
+ # X-t, and ii) the second central moment is standard deviaiton of the of the
+ # invariant X-t
+ for (n in 3:N) {
+ # based on formula 19. ga[1] = mean of invariant X-t, ga[2] = sd,
+ # ga[3] = skew, ga[4] = kurtosis...
+ ga[n] <- mu[n] / (ga[2]^n)
+ }
+ return(list(ga = ga, mu = mu))
}
-#' Calculates the population standard deviation
+#' @title Calculates the population standard deviation
#'
-#' Calculates the population standard deviation dividing by 'n' instead of 'n-1' equivalent to Matlab
+#' @description Calculates the population standard deviation dividing by 'n'
+#' instead of 'n-1' equivalent to Matlab
#'
#' @param x a generic numeric vector
-#' @return std a numeric with the population standard deviaiton of the generic numeric
+#' @return std a numeric with the population standard deviaiton of the generic
+#' numeric
+#' @author Ram Ahluwalia \email{rahluwalia@@gmail.com}
+#'
#' @export
-#' @author Ram Ahluwalia \email{rahluwalia@@gmail.com}
-std = function( x ) { ( sum( ( x - mean( x ) ) ^ 2 ) / length( x ) ) ^.5 }
+std <- function(x) {
+ return(sum((x - mean(x)) ^ 2) / length(x)) ^ .5
+}
+
#' Generate arbitrary distribution of a shifted-lognormal invariant
#'
-#' \deqn{X = a + e^{ m + sZ }} (formula 14)
+#' @details \deqn{X = a + e^{ m + sZ }} (formula 14)
#'
#' @param J a numeric with the number of scenarios
-#' @param a a numeric with the location shift parameter. Mean of distribution will be exp(a)
+#' @param a a numeric with the location shift parameter. Mean of
+#' distribution will be exp(a)
#' @param m log of the mean of the distribution
#' @param s log of the standard deviation of the distribution
#'
-#' @return X a numeric vector with i.i.d. lognormal samples based on parameters J, a, m, and s where X = a + exp( m + s * Z )
+#' @return X a numeric vector with i.i.d. lognormal samples based on
+#' parameters J, a, m, and s where X = a + exp(m + s * Z)
#' @author Ram Ahluwalia \email{rahluwalia@@gmail.com}
#' @export
-GenerateLogNormalDistribution = function( J, a, m, s )
-{
- Z = rnorm( J / 2 , 0 , 1 ) # create J/2 draws from the standard normal
- Z = c( Z , -Z) / std( Z ) # a Jx1 numeric vector
- X = a + exp( m + s * Z ) # a Jx1 numeric vector
-
- return( X = X )
+
+GenerateLogNormalDistribution <- function(J, a, m, s) {
+ Z <- rnorm(J / 2, 0, 1) # create J/2 draws from the standard normal
+ Z <- c(Z, -Z) / std(Z) # a Jx1 numeric vector
+ X <- a + exp(m + s * Z) # a Jx1 numeric vector
+
+ return(X = X)
}
-
Modified: pkg/Meucci/R/LognormalCopulaPdf.R
===================================================================
--- pkg/Meucci/R/LognormalCopulaPdf.R 2015-08-20 09:49:54 UTC (rev 3983)
+++ pkg/Meucci/R/LognormalCopulaPdf.R 2015-08-20 15:52:02 UTC (rev 3984)
@@ -1,7 +1,9 @@
-#' @title Computes the pdf of the copula of the lognormal distribution at the generic point u in the unit hypercube.
+#' @title Computes the pdf of the copula of the lognormal distribution at the
+#' generic point u in the unit hypercube.
#'
-#' @description Computes the pdf of the copula of the lognormal distribution at the generic point u in the unit hypercube,
-#' as described in A. Meucci, "Risk and Asset Allocation", Springer, 2005.
+#' @description Computes the pdf of the copula of the lognormal distribution at
+#' the generic point u in the unit hypercube, as described in A. Meucci,
+#' "Risk and Asset Allocation", Springer, 2005.
#'
#' @param u [vector] (J x 1) grades
#' @param Mu [vector] (N x 1) location parameter
@@ -10,29 +12,29 @@
#' @return F_U [vector] (J x 1) PDF values
#'
#' @references
-#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management" \url{http://symmys.com/node/170},
-#' "E 36 - Pdf of the lognormal copula".
+#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management"
+#' \url{http://symmys.com/node/170}, "E 36 - Pdf of the lognormal copula".
#'
#' See Meucci's script for "LognormalCopulaPdf.m"
#'
#' @author Xavier Valls \email{xaviervallspla@@gmail.com}
#' @export
-LognormalCopulaPdf = function( u, Mu, Sigma )
-{
- N = length( u );
- s = sqrt( diag( Sigma ));
+LognormalCopulaPdf <- function(u, Mu, Sigma) {
+ N <- length(u)
+ s <- sqrt(diag(Sigma))
- x = qlnorm( u, Mu, s );
+ x <- qlnorm(u, Mu, s)
- Numerator = ( 2 * pi ) ^ ( -N / 2 ) * ( (det ( Sigma ) ) ^ ( -0.5 ) ) /
- prod(x) * exp( -0.5 * t(log(x) - Mu) %*% mldivide( Sigma , ( log( x ) - Mu ), pinv=FALSE ) );
+ Numerator <- (2 * pi) ^ (-N / 2) * ((det (Sigma)) ^ (-0.5)) /
+ prod(x) * exp(-0.5 * t(log(x) - Mu) %*%
+ mldivide(Sigma, (log(x) - Mu), pinv=FALSE))
- fs = dlnorm( x, Mu, s);
+ fs <- dlnorm(x, Mu, s)
- Denominator = prod(fs);
+ Denominator <- prod(fs)
- F_U = Numerator / Denominator;
+ F_U <- Numerator / Denominator
- return ( F_U );
-}
\ No newline at end of file
+ return (F_U)
+}
Modified: pkg/Meucci/R/NormalCopulaPdf.R
===================================================================
--- pkg/Meucci/R/NormalCopulaPdf.R 2015-08-20 09:49:54 UTC (rev 3983)
+++ pkg/Meucci/R/NormalCopulaPdf.R 2015-08-20 15:52:02 UTC (rev 3984)
@@ -1,7 +1,9 @@
-#' @title Computes the pdf of the copula of the normal distribution at the generic point u in the unit hypercube
+#' @title Computes the pdf of the copula of the normal distribution at the
+#' generic point u in the unit hypercube
#'
-#' @description Computes the pdf of the copula of the normal distribution at the generic point u in the unit
-#' hypercube, as described in A. Meucci, "Risk and Asset Allocation", Springer, 2005.
+#' @description Computes the pdf of the copula of the normal distribution at the
+#' generic point u in the unit hypercube, as described in A. Meucci, "Risk and
+#' Asset Allocation", Springer, 2005.
#'
#' @param u [vector] (J x 1) grade
#' @param Mu [vector] (N x 1) mean
@@ -10,28 +12,28 @@
#' @return F_U [vector] (J x 1) PDF values
#'
#' @references
-#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management" \url{http://symmys.com/node/170},
-#' "E 33 - Pdf of the normal copula".
+#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management"
+#' \url{http://symmys.com/node/170}, "E 33 - Pdf of the normal copula".
#'
#' See Meucci's script for "NormalCopulaPdf.m"
#'
#' @author Xavier Valls \email{xaviervallspla@@gmail.com}
#' @export
-NormalCopulaPdf = function( u, Mu, Sigma )
-{
- N = length( u );
- s = sqrt( diag( Sigma ));
+NormalCopulaPdf <- function(u, Mu, Sigma) {
+ N <- length(u)
+ s <- sqrt(diag(Sigma))
- x = qnorm( u, Mu, s );
+ x <- qnorm(u, Mu, s)
- Numerator = ( 2 * pi ) ^ ( -N / 2 ) * ( (det ( Sigma ) ) ^ ( -0.5 ) ) * exp( -0.5 * t(x - Mu) %*% mldivide( Sigma , ( x - Mu )));
+ Numerator <- (2 * pi) ^ (-N / 2) * ((det (Sigma)) ^ (-0.5)) *
+ exp(-0.5 * t(x - Mu) %*% mldivide(Sigma, (x - Mu)))
- fs = dnorm( x, Mu, s);
+ fs <- dnorm(x, Mu, s)
- Denominator = prod(fs);
+ Denominator <- prod(fs)
- F_U = Numerator / Denominator;
+ F_U <- Numerator / Denominator
- return ( F_U );
-}
\ No newline at end of file
+ return (F_U)
+}
Modified: pkg/Meucci/R/PerformIidAnalysis.R
===================================================================
--- pkg/Meucci/R/PerformIidAnalysis.R 2015-08-20 09:49:54 UTC (rev 3983)
+++ pkg/Meucci/R/PerformIidAnalysis.R 2015-08-20 15:52:02 UTC (rev 3984)
@@ -1,7 +1,8 @@
#' @title Performs simple invariance (i.i.d.) tests on a time series.
#'
-#' @description This function performs simple invariance (i.i.d.) tests on a time series, as described in
-#' A. Meucci "Risk and Asset Allocation", Springer, 2005
+#' @description This function performs simple invariance (i.i.d.) tests on a
+#' time series, as described in A. Meucci "Risk and Asset Allocation", Springer,
+#' 2005
#'
#' @param Dates : [vector] (T x 1) dates
#' @param Data : [matrix] (T x N) data
@@ -9,54 +10,55 @@
#'
#' @note it checks the evolution over time
#'
-#' it checks that the variables are identically distributed by looking at the histogram of two subsamples
+#' it checks that the variables are identically distributed by looking at the
+#' histogram of two subsamples
#'
-#' it checks that the variables are independent by looking at the 1-lag scatter plot
+#' it checks that the variables are independent by looking at the 1-lag
+#' scatter plot
#'
#' under i.i.d. the location-dispersion ellipsoid should be a circle
#'
#' @references
-#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management" \url{http://symmys.com/node/170}.
+#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management"
+#' \url{http://symmys.com/node/170}. See Meucci's script for
+#' "PerformIidAnalysis.m"
#'
-#' See Meucci's script for "PerformIidAnalysis.m"
-#'
#' @author Xavier Valls \email{xaviervallspla@@gmail.com}
#' @export
-PerformIidAnalysis = function( Dates = dim( Data)[1], Data, Str = "")
-{
+PerformIidAnalysis <- function(Dates = dim(Data)[1], Data, Str = "") {
- ##########################################################################################################
+ ############################################################################
### Time series over time
- dev.new();
- plot( Dates, Data, main = Str );
- #datetick( 'x', 'mmmyy', 'keeplimits', 'keepticks' );
-
+ dev.new()
+ plot(Dates, Data, main = Str)
+ #datetick('x', 'mmmyy', 'keeplimits', 'keepticks')
- ##########################################################################################################
- ### Test "identically distributed hypothesis": split observations into two sub-samples and plot histogram
- Sample_1 = Data[ 1:round(length(Data) / 2) ];
- Sample_2 = Data[(round(length(Data)/2) + 1) : length(Data) ];
- num_bins_1 = round(5 * log(length(Sample_1)));
- num_bins_2 = round(5 * log(length(Sample_2)));
- X_lim = c( min(Data) - .1 * (max(Data) - min(Data)), max(Data) + .1 * (max(Data) - min(Data)));
+ ############################################################################
+ ### Test "identically distributed hypothesis": split observations into two
+ ### sub-samples and plot histogram
+ Sample_1 <- Data[1:round(length(Data) / 2)]
+ Sample_2 <- Data[(round(length(Data)/2) + 1) : length(Data)]
+ num_bins_1 <- round(5 * log(length(Sample_1)))
+ num_bins_2 <- round(5 * log(length(Sample_2)))
+ X_lim <- c(min(Data) - .1 * (max(Data) - min(Data)), max(Data) + .1 *
+ (max(Data) - min(Data)))
- dev.new();
+ dev.new()
- layout( matrix(c(1,1,2,2,0,3,3,0), 2, 4, byrow = TRUE), heights=c(1,1,1));
- hist(Sample_1, num_bins_1, xlab = "", ylab = "", main = "first half" );
- hist(Sample_2, num_bins_2, xlab = "", ylab = "", main = "second half" );
-
- ##########################################################################################################
- ### Test "independently distributed hypothesis": scatter plot of observations at lagged times
-
+ layout(matrix(c(1,1,2,2,0,3,3,0), 2, 4, byrow = TRUE), heights=c(1,1,1))
+ hist(Sample_1, num_bins_1, xlab = "", ylab = "", main = "first half")
+ hist(Sample_2, num_bins_2, xlab = "", ylab = "", main = "second half")
- X = Data[ 1 : length(Data)-1 ];
- Y = Data[ 2 : length(Data) ];
- plot(X, Y, main="changes in implied vol");
+ ############################################################################
+ ### Test "independently distributed hypothesis": scatter plot of
+ ### observations at lagged times
- m = cbind( apply( cbind( X, Y ), 2, mean ));
- S = cov( cbind( X, Y ));
- TwoDimEllipsoid( m, S, 2, FALSE, FALSE);
+ X <- Data[1 : length(Data)-1]
+ Y <- Data[2 : length(Data)]
+ plot(X, Y, main = "changes in implied vol")
-}
\ No newline at end of file
+ m <- cbind(apply(cbind(X, Y), 2, mean))
+ S <- cov(cbind(X, Y))
+ TwoDimEllipsoid(m, S, 2, FALSE, FALSE)
+}
Modified: pkg/Meucci/R/StudentTCopulaPdf.R
===================================================================
--- pkg/Meucci/R/StudentTCopulaPdf.R 2015-08-20 09:49:54 UTC (rev 3983)
+++ pkg/Meucci/R/StudentTCopulaPdf.R 2015-08-20 15:52:02 UTC (rev 3984)
@@ -1,7 +1,9 @@
-#' @title Pdf of the copula of the Student t distribution at the generic point u in the unit hypercube
+#' @title Pdf of the copula of the Student t distribution at the generic point u
+#' in the unit hypercube
#'
-#' @description Pdf of the copula of the Student t distribution at the generic point u in the unit hypercube,
-#' as described in A. Meucci, "Risk and Asset Allocation", Springer, 2005.
+#' @description Pdf of the copula of the Student t distribution at the generic
+#' point u in the unit hypercube, as described in A. Meucci, "Risk and Asset
+#' Allocation", Springer, 2005.
#'
#' @param u : [vector] (J x 1) grade
#' @param nu : [numerical] degrees of freedom
@@ -12,31 +14,31 @@
#' @return F_U : [vector] (J x 1) PDF values
#'
#' @references
-#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management" \url{http://symmys.com/node/170},
-#' "E 88 - Copula vs. Correlation".
+#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management"
+#' \url{http://symmys.com/node/170}, "E 88 - Copula vs. Correlation".
#'
#' See Meucci's script for "StudentTCopulaPdf.m"
#'
#' @author Xavier Valls \email{xaviervallspla@@gmail.com}
#' @export
-StudentTCopulaPdf = function( u, nu, Mu, Sigma )
-{
- N = length( u );
- s = sqrt( diag( Sigma ));
+StudentTCopulaPdf <- function(u, nu, Mu, Sigma) {
+ N <- length(u)
+ s <- sqrt(diag(Sigma))
- x = Mu + s * qt( u, nu);
+ x <- Mu + s * qt(u, nu)
+ #z2 <- t(x - Mu) %*% inv(Sigma) * (x-Mu)
+ z2 <- t(x - Mu) %*% mldivide(Sigma, (x - Mu))
+ K <- (nu * pi) ^ (-N / 2) * gamma((nu + N) / 2) / gamma(nu / 2) *
+ ((det(Sigma)) ^ (-0.5))
+ Numerator <- K * (1 + z2 / nu) ^ (-(nu + N) / 2)
- z2 = t(x - Mu) %*% mldivide( Sigma, (x - Mu)); #z2 = t(x - Mu) %*% inv(Sigma) * (x-Mu);
- K = ( nu * pi )^( -N / 2 ) * gamma( ( nu + N ) / 2 ) / gamma( nu / 2 ) * ( ( det( Sigma ) )^( -0.5 ));
- Numerator = K * (1 + z2 / nu)^(-(nu + N) / 2);
-
- fs = dt((x - Mu) / s , nu);
+ fs <- dt((x - Mu) / s, nu)
- Denominator = prod(fs);
+ Denominator <- prod(fs)
- F_U = Numerator / Denominator;
+ F_U <- Numerator / Denominator
- return ( F_U );
+ return (F_U)
}
Modified: pkg/Meucci/R/TwoDimEllipsoid.R
===================================================================
--- pkg/Meucci/R/TwoDimEllipsoid.R 2015-08-20 09:49:54 UTC (rev 3983)
+++ pkg/Meucci/R/TwoDimEllipsoid.R 2015-08-20 15:52:02 UTC (rev 3984)
@@ -1,20 +1,29 @@
-#'@title Computes the location-dispersion ellipsoid of the normalized first diagonal and off-diagonal elements
-#' of a 2x2 Wishart distribution as a function of the inputs
+#' @title Computes the location-dispersion ellipsoid of the normalized first
+#' diagonal and off-diagonal elements of a 2x2 Wishart distribution as a
+#' function of the inputs
#'
-#' @description This function computes the location-dispersion ellipsoid of the normalized (unit variance,
-#' zero expectation)first diagonal and off-diagonal elements of a 2x2 Wishart distribution as a function
-#' of the inputs, as described in A. Meucci, "Risk and Asset Allocation", Springer, 2005, Chapter 2.
+#' @description This function computes the location-dispersion ellipsoid of the
+#' normalized (unit variance, zero expectation)first diagonal and off-diagonal
+#' elements of a 2x2 Wishart distribution as a function of the inputs, as
+#' described in A. Meucci, "Risk and Asset Allocation", Springer, 2005, Ch 2.
#'
-#' @param Location : [vector] (2 x 1) location vector (typically the expected value
-#' @param Square_Dispersion : [matrix] (2 x 2) scatter matrix Square_Dispersion (typically the covariance matrix)
-#' @param Scale : [scalar] a scalar Scale, that specifies the scale (radius) of the ellipsoid
-#' @param PlotEigVectors : [boolean] true then the eigenvectors (=principal axes) are plotted
-#' @param PlotSquare : [boolean] true then the enshrouding box is plotted. If Square_Dispersion is the covariance
+#' @param Location : [vector] (2x1) location vector (typically the
+#' expected value
+#' @param Square_Dispersion : [matrix] (2x2) scatter matrix Square_Dispersion
+#' (typically the covariance matrix)
+#' @param Scale : [scalar] a scalar Scale, that specifies the
+#' scale (radius) of the ellipsoid
+#' @param PlotEigVectors : [boolean] true then the eigenvectors (=principal
+#' axes) are plotted
+#' @param PlotSquare : [boolean] true then the enshrouding box is
+#' plotted. If Square_Dispersion is the
+#' covariance
#'
#' @return E : [figure handle]
#'
#' @references
-#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management" \url{http://symmys.com/node/170}.
+#' A. Meucci - "Exercises in Advanced Risk and Portfolio Management"
+#' \url{http://symmys.com/node/170}.
#'
#' See Meucci's script for "TwoDimEllipsoid.m"
#'
@@ -22,81 +31,81 @@
#' @export
-TwoDimEllipsoid = function( Location, Square_Dispersion, Scale = 1, PlotEigVectors = FALSE, PlotSquare = FALSE )
-{
+TwoDimEllipsoid <- function(Location, Square_Dispersion, Scale = 1,
+ PlotEigVectors = FALSE, PlotSquare = FALSE) {
- ##########################################################################################################
- ### compute the ellipsoid in the r plane, solution to ((R-Location)' * Dispersion^-1 * (R-Location) ) = Scale^2
-
- Eigen = eigen(Square_Dispersion);
- Centered_Ellipse = c();
- Angle = seq( 0, 2*pi, pi/500 );
- NumSteps = length(Angle);
-
- for( i in 1 : NumSteps )
- {
- # normalized variables (parametric representation of the ellipsoid)
- y = rbind( cos( Angle[ i ] ), sin( Angle[ i ] ) );
- Centered_Ellipse = c( Centered_Ellipse, Eigen$vectors %*% diag(sqrt(Eigen$values), length(Eigen$values)) %*% y ); ##ok<AGROW>
+ ############################################################################
+ ### compute the ellipsoid in the r plane, solution to ((R-Location)' *
+ ### Dispersion^-1 * (R-Location)) = Scale^2
+
+ Eigen <- eigen(Square_Dispersion)
+ Centered_Ellipse <- c()
+ Angle <- seq(0, 2 * pi, pi / 500)
+ NumSteps <- length(Angle)
+
+ for (i in 1 : NumSteps) {
+ # normalized variables (parametric representation of the ellipsoid)
+ y <- rbind(cos(Angle[i]), sin(Angle[i]))
+ Centered_Ellipse <- c(Centered_Ellipse, Eigen$vectors %*%
+ diag(sqrt(Eigen$values), length(Eigen$values)) %*% y)
}
- R = Location %*% array( 1, NumSteps ) + Scale * Centered_Ellipse;
+ R <- Location %*% array(1, NumSteps) + Scale * Centered_Ellipse
- ##########################################################################################################
- ### Plot the ellipsoid
-
- E = lines( R[1, ], R[2, ], col = "red", lwd = 2 );
+ ############################################################################
+ ### Plot the ellipsoid
- ##########################################################################################################
- ### Plot a rectangle centered in Location with semisides of lengths Dispersion[ 1] and Dispersion[ 2 ], respectively
-
- if( PlotSquare )
- {
- Dispersion = sqrt( diag( Square_Dispersion ) );
- Vertex_LowRight_A = Location[ 1 ] + Scale * Dispersion[ 1 ];
- Vertex_LowRight_B = Location[ 2 ] - Scale * Dispersion[ 2 ];
- Vertex_LowLeft_A = Location[ 1 ] - Scale * Dispersion[ 1 ];
- Vertex_LowLeft_B = Location[ 2 ] - Scale * Dispersion[ 2 ];
- Vertex_UpRight_A = Location[ 1 ] + Scale * Dispersion[ 1 ];
- Vertex_UpRight_B = Location[ 2 ] + Scale * Dispersion[ 2 ];
- Vertex_UpLeft_A = Location[ 1 ] - Scale * Dispersion[ 1 ];
- Vertex_UpLeft_B = Location[ 2 ] + Scale * Dispersion[ 2 ];
-
- Square = rbind( c( Vertex_LowRight_A, Vertex_LowRight_B ),
- c( Vertex_LowLeft_A, Vertex_LowLeft_B ),
- c( Vertex_UpLeft_A, Vertex_UpLeft_B ),
- c( Vertex_UpRight_A, Vertex_UpRight_B ),
- c( Vertex_LowRight_A, Vertex_LowRight_B ) );
+ E <- lines(R[1, ], R[2, ], col = "red", lwd = 2)
- h = lines(Square[ , 1 ], Square[ , 2 ], col = "red", lwd = 2 );
-
- }
+ ############################################################################
+ ### Plot a rectangle centered in Location with semisides of lengths
+ ### Dispersion[1] and Dispersion[2], respectively
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/returnanalytics -r 3984
More information about the Returnanalytics-commits
mailing list