Abstract
This paper characterizes the conventional moment-based Schmeiser-Deutsch (S-D) class of distributions through the method of L-moments. The system can be used in a variety of settings such as simulation or modeling various processes. A procedure is also described for simulating S-D distributions with specified L-moments and L-correlations. The Monte Carlo results presented in this study indicate that the estimates of L-skew, L-kurtosis, and L-correlation associated with the S-D class of distributions are substantially superior to their corresponding conventional product-moment estimators in terms of relative bias—most notably when sample sizes are small.
1. Introduction
The conventional moment-based Schmeiser-Deutsch (S-D) [1] class of distributions has demonstrated to be useful for modeling or simulating phenomena in the contexts of operations research and industrial engineering. Some examples include modeling stochastic inventory processes and lead time distributions [2–7], unpaced line efficiency [8, 9], two-stage production systems [10], stochastic activity networks [11, 12], and the newsboy problem [13]. Further, it is also common practice for methodologists to investigate the Type I error and power properties associated with inferential statistics (e.g., [14]). In many cases, these investigations may only require an elementary transformation to produce distributions with specified values of conventional skew, kurtosis, and Pearson correlations (or Gaussian copulas). The S-D class of distributions is particularly well suited for this task as it is computationally efficient because it only requires the knowledge of four parameters and an algorithm that generates zero-one uniform pseudorandom deviates.
Specifically, the quantile function for generating S-D distributions can be succinctly described as in [1]: where is zero-one uniformly distributed. The values of and are the location and scale parameters, while and are the shape parameters that determine the skew and kurtosis. The system of equations for determining the parameters in (1.1) for a S-D distribution with prespecified values of mean, variance, skew, and kurtosis is given in the Appendix. Figure 1 gives an example of an S-D distribution.
There are problems associated with conventional moment-based estimators (e.g., skew and kurtosis in Figure 1) insofar as they can be substantially biased, have high variance, or can be influenced by outliers. For example, inspection of Figure 1 indicates, on average, that the estimates of attenuate 27.51% and 38.82% below their associated population parameters. Note that each estimate of in Figure 1 was calculated based on samples of size and the formulae currently used by most commercial software packages such as SAS, SPSS, and Minitab for computing skew and kurtosis.
However, -moment-based estimators, such as -skew (), -kurtosis (), and the -correlation, have been introduced to address the limitations associated with conventional moment-based estimators [15–18]. Specifically, some of the advantages that -moments have over conventional moments are that they (i) exist whenever the mean of the distribution exists, (ii) are nearly unbiased for all sample sizes and distributions, and (iii) are more robust in the presence of outliers. For example, the estimates of in Figure 1 are relatively much closer to their respective parameters with smaller relative standard errors than their corresponding conventional moment-based analogs . More specifically, the estimates of that were simulated are, on average, 9.19% below and 5.7% above their respective parameters.
Thus, the present aim here is to characterize the S-D class of distributions through the method of -moments. The characterization will enable researchers to model or simulate nonnormal distributions with specified values of -skew, -kurtosis, and -correlation. The rest of the paper is outlined as follows. In Section 2, a summary of univariate -moment theory is provided as well as additional properties associated with S-D distributions. In Section 3, the derivation of the system of equations for specifying values of -skew and -kurtosis for the S-D class of distributions is subsequently provided. In Section 4, the coefficient of -correlation is introduced and the equations are developed for determining intermediate correlations for specified -correlations associated with the S-D class of distributions. In Section 5, the steps for implementing a simulation procedure are described. Numerical examples and the results of a simulation are also provided to confirm the derivations and compare the new methodology with its conventional moment-based counterparts. In Section 6, the results of the simulation are discussed and concluding comments are made.
2. Preliminaries on Univariate L-Moments and the Schmeiser-Deutsch Class of Distributions
2.1. Univariate L-Moments
Let be identically and independently distributed random variables each with continuous pdf , cdf , order statistics denoted as , and -moments defined in terms of either linear combinations of (i) expectations of order statistics or (ii) probability weighted moments (). For the purposes considered herein, the first four -moments associated with are expressed as [16, pages 20–22] where the are determined from where . The coefficients associated with in (2.2) are obtained from shifted orthogonal Legendre polynomials and are computed as shown in [16, page 20].
The -moments and in (2.1) are measures of location and scale and are the arithmetic mean and one-half the coefficient of mean difference, respectively. Higher order -moments are transformed to dimensionless quantities referred to as -moment ratios defined as for , and where and are the analogs to the conventional measures of skew and kurtosis. In general, -moment ratios are bounded in the interval as is the index of -skew () where a symmetric distribution implies that all -moment ratios with odd subscripts are zero. Other smaller boundaries can be found for more specific cases. For example, the index of -kurtosis () has the boundary condition for continuous distributions of [19] .
2.2. The Schmeiser-Deutsch (S-D) Class of Distributions
The cdf and pdf associated with the S-D quantile function in (1.1) are expressed as in [1]: Setting in (2.4) produces symmetric S-D densities with a lower bound of kurtosis of as . Positive (negative) skew is produced for cases where () and , where reverses the direction of skew. For (), the unique mode (antimode) is located at and where produces uniform distributions. For example, the distribution in Figure 1 is bounded in the range of with mode, mean, and variance of 10, 10.042, and 0.0271, respectively. See the Appendix for the formulae for computing the moments associated with S-D distributions. In the next section, the system of -moments for the class of S-D distributions is derived.
3. L-Moments for the Schmeiser-Deutsch (S-D) Class of Distributions
The derivation of the first four -moments associated with the S-D class of distributions begins by defining the probability weighted moments based on (2.2) in terms of (1.1) as where and are the zero-one uniform cdf and pdf. As such, integrating (3.1) for and simplifying using (2.1) give , , , and as Thus, given user-specified values of , , , and , (3.2)–(3.5) can be numerically solved to obtain the parameters for , , , and . Inspection of (3.4) and (3.5) indicates that the solutions to and are independent of the location and scale parameters ( and ). As with the conventional S-D class of symmetric distributions () the lower-bound of -kurtosis is obtained as .
Table 1 gives four examples of S-D distributions. These four distributions are used in the simulation portion of this study in Section 5. Distribution 1 was used by Lau [13] for stochastic modeling related to the newsboy problem and Distribution 2 was used by Aardal et al. [2] in the context of modeling inventory-control systems. The values of conventional skew and kurtosis for the four distributions were determined based on (A.2) in the Appendix. In the next section, we introduce the topic of the -correlation and subsequently develop the methodology for simulating S-D distributions with specified -correlations.
4. L-Correlations for the Schmeiser-Deutsch (S-D) Class of Distributions
The -correlation [17, 18] is introduced by considering two random variables and with distribution functions and , respectively. The second -moments of and can alternatively be expressed as The second -comoments of toward and toward are As such, the -correlations of toward and toward are expressed as The -correlation in (4.5) (or (4.6)) is bounded such that where a value of () indicates a strictly increasing (decreasing) monotone relationship between the two variables. In general, we would also note that .
In the context of the -moment-based S-D class of distributions, suppose it is desired to simulate a -variate distribution based on quantile functions of the forms in (1.1) with a specified -correlation matrix and where each distribution has its own specified values of and . Let denote standard normal variables where the distribution functions and bivariate density function associated with and are expressed as Using (4.7), it follows that the th S-D distribution associated with (1.1) can be expressed as since is zero-one uniformly distributed. As such, using (4.5), the -correlation of toward can be evaluated using the solved values of the parameters for , a specified intermediate correlation (IC) in (4.9), and the following integral generally expressed as where it is required that the location and scale parameters (, ) in (4.10) are to be solved such that will have the values of and in (3.2) and (3.3), that is, set to the values of and associated with the unit normal distribution. Note that this requirement is not a limitation as the -correlation is invariant to linear transformations [17]. Further, we would point out that the purpose of the IC () in (4.9) and (4.10) is to adjust for the effect of the transformation , which is induced by the parameters, such that has its specified -correlation () toward . Analogously, the -correlation of toward is expressed as Note for the special case that if in (4.10) and in (4.11) have the same parameters, that is, ;; ; , then . Provided in Algorithm 1 is source code written in Mathematica [20] that implements the computation of an IC () based on (4.10). The details for simulating S-D distributions with specified values of -skew, -kurtosis, and -correlations are described in the next section.
| ||||||||||||||||||||||||||||||||||||||||||||
5. The Procedure for Simulation and Monte Carlo Study
To implement a method for simulating S-D distributions with specified -moments and -correlations we suggest the following six steps.(1)Specify the -moments for transformations of the forms in (1.1), that is, and obtain the solutions for the parameters of , , , and by solving (3.2)–(3.5) with and using the specified values of -skew () and -kurtosis () for each distribution. Specify a matrix of -correlations () for toward , where .(2)Compute the (Pearson) intermediate correlations (ICs) by substituting the solutions of the parameters from Step (1) into (4.10) and then numerically integrate to solve for (see Algorithm 1 for an example). Repeat this step separately for all pairwise combinations of correlations.(3)Assemble the ICs into a matrix and decompose this matrix using a Cholesky factorization. Note that this step requires the IC matrix to be positive definite.(4)Use the results of the Cholesky factorization from Step (3) to generate standard normal variables correlated at the intermediate levels as follows: where are independent standard normal random variables and where represents the element in the th row and the th column of the matrix associated with the Cholesky factorization performed in Step (3).(5)Substitute from Step (4) into the following Taylor series-based expansion for the standard normal cdf [21]: where denotes the standard normal pdf and where the absolute error associated with (5.2) is less than .(6)Substitute the zero-one uniform deviates, , generated from Step (5) into the equations of the form of to generate the S-D distributions with the specified -moments and -correlations.To demonstrate the steps above, and evaluate the proposed method, a comparison between the proposed -moment and conventional product-moment-based procedures is subsequently described. Specifically, the parameters for the distributions in Table 1 are used as a basis for a comparison using the specified correlation matrix in Table 2. Tables 3 and 4 give the solved IC matrices for the conventional moment and -moment-based methods, respectively. See Algorithm 2 for an example of computing ICs for the conventional method. Tables 5 and 6 give the results of the Cholesky decompositions on the IC matrices, which are then used to create with the specified ICs by making use of the formulae given in (5.1) of Step (4) with . The values of are subsequently transformed to using (5.2) and then substituted into equations of the forms in (1.1) to produce for both methods.
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
In terms of the simulation, a Fortran algorithm was written for both methods to generate 25,000 independent sample estimates for the specified parameters of (i) conventional skew (), kurtosis (), and Pearson correlation (); (ii) -skew (), -kurtosis (), and -correlation (). All estimates were based on sample sizes of and . The formulae used for computing estimates of were based on Fisher’s -statistics, that is, the formulae currently used by most commercial software packages such as SAS, SPSS, and Minitab, for computing indices of skew and kurtosis (where for the standard normal distribution). The formulae used for computing estimates of were Headrick’s equations and [22]. The estimate for was based on the usual formula for the Pearson product-moment of correlation statistic and the estimate for was computed based on (4.5) using the empirical forms of the cdfs in (4.1) and (4.3). The estimates for and were both transformed using Fisher’s transformation. Bias-corrected-accelerated-bootstrapped average (mean) estimates, confidence intervals (C.I.s), and standard errors were subsequently obtained for the estimates associated with the parameters (, , , ) using 10,000 resamples via the commercial software package Spotfire S+ [23]. The bootstrap results for the estimates of the means and C.I.s associated with and were transformed back to their original metrics (i.e., estimates for and ). Further, if a parameter () was outside its associated bootstrap C.I., then an index of relative bias (RB) was computed for the estimate () as: . Note that the small amount of bias associated with any bootstrap C.I. containing a parameter was considered negligible and thus not reported. The results of the simulation are reported in Tables 7–12 and are discussed in the next section.
6. Discussion and Conclusion
One of the primary advantages that -moments have over conventional moment-based estimators is that they can be far less biased when sampling is from distributions with more severe departures from normality (e.g. [16, 21]). Inspection of the simulation results in Tables 7 and 8 of this study clearly indicates that this is also the case for the S-D class of distributions. Specifically, the superiority that estimates of -moment ratios () have over their corresponding conventional moment based counterparts () is obvious. For example, with samples of size the estimates of skew and kurtosis for Distribution 4 (Table 7) were, on average, 85.65% and 74.74% of their associated population parameters, whereas the estimates of -skew and -kurtosis were 97.34% and 92.70% of their respective parameters. Similar results were also obtained for Distributions 1 and 3 and thus not reported. It is also evident from Tables 7 and 8 that -skew and -kurtosis are more efficient as their relative standard errors RSE = (standard error/estimate) × 100 are smaller than the conventional estimators of skew and kurtosis. For example, in terms of Distribution 4, inspection of Tables 7 and 8 () indicates RSE measures of RSE() = 0.1533% and RSE() compared with RSE() and RSE() . This demonstrates that -skew and -kurtosis have more precision because they have less variance around their estimates.
Presented in Tables 9–12 are the results associated with the conventional Pearson and -correlations. Inspection of Tables 9 and 10 indicates that the -correlation is substantially superior to the Pearson correlation in terms of relative bias for small sample sizes. For example, in terms of a moderate correlation (Table 9, , ) the relative bias for Distributions 3 and 4 was 23.95% for the Pearson correlation compared to 5.43% for the -correlation (Table 10, , ). For large sample sizes (Tables 11 and 12, ), both procedures performed adequately as their estimates were in close proximity with their respective population parameters.
In summary, the proposed -moment-based S-D class of distributions is an attractive alternative to the conventional moment-based S-D system. In particular, the -moment-based system has distinct advantages when leptokurtic distributions and small sample sizes are of concern. Finally, we would note that Mathematica Version 8.0 [20] source code is available from the authors for implementing the -moment-based method.
Appendix
System of Conventional Moment-Based Equations for S-D Distributions
The moments () associated with the S-D class of distributions in (1.1) can be determined from where is the zero-one uniform pdf. The mean, variance, skew, and kurtosis are defined in general as in [24]: The moments associated with the location and scale parameters in (A.2) are The moments associated with the shape parameters of skew and kurtosis in (A.2) are