ISRN Applied MathematicsVolume 2012 (2012), Article ID 475781, 16 pagesdoi:10.5402/2012/475781

Research Article

An L-Moment-Based Analog for the Schmeiser-Deutsch Class of Distributions

Todd C. Headrick¹ and Mohan D. Pant²

¹Section on Statistics and Measurement, Department of EPSE, Southern Illinois University Carbondale, 222-J Wham Building , Carbondale, IL 62901-4618, USA²University of Texas at Arlington, 320B Science Hall, Arlington, TX 76019, USA

Received 28 June 2012; Accepted 7 August 2012

Academic Editors: I. Doltsinis and E. Yee

Copyright © 2012 Todd C. Headrick and Mohan D. Pant. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

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.

Figure 1: A graph of an S-D distribution based on the pdf in (2.4). The values of skew and kurtosis were determined based on (A.2) for and in the Appendix. The values of and were determined based on (3.4) and (3.5). The estimates (; ) and bootstrap confidence intervals (C.I.s) were based on resampling 25,000 statistics. Each sample statistic was based on a sample size of .

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.

Table 1: Four S-D distributions and their associated -moments, conventional moments (C-moments), and S-D parameters used in the simulation. The values of -skew () and -kurtosis () are based on (3.4) and (3.5). The values of skew () and kurtosis () are based on the expressions in (A.2) of the Appendix.

Distribution
-moments
C-moments
S-D parameters

-moment
C-moment

(1) Moderate positive skew and platykurtic

(2) Symmetric and mesokurtic

(3) Positive skew and leptokurtic

(4) Negative skew and leptokurtic

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.

(* Intermediate Correlation for Distributions 1 and 2. See Table 4. *)

;

Needs[“MultivariateStatistics”]

DF[MultinormalDistribution[, , ], ];

DF[NormalDistribution[0, 1], ];

DF[NormalDistribution[0, 1], ];

(* S-D -moment Parameters for Distribution in Table 1*)

;

;

;

;

(* Quantile function for (1.1) *)

;

;

Piecewise[, , ]

(* Compute the specified correlation in Table 2 for Distributions 1 and 2 *)

0.70

Algorithm 1: Mathematica source code for computing intermediate correlations for specified -correlations. The example is for distributions and () in Table 1. See Tables 2 and 4.

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.

Table 2: Specified correlation matrix for the distributions in Table 1.

1
2
3
4

1
1

2
0.70
1

3
0.60
0.60
1

4
0.35
0.60
0.40
1

Table 3: Intermediate correlations for the conventional moment procedure.

1
2
3
4

1
1

2
0.761371
1

3
0.945345
0.813672
1

4
0.810720
0.956258
0.910196
1

Table 4: Intermediate correlations for the -moment procedure.

1
2
3
4

1
1

2
0.714645
1

3
0.615758
0.584908
1

4
0.362634
0.584908
0.332108
1

Table 5: Cholesky decompositions for the conventional moment procedure.

0

0
0

0
0
0

Table 6: Cholesky decompositions for the -moment procedure.

0

0
0

0
0
0

(* Intermediate Correlation for Distributions 1 and 2. See Table 3. *)

0.761371;

Needs[“MultivariateStatistics”]

DF;

DF;

DF;

(* S-D C-Moment Parameters Distributions and in Table 1*)

;

;

;

;

;

;

;

;

(* Quantile functions for (1.1) *)

;

;

Piecewise

;

;

Piecewise

(* Compute the specified correlation in Table 2 for Distributions 1 and 2. *)

NIntegrate[, , , MethodMultiDimensional]

0.70

Algorithm 2: Mathematica source code for computing intermediate correlations for specified conventional Pearson correlations. The example is for distributions and () in Table 1. See also Tables 2 and 3.

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.

Table 7: Skew () and kurtosis () results for distributions 2 and 4 in Table 1.

Dist
Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %

2

0.0019
−0.0055, 0.0094
0.00381
—

1.166
1.1496, 1.1821
0.00825
51.59

4

−3.581
−3.5917, −3.5702
0.00549
−14.35

13.61
13.5278, 13.6973
0.04337
−25.26

2

−0.0005
−0.0013, 0.0004
0.00042
—

0.7798
0.7778, 0.7817
0.00100
1.38

4

−4.190
−4.1936, −4.1866
0.00177
0.22

18.44
18.4046, 18.4738
0.01775
1.26

Table 8: -skew () and -kurtosis () results for distributions 2 and 4 in Table 1.

Dist
Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %

2

0.0002
−0.0014, 0.0018
0.00079
—

0.2981
0.2968, 0.2993
0.00063
4.34

4

−0.8638
−0.8647, −0.8628
0.00049
2.66

0.6956
0.6935, 0.6975
0.00102
7.30

2

−0.0001
−0.0003, 0.0001
0.00010
—

0.2860
0.2859, 0.2862
0.00008
0.11

4

−0.8422
−0.8423,
0.00007
0.095

0.6498
0.6495, 0.6501
0.00015
0.23

Table 9: Correlation results for the conventional moment procedure .

Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %

0.7083
0.7071, 0.7095
0.00121
1.19

0.6396
0.6380, 0.6414
0.00146
6.60

0.3837
0.3826, 0.3848
0.00065
9.63

0.6342
0.6327, 0.6356
0.00122
5.70

0.6342
0.6330, 0.6353
0.00097
5.70

0.4958
0.4919, 0.5000
0.00272
23.95

Table 10: Correlation results for the -moment procedure .

Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %

0.7155
0.7140, 0.7170
0.00152
2.21

0.6155
0.6137, 0.6173
0.00148
2.58

0.3588
0.3562, 0.3611
0.00143
2.51

0.6116
0.6099, 0.6136
0.00149
1.93

0.6057
0.6038, 0.6075
0.00151
0.95

0.4217
0.4177, 0.4259
0.00257
5.43

Table 11: Correlation results for the conventional moment procedure .

Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %

0.7001
0.6999, 0.7003
0.00018
—

0.6004
0.6001, 0.6006
0.00021
0.07

0.3502
0.3500, 0.3504
0.00011
—

0.6002
0.5999, 0.6004
0.00019
—

0.6004
0.6002, 0.6006
0.00015
0.07

0.4009
0.4004, 0.4013
0.00027
0.23

Table 12: Correlation results for the -moment procedure .

Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %

0.7002
0.7000, 0.7004
0.00022
—

0.6005
0.6003, 0.6008
0.00021
0.08

0.3501
0.3498, 0.3505
0.00021
—

0.6000
0.5997, 0.6002
0.00022
—

0.6002
0.5999, 0.6005
0.00021
—

0.3984
0.3979, 0.3988
0.00035
−0.40

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

References

B. W. Schmeiser and S. J. Deutsch, “A versatile four parameter family of probability distributions suitable for simulation,” American Institute of Industrial Engineers (AIIE) Transactions, vol. 9, no. 2, pp. 176–182, 1977.

K. Aardal, Ö Jonsson, and H. Jönsson, “Optimal inventory policies with service-level contraints,” Journal of the Operations Research Society, vol. 40, no. 1, pp. 65–73, 1989.

O. Tang and R. W. Grubbström, “On using higher-order moments for stochastic inventory systems,” International Journal of Production Economics, vol. 104, no. 2, pp. 454–461, 2006.

A. Akcay and S. Tayur, “Improved inventory targets in the presence of limited historical demand data,” Manufacturing and Service Operations Management (M&SOM), vol. 13, no. 3, pp. 297–309, 2011.

J. F. Kottas and H. Lau, “A realistic approach for modeling stochastic lead time distributions,” American Institute of Industrial Engineers (AIIE) Transactions, vol. 11, no. 1, pp. 54–60, 1979.

J. F. Kottas and H. S. Lau, “Inventory control with general demand and lead time distributions,” International Journal of Systems Science, vol. 10, no. 5, pp. 485–492, 1979.

J. F. Kottas and H. S. Lau, “The use of versatile distribution families in some stochastic inventory calculations,” Journal of the Operations Research Society, vol. 31, no. 5, pp. 393–403, 1980.

H. Lau and G. E. Martin, “The effects of skewness and kurtosis of processing times in unpaced lines,” International Journal of Production Research, vol. 25, no. 10, pp. 1483–1492, 1987.

G. E. Martin, “Predictive formulae for unpaced line efficiency,” International Journal of Production Research, vol. 31, no. 8, pp. 1981–1990, 1993.

H. Lau, “The production rate of a two-stage system with stochastic processing times,” International Journal of Production Research, vol. 24, no. 2, pp. 401–412, 1986.

H. Lau, A. Lau, and C. Ho, “Improved moment-estimation formulas using more than three subjective fractiles,” Management Science, vol. 44, no. 3, pp. 346–351, 1998.

C. Somarajan and H.-S. Lau, “Stochastic activity networks: a four-parameter approach for completion times,” International Journal of Systems Science, vol. 25, no. 10, pp. 1607–1619, 1994.

H. S. Lau, “The newsboy problem under alternative optimization objectives,” Journal of the Operations Research Society, vol. 31, no. 6, pp. 525–535, 1980.

T. C. Headrick, Statistical Simulation: Power Method Polynomials and Other Transformations, Chapman & Hall/CRC, Boca Raton, Fla, USA, 2010.

J. R. M. Hosking, “-moments: analysis and estimation of distributions using linear combinations of order statistics,” Journal of the Royal Statistical Society, vol. 52, no. 1, pp. 105–124, 1990.

J. R. M. Hosking and J. R. Wallis, Regional Frequency Analysis: An Approach Based on L-Moments, Cambridge University Press, Cambridge, UK, 1997.

R. Serfling and P. Xiao, “A contribution to multivariate -moments: -comoment matrices,” Journal of Multivariate Analysis, vol. 98, no. 9, pp. 1765–1781, 2007.

T. C. Headrick and M. D. Pant, “Simulating non-normal distributions with specified L-moments and L-correlations,” Statistica Neerlandica. In press

M. C. Jones, “On some expressions for variance, covariance, skewness and -moments,” Journal of Statistical Planning and Inference, vol. 126, no. 1, pp. 97–106, 2004.

Wolfram Research, Mathematica, Version 8.0, Wolfram Research, Champaign, Ill, USA, 2010.

G. Marsaglia, “Evaluating the normal distribution,” Journal of Statistical Software, vol. 11, no. 5, pp. 1–11, 2004.

T. C. Headrick, “A characterization of power method transformations through -moments,” Journal of Probability and Statistics, vol. 2011, Article ID 497463, 22 pages, 2011.

TIBCO Software, TIBCO Spotfire S+ 8.1 for Windows, Palo Alto, Calif, USA, 2008.

T. C. Headrick, R. K. Kowalchuk, and Y. Sheng, “Parametric probability densities and distribution functions for Tukey -and- transformations and their use for fitting data,” Applied Mathematical Sciences, vol. 2, no. 9–12, pp. 449–462, 2008.