Skewness and kurtosis unbiased by
Gaussian uncertainties
Abstract
Noise is an unavoidable part of most measurements which can hinder a correct interpretation of the data. Uncertainties propagate in the data analysis and can lead to biased results even in basic descriptive statistics such as the central moments and cumulants. Expressions of noiseunbiased estimates of central moments and cumulants up to the fourth order are presented under the assumption of independent Gaussian uncertainties, for weighted and unweighted statistics. These results are expected to be relevant for applications of the skewness and kurtosis estimators such as outlier detections, normality tests and in automated classification procedures. The comparison of estimators corrected and not corrected for noise biases is illustrated with simulations as a function of signaltonoise ratio, employing different sample sizes and weighting schemes.
Skewness and kurtosis unbiased by
Gaussian uncertainties
Lorenzo Rimoldini
Observatoire astronomique de l’Université de Genève, chemin des Maillettes 51, CH1290 Versoix, Switzerland
ISDC Data Centre for Astrophysics, Université de Genève, chemin d’Ecogia 16, CH1290 Versoix, Switzerland
lorenzo@rimoldini.info
Draft version: April 29, 2013
1 Introduction
Measurements generally provide an approximate description of real phenomena, because data acquisition compounds many processes which contribute, to a different degree, to instrumental errors (e.g., related to sensitivity or systematic biases) and uncertainties of statistical nature. While instrumental effects are addressed before data analysis, statistical uncertainties propagate in subsequent processing and can affect both precision and accuracy of results, especially at low signaltonoise () ratios. Correcting for biases generated by noise can help the characterization and interpretation of weak signals, and in some cases improve a significant fraction of all data (e.g., the number of astronomical sources increases dramatically near the faint detection threshold, since there are many more sources far away than nearby).
In this paper, noiseunbiased estimates of central moments and cumulants up to the fourth order, which are often employed to characterize the shape of the distribution of data, are derived analytically. Some of the advantages of these estimators include the ease of computation and the ability to encapsulate important features in a few numbers. Skewness and kurtosis measure the degree of asymmetry and peakedness or weight of the tails of the distribution, respectively, and they are useful for the detection of outliers, the assessment of departures from normality of the data (D’Agostino, 1986), the classification of light variations of astronomical sources (Rimoldini, 2013a) and many other applications. Various estimators of skewness and kurtosis are available in the literature (e.g., Moors et al., 1996; Hosking, 1990; Groeneveld & Meeden, 1984; Bowley, 1920), some of which aim at mitigating the sensitivity to outliers of the conventional formulations. On the other hand, robust measures might miss important features of signals, especially when these are characterized by outliers (as in astronomical time series where stellar bursts or eclipses from binary systems represent rare events in the light curve) and weighting might help distinguish true outliers from spurious data (employing additional information such as the accuracy of each measurement), so the traditional forms of weighted central moments and cumulants are employed in this work.
Moments are usually computed on random variables. Herein, their application is extended to data generated from deterministic functions and randomized by the uneven sampling of a finite number of measurements and by their uncertainties, whereas the corresponding ‘population’ statistics are defined in the limit of an infinite regular sampling with no random or systematic errors. This scenario is common in astronomical time series, where measurements are typically nonregular due to observational constraints, they are unavoidably affected by noise, and sometimes also not very numerous: all of these aspects introduce some level of randomness in the characterization of the underlying signal of a star.
While the effects of sampling and sample size on time series are studied in Rimoldini (2013a, b), this work addresses the bias, precision and accuracy of estimators when measurements are affected (mostly) by Gaussian uncertainties. Bias is defined as the difference between expectation and population values and thus expresses a systematic deviation from the true value. Precision is described by the dispersion of measurements, while accuracy is related to the distance of an estimator from the true value and thus combines the bias and precision concepts (e.g., accuracy can be measured by the mean square error, defined by the sum of bias and uncertainty in quadrature).
Noiseunbiased expressions are provided for the variance, skewness and kurtosis (central moments and cumulants), weighted and unweighted, assuming Gaussian uncertainties and independent measurements. The dependence of noiseunbiased estimators on S/N is illustrated with simulations employing different sample sizes and two weighting schemes: the common inversesquared uncertainties and interpolationbased weights as described in Rimoldini (2013a). The latter demonstrated a significant improvement in the precision of weighted estimators at the high end.
This paper is organized as follows. The notation employed throughout is defined in Sec. 2, followed by the description of the method to estimate Gaussiannoise unbiased moments in Sec. 3. Noiseunbiased estimates of moments and cumulants (biased and unbiased by samplesize) are presented in Sections 4 and 5, in both weighted and unweighted formulations, and the special case of errorweighted estimators is presented in Sec. 6. The noiseunbiased estimators are compared with the uncorrected (noisebiased) counterparts with simulated signals as a function of ratio in Sec. 7, including weighted and unweighted schemes and two different sample sizes. Conclusions are drawn in Sec. 8, followed by detailed derivations of the noiseunbiased estimators in App. A.
2 Notation
For a set of measurements , the following quantities are defined.

Population central moments with mean , where denotes expectation, and cumulants , , (e.g., Stuart & Ord, 1969).

The sum of the th power of weights is defined as .

The mean of a generic set of elements associated with weights is .

Sample central moments and corresponding cumulants .

Samplesize unbiased estimates of central moments and cumulants , i.e., and .

The standardized skewness and kurtosis are defined as , , , and , with population values and . and satisfy consistency (for ) but are not unbiased in general (e.g., see Heijmans, 1999, for exceptions).

Noiseunbiased estimates of central moments and cumulants are denoted by an asterisk superscript.

No systematic errors are considered herein and random errors are simply referred to as errors or uncertainties.

Statistics weighted by the inversesquared uncertainties are called ‘errorweighted’ for brevity and interpolationbased weights computed in phase (Rimoldini, 2013a) are named ‘phase weights’.
3 Method
The goal is to derive an estimator as a function of observables (measurements with corresponding uncertainties ) which is unbiased by the noise in the data, i.e., such that the expectation equals the estimator in terms of the true (unknown) values aimed at by the measurements.
The noiseunbiased estimator is obtained with the following procedure and assumptions. If independent measurements are associated with independent Gaussian uncertainties , the expected value of the estimator is evaluated from measurements and the joint probability density , for given true values and measurement uncertainties :
(1) 
where
(2) 
As shown in App. A, the expectation of the estimators considered herein can be decomposed as
(3) 
Thus, the noisefree estimator can be estimated in terms of measurements and uncertainties by the noiseunbiased estimator , where and, by definition, . The term is derived first by computing and then by replacing terms depending on in with terms as a function of which satisfy the requirement (see App. A). A property often used in the following sections is that a noiseunbiased linear combination of estimators is equivalent to the linear combination of noiseunbiased estimators:
(4) 
where the coefficients are independent of the measurements .
4 Gaussiannoise unbiased sample moments and cumulants
Weighted sample central moments unbiased by Gaussian uncertainties, such as the variance , skewness , kurtosis and the respective cumulants about the weighted mean are derived assuming independent measurements , uncertainties and weights , as described in full detail in App. A. They are defined as follows:
(5)  
(6)  
(7)  
(8)  
(9) 
By definition, the above expressions satisfy
(10) 
The unweighted forms can be obtained by substituting (for all ) and (for all ) in all terms, leading to:
(11)  
(12)  
(13)  
(14)  
(15) 
5 Gaussiannoise and samplesize unbiased moments and cumulants
The estimates of weighted central moments which are unbiased by both samplesize and Gaussian uncertainties, such as the variance , skewness , kurtosis and the respective cumulants, are defined in terms of the noiseunbiased sample estimators as follows:
(16)  
(17)  
(18)  
(19) 
The derivation of the samplesize unbiased weighted estimators is described in Rimoldini (2013b). The corresponding unweighted forms can be achieved by direct substitution for all , leading to:
(20)  
(21)  
(22)  
(23) 
6 Special cases
If weights are related to measurement errors as , the noiseunbiased weighted sample moments and cumulants reduce to the following expressions:
(24)  
(25)  
(26)  
(27)  
(28) 
In the case of constant errors, i.e., for all , some of the unweighted estimators are equivalent or similar to their noiseunbiased counterparts:
Skewness:  (29)  
Kurtosis:  (30) 
where the approximation holds for large values of or ratios since
(31) 
considering that, for constant errors, and . For sample cumulants up to the fourth order, only the variance depends strongly on noise. However, this is an important estimator because it is often involved in definitions of standardized skewness ( and ) and kurtosis ( and ) as follows:
(32)  
(33) 
For consistency with the above definitions, the noiseunbiased equivalents are defined as
(34)  
(35) 
although the truly noiseunbiased expressions should have been computed on the ratios in Eqs (32)–(33). The application of Eqs (34)–(35) should generally be restricted to larger samples (e.g., ) with ratios greater than a few, in order to avoid nonpositive values of or .
7 Estimators as a function of signaltonoise ratio
Noisebiased and unbiased estimators are compared as a function of signaltonoise ratio with simulated data and different weighting schemes for specific signals, sampling and error laws. The values of the population moments of the continuous simulated periodic ‘true’ signal are computed averaging in phase as follows:
(36) 
7.1 Simulation
Simulated signals are described by a sinusoidal function to the fourth power, which has a nonzero skewness and thus makes it possible to evaluate the precision and accuracy of the skewness standardized by the estimated variance without simply reflecting the accuracy of the variance. The level is evaluated by the ratio of the standard deviation of the true signal and the root of the mean of squared measurement uncertainties (assumed independent of the signal). The signal is sampled and 1000 times at phases randomly drawn from a uniform distribution, while the ratio varies from 1 to 1000 and determines the uncertainties of measurements as follows:
(37)  
(38)  
(39) 
where the th measurement is drawn from a normal distribution of mean and variance . The latter is defined in terms of a variable randomly drawn from a uniform distribution so that measurement uncertainties vary by up to a factor of 3 for a given and ratio. Simulations were repeated times for each ratio (for and 1000).
The dependence of weighted estimators on sample size and the corresponding unbiased expressions were presented in Rimoldini (2013b). Herein, only large sample sizes are employed so that samplesize biases are negligible with respect to the ones resulting from small ratios. A sample signal and simulated data are illustrated in Fig. 1 for and . The reference population values of the mean, variance, skewness and kurtosis of the simulated signal are listed in table 1 of Rimoldini (2013b).
Error weights are defined by , while mixed errorphase weights follow Rimoldini (2013a), assuming phasesorted data:
(40)  
(41)  
(42)  
(43)  
(44) 
Weighting effectively decreases the sample size, since more importance is given to some data at the expense of other ones and results depend mostly on fewer ‘relevant’ measurements (e.g., weighting by the inversesquared uncertainties can worsen precision at high levels). Weighted procedures are desirable when the dispersion and bias of estimators from an effectively reduced sample size are smaller than the improvements in precision and accuracy (e.g., weighting by inversesquared uncertainties can improve both precision and accuracy at low ratios). Also, weighting might exploit correlations in the data to improve precision, as it is shown employing phase weights (Rimoldini, 2013a). Since correlated data do not satisfy the assumptions of the expressions derived herein, their application might return biased results. However, small biases could be justified if improvements in precision are significant and, depending on the extent of the application, larger biases could be mitigated with mixed weighting schemes, such as the one described by Eqs (40)–(44).
Estimators derived herein assume a single weighting scheme and combinations of estimators (like the variance and the mean in the standardized skewness and kurtosis) are expected to apply the same weights to terms associated with the same measurements. The function constitutes just an example to achieve a mixed weighting scheme: tuning parameters offer the possibility to control the transition from errorweighted to phaseweighted estimators (in the limits of low and high , respectively) and thus reach a compromise solution between precision and accuracy for all values of , according to the specific estimators, signals, sampling, errors, sample sizes and their distributions in the data.
7.2 Results
The results of simulations are illustrated for sample estimators, since the conclusions in Rimoldini (2013b) suggested that phaseweighted sample estimators can be more accurate and precise than the samplesize unbiased counterparts in most cases, especially for large sample sizes as considered herein.
Figure 2 illustrates the sample mean in the various scenarios considered in the simulations: sample sizes of and 1000, unweighted and with different weighting schemes (errorweighted, phaseweighted and combined errorphase weighted). While accuracy is the same in all cases, the best precision of the mean is achieved employing phase weights (including the low end, unlike other estimators).
Figures 3–16 compare noisebiased (‘uncorrected’ ) and noiseunbiased (‘corrected’ ) estimators as a function of , evaluating the following deviations from the population values:
(45)  
(46)  
(47)  
(48) 
in both weighted and unweighted cases, for and 1000. The dependence on is described in more details in Rimoldini (2013b). Estimators standardized by both true and estimated variance are presented to help interpret the behaviour of the ratios from their components.
All figures confirm that ‘corrected’ and ‘uncorrected’ estimators have similar precision and accuracy at high levels (typically for ). Noiseunbiased estimators are found to be the most accurate in all cases and over the whole range tested. Their precision is generally similar to the noiseuncorrected counterparts, apart from estimators standardized by the estimated variance, such as , and , for which the uncorrected version can be much more precise (although biased) for , typically. As expected, the precision of estimators employing measurements per sample was greater than the one obtained with sample sizes of .
Weighting by the inverse of squared measurement errors made the estimators slightly less precise at high ratios, but more precise and accurate at low levels (except for the mean).
Weighting by phase intervals led to a significant improvement in precision of all estimators in the limit of large ratios and a reduction of precision at low (apart from the case of the mean). Tuning parameters such as and in Eq. (40) were able to mitigate the imprecision at low reducing to the errorweighted results, which appeared to be the most accurate and precise in the limit of low ratios (in these simulations). This solution might provide a reasonable compromise between precision and accuracy of all estimators, at least for .
Figures 5–8 show that the skewness moment is quite unbiased by noise, while the standardized version is underestimated at high because of the overestimated variance (as shown in Figs 3–4). While the accuracy of deteriorates at low , its precision is much less affected by noise.
The kurtosis moment (Figs 9–12) is less precise and accurate than the noiseunbiased equivalent, and its normalization by the squared variance reduces dramatically its inaccuracy and imprecision (since and exhibit a similar trend as a function of ). The kurtosis cumulant , instead, is much closer to its noiseunbiased counterpart, as shown in Figs 13–16. The normalization of by the squared variance improves its precision at the cost of lower accuracy for : the bias of is similar to (greater than) the precision of for ().
The lower the level is, the less precise estimators are and the noiseunbiased variance can be underestimated (and even become nonpositive). Thus, the skewness and kurtosis estimators standardized by or , as in Eqs (34)–(35), should be avoided in circumstances that combine small sample sizes (up to a few dozens of elements) and low S/N ratios (of the order of a few or less).
Figures related to moments and cumulants of irregularly sampled sinusoidal signals are very similar to the ones presented herein, with the exception of , which would have a similar precision but with no bias, as a consequence of the null skewness of a sinusoidal signal (since the mean of estimates is zero, they are not biased by the standardization with an overestimated noisebiased variance).
From the comparison of noisebiased and unbiased estimators with different weighting schemes, it appears that, for large sample sizes, noiseunbiased phaseweighted estimators are usually the most accurate for (apart from the special cases of standardized skewness and kurtosis when their true value is zero). For noisy signals (e.g., ), error weighting seems the most appropriate, at least with Gaussian uncertainties, thus noiseunbiased errorphase weighted estimators can provide a satisfactory compromise in general. Further improvements might be achieved by tuning parameters better fitted to estimators and signals of interest, in view of specific requirements of precision and accuracy.
8 Conclusions
Exact expressions of noiseunbiased skewness and kurtosis were provided in the unweighted and weighted formulations, under the assumption of independent data and Gaussian uncertainties. Such estimators can be particularly useful in the processing, interpretation and comparison of data characterized by low regimes.
Simulations of an irregularly sampled skewed periodic signal were employed to compare noisebiased and unbiased estimators as a function of in the unweighted, inversesquared error weighted and phaseweighted schemes. While noiseunbiased estimators were found more accurate in general, they were less precise than the uncorrected counterparts at low ratios. The application of a mixed weighting scheme involving phase intervals and uncertainties was able to balance precision and accuracy on a wide range of levels. The effect of noiseunbiased estimators and different weighting schemes on the characterization and classification of astronomical time series is described in Rimoldini (2013a).
Acknowledgments
The author thanks M. Süveges for many discussions and valuable comments on the original manuscript.
References
 Bowley (1920) Bowley A.L., 1920, Elements of Statistics, Charles Scribner’s Sons, New York
 D’Agostino (1986) D’Agostino R.B., 1986, Goodnessoffit techniques, D’Agostino & Stephens eds., Marcel Dekker, New York, p. 367
 Groeneveld & Meeden (1984) Groeneveld R.A., Meeden G., 1984, The Statistician, 33, 391
 Heijmans (1999) Heijmans R., 1999, Statistical Papers, 40, 107
 Hosking (1990) Hosking J.R.M., 1990, J. R. Statist. Soc. B, 52, 105
 Moors et al. (1996) Moors J.J.A., Wagemakers R.Th.A., Coenen V.M.J., Heuts R.M.J., Janssens M.J.B.T., 1996, Statistica Neerlandica, 50, 417
 Rimoldini (2013a) Rimoldini L., 2013a, preprint (arXiv:1304.6616)
 Rimoldini (2013b) Rimoldini L., 2013b, preprint (arXiv:1304.6564)
 Stuart & Ord (1969) Stuart A., Ord J., 1969, Kendall’s Advanced Theory of Statistics, Charles Griffin & Co. Ltd, London
Appendix A Derivation of noiseunbiased moments
The derivations presented in this Appendix involve weighted estimators under the assumption of independent measurements, uncertainties and weights. Definitions and some of the relations often employed herein are listed below.

For brevity, , and and are implied to involve all (from the 1st to the th) terms, unless explicitly stated otherwise.

The following integral solutions are often employed:
(49) 
The expected value of a generic estimator of independent data with Gaussian uncertainties is computed as follows
(50) (51)