# Calculation of Higher-Order Cross-moments

How can I calculate standardized central cross-moments for 2 time-series?

The 4th-order standardized central moment, kurtosis, is;

 kurt = mean((x-mean(x)).^4)./mean(x.^2).^2;


It should be something like the kurt, for x and y.

I also want to know the 3rd order cross-moment (something like skewness).

I am trying to write a code in MATLAB.

• Is this an attempt at LaTeX? Or are you trying to enter code? If the latter, what language are you using? Oct 18, 2014 at 16:16
• @wolfies I'm trying to write a code in Matlab.
– Rash
Oct 18, 2014 at 16:22
• You included the tag for [discrete-data]. Are your data discrete? (You don't mention anything about that in the question.) Is that important for your situation in some way? Oct 18, 2014 at 17:21
• @gung, My mistake, My data is discrete, in matrix x, and I want to write the codes for calculation.
– Rash
Oct 18, 2014 at 17:34

## 1 Answer

You are in essence looking for a multivariate measure of skew and kurtosis. There are many. I would start with the most establish ones, which are the multivariate skew and kurtosis measures of Mardia 1977 [0].

It seems to me you are more asking for an implementation than about the measures themselves. I don t know of any matlab implementation, but the R code below (from the R library psych) should be fairly easy to translate in matlab:

mardia <- function(x, na.rm=TRUE, plot=TRUE) {
cl <- match.call()
x <- as.matrix(x)     # in case it was a dataframe
if(na.rm) x <- na.omit(x)
n <- dim(x)[1]
p <- dim(x)[2]
x <- scale(x,scale=FALSE)  # zero center
S <- cov(x)
S.inv <- solve(S)
D <- x %*% S.inv %*% t(x)
b1p <- sum(D^3)/n^2
b2p <- tr(D^2)/n
chi.df <- p*(p+1)*(p+2)/6
k <- (p+1)*(n+1)*(n+3)/(n*((n+1)*(p+1) -6))

small.skew <- n*k*b1p/6
M.skew <- n*b1p/6
M.kurt <- (b2p - p * (p+2))*sqrt(n/(8*p*(p+2)))
p.skew <- 1-pchisq(M.skew,chi.df)
p.small <- 1 - pchisq(small.skew,chi.df)
p.kurt <- 2*(1- pnorm(abs(M.kurt)))
d = sqrt(diag(D))
if(plot) {qqnorm(d)
qqline(d)}
results <- list(n.obs=n, n.var=p, b1p=b1p,b2p= b2p, skew=M.skew,
small.skew=small.skew, p.skew=p.skew, p.small=p.small,
kurtosis=M.kurt, p.kurt=p.kurt, d=d, Call=cl)
class(results) <- c("psych", "mardia")
return(results)
}

• [0] K.V. Mardia (1970). Measures of multivariate skewness and kurtosis with applications. Biometrika, 57(3):pp. 519-30, 1970.
• I couldn't understand much about R code, but the paper was very useful, Thanks.
– Rash
Oct 18, 2014 at 16:50