I think it is time to start paying attention to the answers you receive: https://stats.stackexchange.com/a/54734
It has been largely discussed the validity of using a measure of skewness in multimodal distributions, since its interpretation becomes unclear. This is the case of finite mixtures. If your mixture looks (or is) unimodal, then you can use this value to understand a bit how asymmetric it is.
In addition, note that MLE is based on maximising the likelihood, not on matching moments.
Check the following example using a simulated sample from a skew-normal (the sample skewness and the skewness of the mixture are different).
library(sn)
library(moments)
library(mixtools)
set.seed(4321)
samp <- rsn(1000,0,1,10)
skewness(samp)
mix<-normalmixEM(samp,k=2,fast=TRUE)
# Sampling from a 2-gaussian mixture
gaussmix <- function(n,m1,m2,s1,s2,alpha) {
I <- runif(n)<alpha
rnorm(n,mean=ifelse(I,m1,m2),sd=ifelse(I,s1,s2))
}
# A simulated sample
samp2 <- gaussmix(100000,mix$mu[1],mix$mu[2],mix$sigma[1],mix$sigma[2],mix$lambda[1])
skewness(samp2)
Update
The measure of skewness based on the third moment is driven by the tails, rather than the level of asymmetry of the distribution. For example, a Student-$t$ distribution with $\nu<3$ degrees of freedom has undefined skewness in spite of being clearly symmetric. For this reason, it is preferred to use a quantile-based measure of skewness. An example of this of this is the AG measure of skewness which can be defined for any unimodal distribution $F$ as
$$AG=1-2F(\mbox{mode}).$$
This can be estimated by using a nonparametric estimator of $F$, such as a kernel estimator as follows
rm(list=ls())
library(sn)
library(moments)
library(mixtools)
# AG measure of skewness: insert the data and the interval c(minim,maxim) where the mode is located
AG <- function(data,minim,maxim){
n = length(data)
hb = (4*sqrt(var(data))^5/(3*n))^(1/5)
kn = function(x){
k = r = length(x)
for(i in 1:k) r[i] = mean(dnorm((x[i]-data)/hb))/hb
return(r)
}
mode = optimise(f = kn, interval = c(minim,maxim),maximum = TRUE)$maximum
KN = function(x){
k = r = length(x)
for(i in 1:k) r[i] = mean(pnorm((x[i]-data)/hb))
return(r)
}
return(1-2*KN(mode))
}
# Simulated sample
set.seed(707)
samp <- rsn(1000,0,1,10)
# AG and moment-based skewness measure of the sample
AG(samp,0,2)
skewness(samp)
mix<-normalmixEM(samp,k=2,fast=TRUE)
# Sampling from a 2-gaussian mixture
gaussmix <- function(n,m1,m2,s1,s2,alpha) {
I <- runif(n)<alpha
rnorm(n,mean=ifelse(I,m1,m2),sd=ifelse(I,s1,s2))
}
# A simulated sample
mix<-normalmixEM(samp,k=2,fast=TRUE)
samp2 <- gaussmix(100000,mix$mu[1],mix$mu[2],mix$sigma[1],mix$sigma[2],mix$lambda[1])
# AG and moment-based skewness measure of the fitted mixture
AG(samp2,0,2)
skewness(samp2)
As you can see, the AG measure of skewness is similar for both distributions, as expected.
The moral of the story is: The moment-based measure of skewness is NOT the best choice.
The same happens if this measure is applied to your data:
library(sn)
library(moments)
library(mixtools)
# AG measure of skewness: insert the data and the interval c(minim,maxim) where the mode is located
AG <- function(data,minim,maxim){
n = length(data)
hb = (4*sqrt(var(data))^5/(3*n))^(1/5)
kn = function(x){
k = r = length(x)
for(i in 1:k) r[i] = mean(dnorm((x[i]-data)/hb))/hb
return(r)
}
mode = optimise(f = kn, interval = c(minim,maxim),maximum = TRUE)$maximum
KN = function(x){
k = r = length(x)
for(i in 1:k) r[i] = mean(pnorm((x[i]-data)/hb))
return(r)
}
return(1-2*KN(mode))
}
# AG and moment-based skewness measure of the sample
AG(dat,-.1,.1)
skewness(dat)
# Sampling from a 2-gaussian mixture
gaussmix <- function(n,m1,m2,s1,s2,alpha) {
I <- runif(n)<alpha
rnorm(n,mean=ifelse(I,m1,m2),sd=ifelse(I,s1,s2))
}
# A simulated sample
mix<-normalmixEM(dat,k=2,fast=TRUE)
samp2 <- gaussmix(100000,mix$mu[1],mix$mu[2],mix$sigma[1],mix$sigma[2],mix$lambda[1])
# AG and moment-based skewness measure of the fitted mixture
AG(samp2,-.1,.1)
skewness(samp2)