The effect of the shape (skewness) parameter in skew-symmetric models, Part II
Francisco Javier Rubio
August 22, 2025
The effect of the parameter \(\lambda\) in skew-symmetric models
The parameter \(\lambda\) in skew-symmetric models of the type (Azzalini 1985)
\[s(x\mid \mu, \sigma ,\lambda) = 2 f\left(\dfrac{x-\mu}{\sigma}\right)G\left(\lambda\dfrac{x-\mu}{\sigma}\right),\] controls the asymmetry of the probability density function (pdf), however, it has a different role for different choices of \(f\) and \(G\). In the case \(f=\phi\) and \(G=\Phi\), we obtain the skew-normal distribution (Azzalini 1985). In this case, Rubio and Genton (2016) pointed out that the skew normal distribution is virtually symmetric for \(\vert \lambda \vert < 1.25\). In this region, the shape parameter \(\lambda\) affects mainly the location of the mode and the scale of the density function, but the effect on asymmetry is negligible.
I previously analysed the effect of \(\lambda\) in terms of the minimum total variation distance, see: The effect of the shape (skewness) parameter in skew-symmetric models. In this document, we use instead the recently proposed discrepancy measure (Wagenmakers and Grasman 2025):
\[ d(f \mid\mid g) = \int_{\mathcal{X}} \dfrac{f(x)}{f(x)+g(x)} f(x)dx. \] Wagenmakers and Grasman (2025) showed that \(d(f \mid\mid g) = d(g \mid\mid f)\).
In the context of skew-symmetric models, we might be interested on identifying whether the skewness parameter has a non-negligible effect on the pdf or not, even if it is theoretically non-zero (Dette, Ley, and Rubio 2018). Thus, we can think of comparing a skew-symmetric model with parameters \((0,1,\lambda)\) and the corresponding ``closest’’ nested symmetric density \(f\) with location \(\mu\) and scale \(\sigma\). In the skew-normal case, we would like to compare the skew normal density with parameters \((0,1,\lambda)\) and the closest Normal density with location \(\mu\) and scale \(\sigma\). This is, we would like to find the closest normal model to the skew normal model. For this purpose, we can employ the discrepancy measure above and define the minimum discrepancy (translated to \(0\))
\[\begin{align} m(\lambda) &= \min_{\mu,\sigma} d\left(s(\mu,\sigma,0) \mid\mid s(0,1,\lambda)\right) - \dfrac{1}{2} \\ &= \min_{\mu,\sigma} \int_{-\infty}^{\infty} \dfrac{\dfrac{1}{\sigma^2}f\left(\dfrac{x-\mu}{\sigma}\right)^2}{\dfrac{1}{\sigma}f\left(\dfrac{x-\mu}{\sigma}\right)+2f(x)G(\lambda x)} dx - \dfrac{1}{2}. \end{align}\]
The following R code shows the minimum discrepancy between several skew-symmetric models and the corresponding nested symmetric models, as a function of \(\lambda\). We can see that the effect of the skewness parameter differ for different skew-symmetric distributions. This also provides an intuition about the nature of some inferential problems which appear when estimating \(\lambda\), such as those in the skew-normal distribution (Pewsey 2000). That is, in a neighbourhood of \(\lambda=0\), it is difficult to distinguish the skew-symmetric model from the nested symmetric model (Rubio and Genton 2016).
The Skew Normal distribution
##################################################################################################
# Skew Normal
##################################################################################################
rm(list=ls())
library(sn)## Loading required package: stats4##
## Attaching package: 'sn'## The following object is masked from 'package:stats':
##
## sd# The Skew Normal is virtually symmetric for |lambda| < 1.25
tempf1 <- Vectorize(function(x) dsn(x,0,1,1))
tempf2 <- Vectorize(function(x) dnorm(x,0.534,exp(-0.2)))
curve(tempf1,-5,5, ylab = "Density")
curve(tempf2,-5,5,add=T,lty=2)
# Minimum discrepancy between the SN distribution and the closest Normal model
discrepancy_min <- Vectorize(function(lambda){
disc <- function(par){
tempf <- Vectorize(function(x) dnorm(x,par[1],exp(par[2]))*dnorm(x,par[1],exp(par[2]))/(dnorm(x,par[1],exp(par[2]))+2*dnorm(x)*pnorm(lambda*x)) )
int <- integrate(tempf,-5,5)$value
return(int)
}
val <- optim(c(0,0),disc)$value
return(abs(val-0.5))
})
discrepancy_min(0)## [1] 2.866516e-07discrepancy_min(1)## [1] 0.0003774259discrepancy_min(1000)## [1] 0.04169079curve(discrepancy_min, -10, 10, n = 300, lwd = 2, ylab = "m", xlab = expression(lambda))
abline(h=0.005, lwd = 2, col = "red")
The Skew Logistic distribution
##################################################################################################
# Skew Logistic
##################################################################################################
rm(list=ls())
# Minimum discrepancy between the skew-logistic distribution and the closest logistic model
discrepancy_min <- Vectorize(function(lambda){
disc <- function(par){
tempf <- Vectorize(function(x) dlogis(x,par[1],exp(par[2]))*dlogis(x,par[1],exp(par[2]))/(dlogis(x,par[1],exp(par[2]))+2*dlogis(x)*plogis(lambda*x)) )
int <- integrate(tempf,-15,15)$value
return(int)
}
val <- optim(c(0,0),disc)$value
return(abs(val-0.5))
})
discrepancy_min(0)## [1] 3.058412e-07discrepancy_min(1)## [1] 0.003111176discrepancy_min(1000)## [1] 0.05050229curve(discrepancy_min, -10, 10, n = 300, lwd = 2, ylab = "m", xlab = expression(lambda))
abline(h=0.005, lwd = 2, col = "red")
The Skew t(2) distribution
##################################################################################################
# Skew Logistic
##################################################################################################
rm(list=ls())
dts <- function(x, mu, sigma, nu) dt( (x - mu)/sigma, df = nu)/sigma
# Minimum discrepancy between the skew-t(2) distribution and the closest t(2) model
discrepancy_min <- Vectorize(function(lambda){
disc <- function(par){
tempf <- Vectorize(function(x) dts(x,par[1],exp(par[2]),2)*dts(x,par[1],exp(par[2]),2)/(dts(x,par[1],exp(par[2]),2)+2*dt(x,df=2)*pt(lambda*x,df=2)) )
int <- integrate(tempf,-Inf,Inf)$value
return(int)
}
val <- optim(c(0,0),disc)$value
return(abs(val-0.5))
})
discrepancy_min(0)## [1] 1.110223e-16discrepancy_min(1)## [1] 0.01294834discrepancy_min(1000)## [1] 0.06603632curve(discrepancy_min, -10, 10, n = 300, lwd = 2, ylab = "m", xlab = expression(lambda))
abline(h=0.005, lwd = 2, col = "red")