Expected value of the SHASH distribution?

Question

The sinh-arcsinh (SHASH) distribution has a pdf as follows:

$f(x) = {\delta cosh(\omega)\over \sqrt{1+({x-\theta \over \sigma})^2}}\phi[sinh(\omega)]$ where $\omega=\gamma+\delta sinh^{-1}({x-\theta \over \sigma})$

and $\phi()$ is the standard normal pdf.

Note: $-\infty<\gamma,x,\theta<\infty; 0<\delta, \sigma$

Does anyone know what the expected value and variance of it is? I can see it is closely related to the Johnson SU distribution which has the mean, variance, median, and quantiles analytically defined. It also has the Normal distribution as a special case when $\gamma = 0$ and $ \delta = 1$. I know some distributions don't have these available in closed form analytic expressions....but I don't know if that's the case here. If available, it would be very helpful if someone knows what they are!

Update: It appears the moments are published, but with the caveat that they depend on the modified Bessel function of the second kind ($P_\nu)$ below. E.g.

$E(X_{\gamma, \delta})=-sinh(\gamma/ \delta)P_{1/\delta}$

$Var(X_{\gamma, \delta})={1 \over 2}(cosh({2\gamma \over \delta})P_{2/\delta} - 1) - \mu_{\gamma, \delta}^2$

However -- this doesn't take into account the center, scale in the full version ($\theta, \sigma$). Can someone help modify the mean and var with center and scale included?

Answer 1

Indeed, the moments of this distribution are already calculated in the original paper

https://www.jstor.org/stable/27798865

For the case with location and scale parameters, you just need to use the usual properties of the location-scale family. Let $Z = \mu + \sigma X$, where $X \sim SHASH(\gamma,\delta)$, then

$$E[Z] = \mu + \sigma E[X].$$

Also,

$$Var[Z] = \sigma^2 Var[X].$$

[Here] is an implementation of this distribution in R. Using their command rsas, you can simulate from a distribution with specific parameter values and approximate its mean and variance using mean() and var().

sim <- rsas(1e5,3,2,-1,1.5)
mean(sim)
var(sim)

Alternatively, you can use their function dsas to approximate the mean using numerical integration:

mean_sas <- function(mu,sigma,epsilon,delta){
  tempf <- Vectorize(function(x) x*dsas(x,mu,sigma,epsilon,delta))
  val <- integrate(tempf,-Inf,Inf)$value
  return(val)
}

mean_sas(3,2,-1,1.5) 

Similarly for the variance:

var_sas<- function(mu,sigma,epsilon,delta){
  tempf1 <- Vectorize(function(x) x*dsas(x,mu,sigma,epsilon,delta))
  tempf2 <- Vectorize(function(x) x^2*dsas(x,mu,sigma,epsilon,delta))
  val <- integrate(tempf2,-Inf,Inf)$value - (integrate(tempf1,-Inf,Inf)$value)^2
  return(val)
}


var_sas(3,2,-1,1.5)

Finally, the same ideas can be applied to [their alternative version] of the SHASH distribution.