Ожидаемое возвращение и ковариация из временного ряда возврата

Question

Я пытаюсь смоделировать определенную здесь функцию Matlab ewstats:

https://it.mathworks.com/help/finance/ewstats.html

Результаты, предоставленные Matlab, следующие:

> ExpReturn = 1×2
0.1995    0.1002

> ExpCovariance = 2×2
0.0032   -0.0017
-0.0017    0.0010

Я пытаюсь скопировать пример с пакетом RiskPortfolios R:

https://cran.r -project.org / web / packages / RiskPortfolios / RiskPortfolios.pdf

Я использую код R:

library(RiskPortfolios)

rets <- as.matrix(cbind(c(0.24, 0.15, 0.27, 0.14), c(0.08, 0.13, 0.06, 0.13)))
w <- 0.98

rets
w
meanEstimation(rets, control = list(type = 'ewma', lambda = w))
covEstimation(rets, control = list(type = 'ewma', lambda = w))

Средняя оценка такая же, как в примере, но ковариационная матрица отличается:

> rets
     [,1] [,2]
[1,] 0.24 0.08
[2,] 0.15 0.13
[3,] 0.27 0.06
[4,] 0.14 0.13
> w
[1] 0.98
> 
> meanEstimation(rets, control = list(type = 'ewma', lambda = w))
[1] 0.1995434 0.1002031
> 
> covEstimation(rets, control = list(type = 'ewma', lambda = w))
             [,1]         [,2]
[1,]  0.007045044 -0.003857217
[2,] -0.003857217  0.002123827

Я что-то пропустил?Спасибо

Answer 1

Они дают один и тот же ответ, если используется type = "lw":

round(covEstimation(rets, control = list(type = 'lw')), 4)

##   0.0032 -0.0017
##  -0.0017  0.0010

Answer 2

Для тех, кому в конечном итоге требуется эквивалентная ewstats функция в R, вот код, который я написал:

ewstats <- function(RetSeries, DecayFactor=NULL, WindowLength=NULL){

  #EWSTATS Expected return and covariance from return time series.  
  #   Optional exponential weighting emphasizes more recent data.
  # 
  #   [ExpReturn, ExpCovariance, NumEffObs] = ewstats(RetSeries, ...
  #                                           DecayFactor, WindowLength)
  #  
  #   Inputs:
  #     RetSeries : NUMOBS by NASSETS matrix of equally spaced incremental 
  #     return observations.  The first row is the oldest observation, and the
  #     last row is the most recent.
  #
  #     DecayFactor : Controls how much less each observation is weighted than its
  #     successor.  The k'th observation back in time has weight DecayFactor^k.
  #     DecayFactor must lie in the range: 0 < DecayFactor <= 1. 
  #     The default is DecayFactor = 1, which is the equally weighted linear
  #     moving average Model (BIS).  
  #
  #     WindowLength: The number of recent observations used in
  #     the computation.  The default is all NUMOBS observations.
  #
  #   Outputs:
  #     ExpReturn : 1 by NASSETS estimated expected returns.
  # 
  #     ExpCovariance : NASSETS by NASSETS estimated covariance matrix.  
  #
  #     NumEffObs: The number of effective observations is given by the formula:
  #     NumEffObs = (1-DecayFactor^WindowLength)/(1-DecayFactor).  Smaller
  #     DecayFactors or WindowLengths emphasize recent data more strongly, but
  #     use less of the available data set.
  #
  #   The standard deviations of the asset return processes are given by:
  #   STDVec = sqrt(diag(ECov)).  The correlation matrix is :
  #   CorrMat = VarMat./( STDVec*STDVec' )
  #
  #   See also MEAN, COV, COV2CORR.


  NumObs <- dim(RetSeries)[1]
  NumSeries <- dim(RetSeries)[2]

  # size the series and the window
  if (is.null(WindowLength)) {
    WindowLength <- NumObs
  }

  if (is.null(DecayFactor)) {
    DecayFactor = 1
  }


  if (DecayFactor <= 0 | DecayFactor > 1) {
    stop('Must have 0< decay factor <= 1.')
  }


  if (WindowLength > NumObs){
    stop(sprintf('Window Length #d must be <= number of observations #d',
                 WindowLength, NumObs))
  }


  # ------------------------------------------------------------------------
  # size the data to the window
  RetSeries <- RetSeries[NumObs-WindowLength+1:NumObs, ]

  # Calculate decay coefficients
  DecayPowers <- seq(WindowLength-1, 0, by = -1)
  VarWts <- sqrt(DecayFactor)^DecayPowers
  RetWts <- (DecayFactor)^DecayPowers

  NEff = sum(RetWts) # number of equivalent values in computation

  # Compute the exponentially weighted mean return
  WtSeries <- matrix(rep(RetWts, times = NumSeries),
                     nrow = length(RetWts), ncol = NumSeries) * RetSeries

  ERet <- colSums(WtSeries)/NEff;

  # Subtract the weighted mean from the original Series
  CenteredSeries <- RetSeries - matrix(rep(ERet, each = WindowLength),
                                       nrow = WindowLength, ncol = length(ERet))

  # Compute the weighted variance
  WtSeries <- matrix(rep(VarWts, times = NumSeries),
                     nrow = length(VarWts), ncol = NumSeries) * CenteredSeries

  ECov <- t(WtSeries) %*% WtSeries / NEff

  list(ExpReturn = ERet, ExpCovariance = ECov, NumEffObs = NEff)
}

Answer 3

Они используют разные алгоритмы.Из руководства RiskPortfolio:

ewma ... См. RiskMetrics (1996)

На странице справки Matlab:

Тамнет никакой связи между функцией ewstats и подходом RiskMetrics® для определения ожидаемой доходности и ковариации из временного ряда возврата.

К сожалению, Matlab не сообщает нам, какой алгоритм используется.