1# Compute co-moment matrices
2
3#' calculate centered Returns
4#'
5#' the \eqn{n}-th centered moment is calculated as \deqn{ }{moment^n(R) =
6#' E[R-E(R)^n]}\deqn{ \mu^{(n)}(R) = E\lbrack(R-E(R))^n\rbrack }{moment^n(R) =
7#' E[R-E(R)^n]}
8#'
9#' These functions are used internally by PerformanceAnalytics to calculate
10#' centered moments for a multivariate distribution as well as the standardized
11#' moments of a portfolio distribution.  They are exposed here for users who
12#' wish to use them directly, and we'll get more documentation written when we
13#' can.
14#'
15#' These functions were first utilized in Boudt, Peterson, and Croux (2008),
16#' and have been subsequently used in our other research.
17#'
18#' ~~ Additional Details will be added to documentation as soon as we have time
19#' to write them. Documentation Patches Welcome. ~~
20#'
21#' @aliases centeredcomoment centeredmoment Return.centered
22#' @param R an xts, vector, matrix, data frame, timeSeries or zoo object of
23#' asset returns
24#' @param Ra an xts, vector, matrix, data frame, timeSeries or zoo object of
25#' asset returns
26#' @param Rb an xts, vector, matrix, data frame, timeSeries or zoo object of
27#' index, benchmark, portfolio, or secondary asset returns to compare against
28#' @param power power or moment to calculate
29#' @param p1 first power of the comoment
30#' @param p2 second power of the comoment
31#' @param normalize whether to standardize the calculation to agree with common
32#' usage, or leave the default mathematical meaning
33#' @param \dots any other passthru parameters
34#' @author Kris Boudt and Brian Peterson
35#' @references Boudt, Kris, Brian G. Peterson, and Christophe Croux. 2008.
36#' Estimation and Decomposition of Downside Risk for Portfolios with Non-Normal
37#' Returns. Journal of Risk. Winter.
38#'
39#' Martellini, L. and Ziemann, V., 2010. Improved estimates of higher-order
40#' comoments and implications for portfolio selection. Review of Financial
41#' Studies, 23(4):1467-1502.
42#'
43#' Ranaldo, Angelo, and Laurent Favre Sr. 2005. How to Price Hedge Funds: From
44#' Two- to Four-Moment CAPM. SSRN eLibrary.
45#'
46#' Scott, Robert C., and Philip A. Horvath. 1980. On the Direction of
47#' Preference for Moments of Higher Order than the Variance. Journal of Finance
48#' 35(4):915-919.
49###keywords ts multivariate distribution models
50#' @examples
51#'
52#'
53#' data(managers)
54#' Return.centered(managers[,1:3,drop=FALSE])
55#'
56#' @rdname centeredmoments
57#' @export
58Return.centered <-
59function (R,...)
60{ # @author Peter Carl and Kris Boudt
61
62    # DESCRIPTION:
63    # Calculates the returns less the mean return of the asset
64
65    # Inputs:
66    # R: a matrix, data frame, or timeSeries of returns
67
68    # Outputs:
69    # A timeseries of the calculated series
70
71    # FUNCTION:
72
73    # Transform input data to a timeseries (zoo) object
74    R = checkData(R, method="zoo")
75
76    # Get dimensions and labels
77    columns.a = ncol(R)
78    rows.a = nrow(R)
79
80    if(columns.a==1){
81       R.centered = zoo(NA);
82       R.mean = zoo(NA);
83       R.mean = mean(R[, drop=FALSE])
84       R.centered = R[ , drop=FALSE] - R.mean
85    }else{
86       R.mean = apply(R,2,'mean', na.rm=TRUE)
87       # returns a vector holding the mean return for each asset
88
89       R.centered = R - matrix( rep(R.mean,rows.a), ncol= columns.a, byrow=TRUE)
90       # return the matrix of centered returns
91   }
92
93
94   # RESULTS:
95    return(R.centered)
96}
97
98###############################################################################
99
100#' @rdname CoMoments
101#' @name CoMoments
102NULL
103
104#' @rdname CoMoments
105#' @export
106CoSkewnessMatrix <-
107function (R, ...)
108{ # @author Kris Boudt
109    return(M3.MM(R))
110}
111
112###############################################################################
113
114#' @rdname CoMoments
115#' @export
116CoKurtosisMatrix <-
117function (R, ...)
118{ # @author Kris Boudt
119    return(M4.MM(R))
120}
121
122###############################################################################
123
124
125#' @rdname centeredmoments
126#' @export
127centeredmoment = function(R,power)
128{# @author Kris Boudt, Peter Carl
129    R = checkData(R)
130    out =  apply(Return.centered(R)^power,2,FUN=mean, na.rm=TRUE)
131    return(out);
132}
133
134###############################################################################
135
136#' @rdname centeredmoments
137#' @export
138centeredcomoment = function(Ra,Rb,p1,p2,normalize=FALSE)
139{# @author Kris Boudt, Peter Carl, and Brian G. Peterson
140
141    Ra = checkData(Ra); Rb = checkData(Rb);
142
143    out = mean( na.omit( Return.centered(Ra)^p1 * Return.centered(Rb)^p2))
144
145    if(normalize) {
146        out=out/ as.numeric(centeredmoment(Rb,power=(p1+p2))) #
147    }
148    return(out);
149}
150
151
152###############################################################################
153
154#' Functions for calculating comoments of financial time series
155#'
156#' calculates coskewness and cokurtosis as the skewness and kurtosis of two
157#' assets with reference to one another.
158#'
159#' Ranaldo and Favre (2005) define coskewness and cokurtosis as the skewness
160#' and kurtosis of a given asset analysed with the skewness and kurtosis of the
161#' reference asset or portfolio.  Adding an asset to a portfolio, such as a
162#' hedge fund with a significant level of coskewness to the portfolio, can
163#' increase or decrease the resulting portfolio's skewness. Similarly, adding a
164#' hedge fund with a positive cokurtosis coefficient will add kurtosis to the
165#' portfolio.
166#'
167#' The co-moments are useful for measuring the marginal contribution of each
168#' asset to the portfolio's resulting risk.  As such, comoments of asset return
169#' distribution should be useful as inputs for portfolio optimization in
170#' addition to the covariance matrix.  Martellini and Ziemann (2007) point out
171#' that the problem of portfolio selection becomes one of selecting tangency
172#' points in four dimensions, incorporating expected return, second, third and
173#' fourth centered moments of asset returns.
174#'
175#' Even outside of the optimization problem, measuring the co-moments should be
176#' a useful tool for evaluating whether or not an asset is likely to provide
177#' diversification potential to a portfolio, not only in terms of normal risk
178#' (i.e. volatility) but also the risk of assymetry (skewness) and extreme
179#' events (kurtosis).
180#' @name CoMoments
181#' @concept co-moments
182#' @concept moments
183#' @aliases CoMoments CoVariance CoSkewness CoKurtosis
184#' @param R an xts, vector, matrix, data frame, timeSeries or zoo object of
185#' asset returns
186#' @param Ra an xts, vector, matrix, data frame, timeSeries or zoo object of
187#' asset returns
188#' @param Rb an xts, vector, matrix, data frame, timeSeries or zoo object of
189#' index, benchmark, portfolio, or secondary asset returns to compare against
190#' @param unbiased TRUE/FALSE whether to use a correction to have an unbiased
191#' estimator, default FALSE
192#' @param as.mat TRUE/FALSE whether to return the full moment matrix or only
193#' the vector with the unique elements (the latter is advised for speed), default
194#' TRUE
195#' @param \dots any other passthru parameters
196#' @author Kris Boudt, Peter Carl, Dries Cornilly, Brian Peterson
197#' @seealso \code{\link{BetaCoSkewness}} \cr \code{\link{BetaCoKurtosis}} \cr
198#' \code{\link{BetaCoMoments}} \cr \code{\link{ShrinkageMoments}} \cr \code{\link{EWMAMoments}}
199#' \cr \code{\link{StructuredMoments}} \cr \code{\link{MCA}}
200#' @references Boudt, Kris, Brian G. Peterson, and Christophe Croux. 2008.
201#' Estimation and Decomposition of Downside Risk for Portfolios with Non-Normal
202#' Returns. Journal of Risk. Winter.
203#'
204#' Boudt, Kris, Cornilly, Dries and Verdonck, Tim. 2017. A Coskewness Shrinkage
205#' Approach for Estimating the Skewness of Linear Combinations of Random Variables.
206#' Submitted. Available at SSRN: https://ssrn.com/abstract=2839781
207#'
208#' Martellini, L., & Ziemann, V. 2010. Improved estimates of higher-order
209#' comoments and implications for portfolio selection. Review of Financial
210#' Studies, 23(4), 1467-1502.
211#'
212#' Ranaldo, Angelo, and Laurent Favre Sr. 2005. How to Price Hedge Funds: From
213#' Two- to Four-Moment CAPM. SSRN eLibrary.
214#'
215#' Scott, Robert C., and Philip A. Horvath. 1980. On the Direction of
216#' Preference for Moments of Higher Order than the Variance. Journal of Finance
217#' 35(4):915-919.
218###keywords ts multivariate distribution models
219#' @examples
220#'
221#' data(managers)
222#' CoVariance(managers[, "HAM2", drop=FALSE], managers[, "SP500 TR", drop=FALSE])
223#' CoSkewness(managers[, "HAM2", drop=FALSE], managers[, "SP500 TR", drop=FALSE])
224#' CoKurtosis(managers[, "HAM2", drop=FALSE], managers[, "SP500 TR", drop=FALSE])
225#'
226#' @export CoVariance
227CoVariance<- function(Ra,Rb)
228{# @author Kris Boudt, Peter Carl
229    Ra= checkData(Ra)
230    Rb= checkData(Rb)
231
232    Ra.ncols = NCOL(Ra)
233    Rb.ncols = NCOL(Rb)
234
235    pairs = expand.grid(1:Ra.ncols, 1:Rb.ncols)
236
237    covar <-function (Ra, Rb)
238    {
239        R = na.omit(cbind(Ra, Rb)) # remove NA's
240        return(centeredcomoment(R[,1],R[,2],p1=1,p2=1,normalize=FALSE))
241    }
242
243    result = apply(pairs, 1, FUN = function(n, Ra, Rb) covar(Ra[,n[1]], Rb[,n[2]]), Ra = Ra, Rb = Rb)
244
245    if(length(result) ==1)
246        return(result)
247    else {
248        dim(result) = c(Ra.ncols, Rb.ncols)
249        colnames(result) = paste("Covariance:", colnames(Rb))
250        rownames(result) = colnames(Ra)
251        return(t(result))
252    }
253}
254
255#' Functions to calculate systematic or beta co-moments of return series
256#'
257#' calculate higher co-moment betas, or 'systematic' variance, skewness, and
258#' kurtosis
259#'
260#' The co-moments, including covariance, coskewness, and cokurtosis, do not
261#' allow the marginal impact of an asset on a portfolio to be directly
262#' measured.  Instead, Martellini and Zieman (2007) develop a framework that
263#' assesses the potential diversification of an asset relative to a portfolio.
264#' They use higher moment betas to estimate how much portfolio risk will be
265#' impacted by adding an asset, in terms of symmetric risk (i.e., volatility),
266#' in asymmetry risk (i.e., skewness), and extreme risks (i.e. kurtosis). That
267#' allows them to show that adding an asset to a portfolio (or benchmark) will
268#' reduce the portfolio's variance to be reduced if the second-order beta of
269#' the asset with respect to the portfolio is less than one.  They develop the
270#' same concepts for the third and fourth order moments.  The authors offer
271#' these higher moment betas as a measure of the diversification potential of
272#' an asset.
273#'
274#' Higher moment betas are defined as proportional to the derivative of the
275#' covariance, coskewness and cokurtosis of the second, third and fourth
276#' portfolio moment with respect to the portfolio weights. The beta co-variance
277#' is calculated as:
278#'
279#' \deqn{ BetaCoV(Ra,Rb) = \beta^{(2)}_{a,b} =
280#' \frac{CoV(R_a,R_b)}{\mu^{(2)}(R_b)} }{BetaCoV(Ra,Rb) =
281#' CoV(Ra,Rb)/centeredmoment(Rb,2)}
282#'
283#' Beta co-skewness is given as:
284#'
285#' \deqn{ BetaCoS(Ra,Rb) = \beta^{(3)}_{a,b}= \frac{CoS(R_a,R_b)}{\mu^{(3)}(R_b)} }{BetaCoS(Ra,Rb) =
286#' CoS(Ra,Rb)/centeredmoment(Rb,3)}
287#'
288#' Beta co-kurtosis is:
289#'
290#' \deqn{ BetaCoK(Ra,Rb)=\beta^{(4)}_{a,b}
291#' = \frac{CoK(R_a,R_b)}{\mu^{(4)}(R_b)} }{BetaCoK(Ra,Rb) =
292#' CoK(Ra,Rb)/centeredmoment(Rb,4)}
293#'
294#' where the \eqn{n}-th centered moment is
295#' calculated as
296#'
297#' \deqn{ \mu^{(n)}(R) = E\lbrack(R-E(R))^n\rbrack }{moment^n(R) = E[R-E(R)^n]}
298#'
299#' A beta is greater than one indicates that no diversification benefits should
300#' be expected from the introduction of that asset into the portfolio.
301#' Conversely, a beta that is less than one indicates that adding the new asset
302#' should reduce the resulting portfolio's volatility and kurtosis, and to an
303#' increase in skewness. More specifically, the lower the beta the higher the
304#' diversification effect on normal risk (i.e. volatility). Similarly, since
305#' extreme risks are generally characterised by negative skewness and positive
306#' kurtosis, the lower the beta, the higher the diversification effect on
307#' extreme risks (as reflected in Modified Value-at-Risk or ER).
308#'
309#' The addition of a small fraction of a new asset to a portfolio leads to a
310#' decrease in the portfolio's second moment (respectively, an increase in the
311#' portfolio's third moment and a decrease in the portfolio's fourth moment) if
312#' and only if the second moment (respectively, the third moment and fourth
313#' moment) beta is less than one (see Martellini and Ziemann (2007) for more
314#' details).
315#'
316#' For skewness, the interpretation is slightly more involved.  If the skewness
317#' of the portfolio is negative, we would expect an increase in portfolio
318#' skewness when the third moment beta is lower than one. When the skewness of
319#' the portfolio is positive, then the condition is that the third moment beta
320#' is greater than, as opposed to lower than, one.
321#'
322#' Because the interpretation of beta coskewness is made difficult by the need
323#' to condition on it's skewness, we deviate from the more widely used measure
324#' slightly.  To make the interpretation consistent across all three measures,
325#' the beta coskewness function tests the skewness and multiplies the result by
326#' the sign of the skewness.  That allows an analyst to review the metric and
327#' interpret it without needing additional information.  To use the more widely
328#' used metric, simply set the parameter \code{test = FALSE}.
329#' @name BetaCoMoments
330#' @concept beta co-moments
331#' @concept moments
332#' @aliases BetaCoMoments BetaCoVariance BetaCoSkewness BetaCoKurtosis SystematicSkewness SystematicKurtosis
333#' @param Ra an xts, vector, matrix, data frame, timeSeries or zoo object of
334#' asset returns
335#' @param Rb an xts, vector, matrix, data frame, timeSeries or zoo object of
336#' index, benchmark, or secondary asset returns to compare against
337#' @param test condition not implemented yet
338#' @author Kris Boudt, Peter Carl, Brian Peterson
339#' @seealso \code{\link{CoMoments}}
340#' @references
341#'
342#' Boudt, Kris, Brian G. Peterson, and Christophe Croux. 2008. Estimation and
343#' Decomposition of Downside Risk for Portfolios with Non-Normal Returns.
344#' Journal of Risk. Winter.
345#'
346#' Martellini, Lionel, and Volker Ziemann. 2007. Improved Forecasts of
347#' Higher-Order Comoments and Implications for Portfolio Selection. EDHEC Risk
348#' and Asset Management Research Centre working paper.
349###keywords ts multivariate distribution models
350#' @examples
351#'
352#' data(managers)
353#'
354#' BetaCoVariance(managers[, "HAM2", drop=FALSE], managers[, "SP500 TR", drop=FALSE])
355#' BetaCoSkewness(managers[, "HAM2", drop=FALSE], managers[, "SP500 TR", drop=FALSE])
356#' BetaCoKurtosis(managers[, "HAM2", drop=FALSE], managers[, "SP500 TR", drop=FALSE])
357#' BetaCoKurtosis(managers[,1:6], managers[,8,drop=FALSE])
358#' BetaCoKurtosis(managers[,1:6], managers[,8:7])
359#'
360#' @export  BetaCoVariance
361BetaCoVariance <- function(Ra,Rb)
362{# @author Kris Boudt, Peter Carl
363    Ra= checkData(Ra)
364    Rb= checkData(Rb)
365
366    Ra.ncols = NCOL(Ra)
367    Rb.ncols = NCOL(Rb)
368
369    pairs = expand.grid(1:Ra.ncols, 1:Rb.ncols)
370
371    bcovar <-function (Ra, Rb)
372    {
373        R = na.omit(cbind(Ra, Rb)) # remove NA's
374        return(centeredcomoment(R[,1],R[,2],p1=1,p2=1,normalize=TRUE))
375    }
376
377    result = apply(pairs, 1, FUN = function(n, Ra, Rb) bcovar(Ra[,n[1]], Rb[,n[2]]), Ra = Ra, Rb = Rb)
378
379    if(length(result) ==1)
380        return(result)
381    else {
382        dim(result) = c(Ra.ncols, Rb.ncols)
383        colnames(result) = paste("Beta Covariance:", colnames(Rb))
384        rownames(result) = colnames(Ra)
385        return(t(result))
386    }
387}
388
389#' @rdname CoMoments
390#' @export
391CoSkewness <- function(Ra,Rb)
392{# @author Kris Boudt, Peter Carl
393    Ra= checkData(Ra)
394    Rb= checkData(Rb)
395
396    Ra.ncols = NCOL(Ra)
397    Rb.ncols = NCOL(Rb)
398
399    pairs = expand.grid(1:Ra.ncols, 1:Rb.ncols)
400
401    cosk <-function (Ra, Rb)
402    {
403        R = na.omit(cbind(Ra, Rb)) # remove NA's
404        return(centeredcomoment(R[,1],R[,2],p1=1,p2=2,normalize=FALSE))
405    }
406
407    result = apply(pairs, 1, FUN = function(n, Ra, Rb) cosk(Ra[,n[1]], Rb[,n[2]]), Ra = Ra, Rb = Rb)
408
409    if(length(result) ==1)
410        return(result)
411    else {
412        dim(result) = c(Ra.ncols, Rb.ncols)
413        colnames(result) = paste("Coskewness:", colnames(Rb))
414        rownames(result) = colnames(Ra)
415        return(t(result))
416    }
417}
418
419#' @rdname BetaCoMoments
420#' @export
421BetaCoSkewness <- function(Ra, Rb, test=FALSE)
422{# @author Kris Boudt, Peter Carl
423    Ra= checkData(Ra)
424    Rb= checkData(Rb)
425
426    Ra.ncols = NCOL(Ra)
427    Rb.ncols = NCOL(Rb)
428
429    pairs = expand.grid(1:Ra.ncols, 1:Rb.ncols)
430
431    bcosk <-function (Ra, Rb)
432    {
433        R = na.omit(cbind(Ra, Rb)) # remove NA's
434        skew = skewness(Rb)
435        # Kris notes: the output should be conditional on the value of the market skewness.
436        if(skew > -0.05 && skew < 0.05 ){
437            warning("skewness is close to zero. The classical definition of the coskewness statistic is not applicable and one should normalize using the comoment without standardization.")
438        }
439        if(test==TRUE){
440#             if(skew < 0)
441#                 multiplier = -1
442#             else
443#                 multiplier = 1
444            stop("Not implemented yet")
445        }
446        else
447            multiplier = 1
448
449        return(multiplier * centeredcomoment(R[,1],R[,2],p1=1,p2=2,normalize=TRUE))
450    }
451
452    result = apply(pairs, 1, FUN = function(n, Ra, Rb) bcosk(Ra[,n[1]], Rb[,n[2]]), Ra = Ra, Rb = Rb)
453
454    if(length(result) ==1)
455        return(result)
456    else {
457        dim(result) = c(Ra.ncols, Rb.ncols)
458        colnames(result) = paste("Beta Coskewness:", colnames(Rb))
459        rownames(result) = colnames(Ra)
460        return(t(result))
461    }
462}
463
464#' @rdname CoMoments
465#' @export
466CoKurtosis <- function(Ra,Rb)
467{# @author Kris Boudt, Peter Carl
468    Ra= checkData(Ra)
469    Rb= checkData(Rb)
470
471    Ra.ncols = NCOL(Ra)
472    Rb.ncols = NCOL(Rb)
473
474    pairs = expand.grid(1:Ra.ncols, 1:Rb.ncols)
475
476    cokurt <-function (Ra, Rb)
477    {
478        R = na.omit(cbind(Ra, Rb)) # remove NA's
479        return(centeredcomoment(R[,1],R[,2],p1=1,p2=3,normalize=FALSE))
480    }
481
482    result = apply(pairs, 1, FUN = function(n, Ra, Rb) cokurt(Ra[,n[1]], Rb[,n[2]]), Ra = Ra, Rb = Rb)
483
484    if(length(result) ==1)
485        return(result)
486    else {
487        dim(result) = c(Ra.ncols, Rb.ncols)
488        colnames(result) = paste("Cokurtosis:", colnames(Rb))
489        rownames(result) = colnames(Ra)
490        return(t(result))
491    }
492}
493
494#' @rdname BetaCoMoments
495#' @export
496BetaCoKurtosis <- function(Ra,Rb)
497{# @author Kris Boudt, Peter Carl
498    Ra= checkData(Ra)
499    Rb= checkData(Rb)
500
501    Ra.ncols = NCOL(Ra)
502    Rb.ncols = NCOL(Rb)
503
504    pairs = expand.grid(1:Ra.ncols, 1:Rb.ncols)
505
506    bcosk <-function (Ra, Rb)
507    {
508        R = na.omit(cbind(Ra, Rb)) # remove NA's
509        return(centeredcomoment(R[,1],R[,2],p1=1,p2=3,normalize=TRUE))
510    }
511
512    result = apply(pairs, 1, FUN = function(n, Ra, Rb) bcosk(Ra[,n[1]], Rb[,n[2]]), Ra = Ra, Rb = Rb)
513
514    if(length(result) ==1)
515        return(result)
516    else {
517        dim(result) = c(Ra.ncols, Rb.ncols)
518        colnames(result) = paste("Beta Cokurtosis:", colnames(Rb))
519        rownames(result) = colnames(Ra)
520        return(t(result))
521    }
522}
523
524
525###############################################################################
526# R (http://r-project.org/) Econometrics for Performance and Risk Analysis
527#
528# Copyright (c) 2004-2020 Peter Carl and Brian G. Peterson
529#
530# This R package is distributed under the terms of the GNU Public License (GPL)
531# for full details see the file COPYING
532#
533# $Id$
534#
535###############################################################################
536