Estimating parameters for univariate skew t

Question

How can I solve the MLE for the skew-t distribution via EM? I am comfortable with the EM methods for t, so could someone show it for the skew-t?

Answer 1

I will base my answer on the skew-$t$-normal distribution of Gomez et al (2007) (see the end of the answer for full references) for reasons discussed shortly. Let $z \sim StN(\mu,\sigma^2,\lambda,\nu)$ denote that $z$ follows a skew-$t$-normal distribution with location parameter $\mu$, scale parameter $\sigma^2$, skewness parameter $\lambda$ and degrees of freedom parameter $\nu$. Then the density of $z$ is given by:

$$ f(z|\mu,\sigma^2,\lambda,\nu) = 2t_\nu(z|\mu,\sigma^2)\Phi\left(\frac{z-\mu}{\sigma}\lambda\right) $$

where $t_\nu(z|\mu,\sigma^2)$ denotes the $t$ density with $\nu$ degrees of freedom, location parameter $\mu$ and scale parameter $\sigma^2$. $\Phi$ denotes the Gaussian cumulative distribution function.

This form is identical to the skew-$t$ distribution proposed by Azzalini and Capitanio (2003) except that the Gaussian CDF controlling the skewness of the density is replaced by the Student $t$ CDF, which depends not just on $\lambda$ but also the degrees of freedom parameter $\nu$. The form given here has the advantage of being computationally simpler to evaluate and I think it's better to decouple the skewness from the degrees of freedom parameter.

To proceed with the EM algorithm, note the following hierarchical decomposition:

$$ \tau \sim \Gamma\left(\frac{\nu}{2},\frac{\nu}{2}\right)\\ \gamma|\tau \sim TN\left(0,\frac{\tau+\lambda^2}{\tau};(0,\infty)\right)\\ z|\tau,\gamma \sim N\left(\mu + \frac{\sigma\lambda}{\tau + \lambda^2}\gamma,\frac{\sigma^2}{\tau + \lambda^2}\right) $$

Where $TN(\mu,\sigma^2;(a,b))$ denotes the truncated normal distribution lying within the interval $(a,b)$. So the idea is to treat $\tau$ and $\gamma$ as latent variables. These are the maximum likelihood estimators of the complete log likelihood:

$$ \widehat{\mu} = \frac{\sum_{i=1}^n[\tau_iz_i] + \lambda^2\sum_{i=1}^nz_i-\lambda\sigma\sum_{i=1}^n\gamma_i}{n\lambda^2+\sum_{i=1}^n\tau_i}\\ \widehat{\sigma^2}=\frac{1}{n}\sum_{i=1}^n\tau_i(z_i-\mu)^2\\ \widehat{\delta}=\frac{\sum_{i=1}^n(z_i-\mu)\gamma_i}{\sum_{i=1}^n(z_i-\mu)^2}\\ \widehat{\lambda}=\widehat{\delta}\widehat{\sigma} $$

there is no closed form estimator $\nu$ so instead you will have to numerically maximize:

$$ \frac{\nu}{2}\log\frac{\nu}{2} - \log\Gamma(\frac{\nu}{2}) + \frac{\nu-1}{2}\log\tau_i-\frac{\nu}{2}\tau_i $$

These are the conditional expectations of the latent variables:

$$ E[\gamma_i|z_i] = (z_i - \mu)\frac{\lambda}{\sigma} + \frac{\phi\left((z_i - \mu)\frac{\lambda}{\sigma}\right)}{\Phi\left((z_i - \mu)\frac{\lambda}{\sigma}\right)}\\ E[\tau_i|z_i] = \frac{\nu+1}{\nu + \frac{(z_i - \mu)^2}{\sigma^2}} $$

I apologize in advance for any potential mistakes in the formulae and suggest you verify them by deriving them yourself.

Also, as I said in my comment, I believe you are better off foregoing the EM algorithm and simply maximizing the observed likelihood with, for example, Nelder-Mead.

---

References:

Azzalini, Adelchi, and Antonella Capitanio. "Distributions generated by perturbation of symmetry with emphasis on a multivariate skew t‐distribution." Journal of the Royal Statistical Society: Series B (Statistical Methodology) 65, no. 2 (2003): 367-389.

Gómez, Héctor W., Osvaldo Venegas, and Heleno Bolfarine. "Skew‐symmetric distributions generated by the distribution function of the normal distribution." Environmetrics 18, no. 4 (2007): 395-407.