In short, I would like to find the parameters of the skew normal so that given first three moments are matched. Find location, scale and shape of skew normal, so that mean = 0, variance = 1 and skewness = arbitrary value.
I would like to simulate from a skew-normal with the following properties
$E[Y] = 0 \\ E[Y^2] = 1\\ E[Y^3] = a \qquad a \quad \text{being any (possible) real number}$
I would say this is the same as a r.v. with mean = 0, variance = 1 and skewness = a. I am trying to use a skew normal Source1, wikipedia
There are formulas for the mean, variance and skewness. My approach was solving this with nleqlsv in R, but my optimization stops far away from a good result.
This is my whole R Code
library(sn)
sn_mean <- function(loc, scale, shape) {
delta <- shape / sqrt(1 + shape^2)
return ( loc + scale * delta * sqrt(2/ pi) )
}
sn_var <- function(loc, scale, shape) {
delta <- shape / sqrt(1 + shape^2)
return ( scale^2 * (1 - 2 * delta^2 / pi) )
}
sn_skew <- function(loc, scale, shape){
delta <- shape / sqrt(1 + shape^2)
return ( (4-pi) / 2 ) * ( ( delta * sqrt(2/pi))^3 / (1 - 2* delta^2 / pi)^(3/2) )
}
skew <- function(X){
return (sum(1/length(X) * ((X-mean(X))/sqrt(var(X)))^3))
}
objfunc <- function(x) {
shape <- x[1]
scale <- x[2]
loc <- x[3]
y <- numeric(3)
y[1] <- 0 - sn_mean(loc, scale, shape)
y[2] <- 1 - sn_var(loc, scale, shape)
y[3] <- 1 - sn_skew(loc, scale, shape)
return(y)
}
a <- nleqslv(c(100,100,100), objfunc, control=list(allowSingular=TRUE))
- Is there a better approach to get a skew normal with specific mean, var and skew?
- How can I improve my optimization to yield a better result?