MKL Summary Statistics computation of skewness and kurtosis seems to be broken when the standard deviation is more than four orders of magnitude smaller than the mean.

Can anyone tell me how to fix this, please?

For example, given normal random samples around 2, I should get skewness and excess kurtosis close to zero. However, at a standard deviation of 1.0e-4, Excess kurtosis jumps to -311. At 1.0e-5, skewness jumps to 11. At 1.0e-6, the absolute values of both skewness and kurtosis are gigantic.

In R, I get the correct values much longer.

I have

confirmed that I can workaround this problem by copying the data,

subtracting the mean and computing the moments with another MKL task.

But, man, this is ugly.

Here is my C++ code for computing summary stats with MKL

------------------------------------------------------------------

// Setup a Summary Statistic task pointing at the data

int status;

VSLSSTaskPtr task;

MKL_INT numRows = m_count;

MKL_INT numCols = 1;

MKL_INT xstorage = VSL_SS_MATRIX_STORAGE_ROWS;

int indices[1] = {1};

status = vsldSSNewTask( &task, &numCols, &numRows, &xstorage, &(m_data[0]), 0, indices );

// Setup 1st, 2nd and 3rd quartile computation

double quants[3];

double q_order[3] = {.25,.5,.75};

MKL_INT q_order_n = 3;

status = vsldSSEditQuantiles( task, &q_order_n, q_order, quants, NULL, NULL );

// Setup max and min computation

double min_est = m_data[0];

double max_est = m_data[0];

status = vsldSSEditTask(task, VSL_SS_ED_MIN, &min_est);

status = vsldSSEditTask(task, VSL_SS_ED_MAX, &max_est);

// Setup moment computation. Computing kurtosis requires all the other moments.

double mean_est, skewness_est, kurtosis_est;

double raw2mom, raw3mom, raw4mom;

double cen2mom, cen3mom, cen4mom;

status = vsldSSEditTask(task, VSL_SS_ED_SKEWNESS, &skewness_est);

status = vsldSSEditTask(task, VSL_SS_ED_KURTOSIS, &kurtosis_est);

status = vsldSSEditMoments( task, &mean_est, &raw2mom, &raw3mom, &raw4mom, &cen2mom, &cen3mom, &cen4mom );

// Request five number summary and moment computation

MKL_UINT64 mask = VSL_SS_MEAN | VSL_SS_KURTOSIS | VSL_SS_SKEWNESS |

VSL_SS_2R_MOM | VSL_SS_3R_MOM | VSL_SS_4R_MOM |

VSL_SS_2C_MOM | VSL_SS_3C_MOM | VSL_SS_4C_MOM |

VSL_SS_QUANTS | VSL_SS_MIN | VSL_SS_MAX;

status = vsldSSCompute(task, mask, VSL_SS_METHOD_FAST);

--------------------------------------------------------------------------

Here is the R code that I used to generate the test data as well and check skewness and kurtosis:

-----------------------------------------------------------------------

library(moments)

setwd("C:/Users/stennican/Documents/Rcode")

N = 500

normDist = rnorm(N)

trial = seq(1,N)

twos = rep(2,N)

twoSd1 = twos + 1e-1*normDist

twoSd2 = twos + 1e-2*normDist

twoSd3 = twos + 1e-3*normDist

twoSd4 = twos + 1e-4*normDist

twoSd5 = twos + 1e-5*normDist

twoSd7 = twos + 1e-7*normDist

twoSd10 = twos + 1e-10*normDist

twoSd13 = twos + 1e-13*normDist

zeros = rep(0,N)

zeroSd1 = 1e-1*normDist

zeroSd2 = 1e-2*normDist

zeroSd3 = 1e-3*normDist

zeroSd4 = 1e-4*normDist

zeroSd5 = 1e-5*normDist

zeroSd7 = 1e-7*normDist

zeroSd10 = 1e-10*normDist

zeroSd13 = 1e-13*normDist

out = data.frame(trial,zeros,twos,

zeroSd1,zeroSd2,zeroSd3,zeroSd4,zeroSd5,

zeroSd7,zeroSd10,zeroSd13,

twoSd1,twoSd2,twoSd3,twoSd4,twoSd5,

twoSd7,twoSd10,twoSd13)

apply(out,2,sd)

apply(out,2,skewness)

apply(out,2,kurtosis)

options(digits = 22)

write.csv(out,"testB42382.csv",row.names=FALSE)