Multidimensional normal in R
I am trying to generate multivariate normal random numbers using a command rmsn
from the package sn
in R. I would like, ideally, to be able to get three columns of numbers with variances and covariance specified, with one column highly skewed. But I try my best to achieve both goals.
The post in malformed normal distribution was related and helpful (and source some code below), but did not fully clarify the issue for me.
I have tried:
a <- c(5, 0, 0) # set shape parameter
s <- diag(3) # create variance-covariance matrix
w <- sqrt(1/(1-((2*(a^2)/(1 + a^2))/pi))) # determine scale parameter to get sd of 1
xi <- w*a/sqrt(1 + a^2)*sqrt(2/pi) # determine location parameter to get mean of 0
apply(rmsn(n=1000, xi=c(xi), Omega=s, alpha=a), 2, sd)
colMeans(rmsn(n=1000, xi=c(xi), Omega=s, alpha=a))
Column and SD values are correct for the second and third columns (which are not skewed), but not the first (which does). Can anyone clarify where my code above is, or my thinking went wrong? I may not understand how to use rmsn
or output. Any help would be appreciated.
source to share
The location is not average (unless there is no bias). From the documentation:
Note that the location vector "xi" does not mean the mean of the distribution vector (which may in fact not even exist if 'df <= 1), and similarly "Omega" is not a covariance matrix distribution
And you can replace Omega=s
with Omega=w
. And that should be a variance matrix: there shouldn't be a square root.
source to share