Good morning.
I am writing a thesis about testing multivariate normality. I would like to do a comparison of power of some tests against given alternatives based on Monte Carlo simulations.
I have a following problem: Some multivariate normality tests are based on multivariate skewness $\beta_{1,k}$ and multivariate kurtosis $\beta_{2,k}$ defined by Mardia this way:
\begin{equation*} \beta_{1,k}(\mathbf{X}) := E [ \left\{ (\mathbf{X} - \boldsymbol \mu)^{T} \boldsymbol \Sigma^{-1} (\mathbf{Y} - \boldsymbol \mu) \right\}^{3} ], \end{equation*}
\begin{equation*} \beta_{2,k}(\mathbf{X}) := E [ \left\{ (\mathbf{X} - \boldsymbol \mu)^{T} \boldsymbol \Sigma^{-1} (\mathbf{X} - \boldsymbol \mu) \right\}^{2} ] \end{equation*} for iid. random vectors $\mathbf{X}$ and $\mathbf{Y}$. It is widely known that for $k$-dimensional normal distribution holds $\beta_{1,k} = 0$ and $\beta_{2,k} = k(k + 2)$.
There are some studies that compare power of normality test against alternatives with $\beta_{1,k} = 0$ but $\beta_{2,k} \neq k(k + 2)$ (such as Students $t$ distribution). I would like to do it in the opposite direction and compare power of my normality tests against alternative with $\beta_{1,k} \neq 0$ but $\beta_{2,k} = k(k + 2)$. My problem is that I am unable to find such a distribution.
Could anybody tell me, how to generate sample from multivariate distribution with $\beta_{1,k} \neq 0$ but $\beta_{2,k} = k(k + 2)$ in $R$ programme? I will be grateful for all useful ideas.
I tried to solve my problem by using bivariate Skew-Normal Distribution from package $sn$ but it is probably not possible for this distribution (with parametrization from $sn$ package) to have $\beta_{1,k} \neq 0$ but $\beta_{2,k} = k(k + 2)$.
Thank You!