# AnL-Moment-Based Analog for the Schmeiser-Deutsch Class of Distributions

Authors: Todd C. Headrick; Mohan D. Pant
Journal: ISRN Applied Mathematics
 (2012)

---

## Abstract

This paper characterizes the conventional moment-based Schmeiser-Deutsch (S-D) class of distributions through the method ofL-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.


---

## Body

## 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]:
(1.1)X=q(U)={γ1-γ2(γ4-U)γ3,forU≤γ4,γ1+γ2(U-γ4)γ3,forU>γ4,
where U is zero-one uniformly distributed. The values of γ1∈(-∞,+∞) and γ2≥0 are the location and scale parameters, while γ3≥0 and γ4∈[0,1] 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 α3 and α4 in the Appendix. The values of τ3 and τ4 were determined based on (3.4) and (3.5). The estimates (α^3,4; τ^3,4) and bootstrap confidence intervals (C.I.s) were based on resampling 25,000 statistics. Each sample statistic was based on a sample size of n=25.There are problems associated with conventional moment-based estimators (e.g., skew and kurtosisα^3,4 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 α^3,4 attenuate 27.51% and 38.82% below their associated population parameters. Note that each estimate of α^3,4 in Figure 1 was calculated based on samples of size n=25 and the formulae currently used by most commercial software packages such as SAS, SPSS, and Minitab for computing skew and kurtosis.However,L-moment-based estimators, such as L-skew (τ^3), L-kurtosis (τ^4), and the L-correlation, have been introduced to address the limitations associated with conventional moment-based estimators [15–18]. Specifically, some of the advantages that L-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 τ^3,4 in Figure 1 are relatively much closer to their respective parameters τ3,4 with smaller relative standard errors than their corresponding conventional moment-based analogs α^3,4. More specifically, the estimates of τ^3,4 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 ofL-moments. The characterization will enable researchers to model or simulate nonnormal distributions with specified values of L-skew, L-kurtosis, and L-correlation. The rest of the paper is outlined as follows. In Section 2, a summary of univariate L-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 L-skew and L-kurtosis for the S-D class of distributions is subsequently provided. In Section 4, the coefficient of L-correlation is introduced and the equations are developed for determining intermediate correlations for specified L-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 UnivariateL-Moments and the Schmeiser-Deutsch Class of Distributions

### 2.1. UnivariateL-Moments

LetY1,…,Yj,…,Yn be identically and independently distributed random variables each with continuous pdf f(y), cdf F(y), order statistics denoted as Y1:n≤⋯≤Yj:n≤⋯≤Yn:n, and L-moments defined in terms of either linear combinations of (i) expectations of order statistics or (ii) probability weighted moments (βi). For the purposes considered herein, the first four L-moments associated with Yj:n are expressed as [16, pages 20–22]
(2.1)λ1=E[Y1:1]=β0,λ2=12E[Y2:2-Y1:2]=2β1-β0,λ3=13E[Y3:3-2Y2:3+Y1:3]=6β2-6β1+β0,λ4=14E[Y4:4-3Y3:4+3Y2:4-Y1:4]=20β3-30β2+12β1-β0,
where the βi are determined from
(2.2)βi=∫y{F(y)}if(y)dy,
where i=0,…,3. The coefficients associated with βi in (2.2) are obtained from shifted orthogonal Legendre polynomials and are computed as shown in [16, page 20].TheL-moments λ1 and λ2 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 L-moments are transformed to dimensionless quantities referred to as L-moment ratios defined as τr=λr/λ2 for r≥3, and where τ3 and τ4 are the analogs to the conventional measures of skew and kurtosis. In general, L-moment ratios are bounded in the interval -1<τr<1  as is the index of L-skew (τ3) where a symmetric distribution implies that all L-moment ratios with odd subscripts are zero. Other smaller boundaries can be found for more specific cases. For example, the index of L-kurtosis (τ4) has the boundary condition for continuous distributions of [19] (5τ32-1)/4<τ4<1.


### 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]:
(2.3)FX(x)={γ4-(γ1-xγ2)1/γ3,ifγ1-γ2γ4γ3≤x≤γ1γ4+(x-γ1γ2)1/γ3,ifγ1≤x≤γ1+γ2(1-γ4)γ3,(2.4)fX(x)=1γ2γ3|γ1-xγ2|(1-γ3)/γ3,forγ1-γ2γ4γ3≤x≤γ1+γ2(1-γ4)γ3.
Setting γ4=0.50 in (2.4) produces symmetric S-D densities with a lower bound of kurtosis of α4=-2 as γ3→0. Positive (negative) skew is produced for cases where 0≤γ4<0.50 (0.50<γ4≤1) and γ3>1, where γ3<1 reverses the direction of skew. For γ3>1 (0<γ3<1), the unique mode (antimode) is located at γ1 and where γ3=1 produces uniform distributions. For example, the distribution in Figure 1 is bounded in the range of 9.9372≤x≤11.2748 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 L-moments for the class of S-D distributions is derived.


## 2.1. UnivariateL-Moments

LetY1,…,Yj,…,Yn be identically and independently distributed random variables each with continuous pdf f(y), cdf F(y), order statistics denoted as Y1:n≤⋯≤Yj:n≤⋯≤Yn:n, and L-moments defined in terms of either linear combinations of (i) expectations of order statistics or (ii) probability weighted moments (βi). For the purposes considered herein, the first four L-moments associated with Yj:n are expressed as [16, pages 20–22]
(2.1)λ1=E[Y1:1]=β0,λ2=12E[Y2:2-Y1:2]=2β1-β0,λ3=13E[Y3:3-2Y2:3+Y1:3]=6β2-6β1+β0,λ4=14E[Y4:4-3Y3:4+3Y2:4-Y1:4]=20β3-30β2+12β1-β0,
where the βi are determined from
(2.2)βi=∫y{F(y)}if(y)dy,
where i=0,…,3. The coefficients associated with βi in (2.2) are obtained from shifted orthogonal Legendre polynomials and are computed as shown in [16, page 20].TheL-moments λ1 and λ2 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 L-moments are transformed to dimensionless quantities referred to as L-moment ratios defined as τr=λr/λ2 for r≥3, and where τ3 and τ4 are the analogs to the conventional measures of skew and kurtosis. In general, L-moment ratios are bounded in the interval -1<τr<1  as is the index of L-skew (τ3) where a symmetric distribution implies that all L-moment ratios with odd subscripts are zero. Other smaller boundaries can be found for more specific cases. For example, the index of L-kurtosis (τ4) has the boundary condition for continuous distributions of [19] (5τ32-1)/4<τ4<1.


## 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]:
(2.3)FX(x)={γ4-(γ1-xγ2)1/γ3,ifγ1-γ2γ4γ3≤x≤γ1γ4+(x-γ1γ2)1/γ3,ifγ1≤x≤γ1+γ2(1-γ4)γ3,(2.4)fX(x)=1γ2γ3|γ1-xγ2|(1-γ3)/γ3,forγ1-γ2γ4γ3≤x≤γ1+γ2(1-γ4)γ3.
Setting γ4=0.50 in (2.4) produces symmetric S-D densities with a lower bound of kurtosis of α4=-2 as γ3→0. Positive (negative) skew is produced for cases where 0≤γ4<0.50 (0.50<γ4≤1) and γ3>1, where γ3<1 reverses the direction of skew. For γ3>1 (0<γ3<1), the unique mode (antimode) is located at γ1 and where γ3=1 produces uniform distributions. For example, the distribution in Figure 1 is bounded in the range of 9.9372≤x≤11.2748 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 L-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 fourL-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
(3.1)βi=∫0γ4q(u,γ1,γ2,γ3,γ4){F(u)}if(u)du+∫γ41q(u,γ1,γ2,γ3,γ4){F(u)}if(u)du,
where F(u)=u and f(u)=1 are the zero-one uniform cdf and pdf. As such, integrating (3.1) for i=0,1,2,3 and simplifying using (2.1) give λ1, λ2, τ3, and τ4 as
(3.2)λ1=γ1+γ2(1-γ4)1+γ3(1+γ3)-γ1γ4+γ4(γ1-γ2γ4γ3(1+γ3)),(3.3)λ2=(-γ2γ41+γ3(2γ4-γ3-2)+γ2(1-γ4)1+γ3(γ3+2γ4))((1+γ3)(2+γ3)),(3.4)τ3=-{((1-γ4)1+γ3(γ3(1-6γ4)-γ32+6(1-2γ4)γ4)+γ41+γ3(6+5γ3+γ32-6(3+γ3)γ4+12γ42))/((3+γ3)((2+γ3-2γ4)γ41+γ3+(1-γ4)1+γ3(γ3+2γ4)))},(3.5)τ4={-(1-γ4)1+γ3((γ3-2)(γ3-1)γ3+12(γ3-2)(γ3-1)γ4+60(γ3-2)γ42+120γ43)+γ41+γ3(-24-26γ3-9γ32-γ33+12(3+γ3)(4+γ3)γ4-60(4+γ3)γ42+120γ43)}/{(3+γ3)(4+γ3)(γ41+γ3(2γ4-γ3-2)-(1-γ4)1+γ3(γ3+2γ4))}.
Thus, given user-specified values of λ1, λ2, τ3, and τ4, (3.2)–(3.5) can be numerically solved to obtain the parameters for γ1, γ2, γ3, and γ4. Inspection of (3.4) and (3.5) indicates that the solutions to τ3 and τ4 are independent of the location and scale parameters (λ1 and λ2). As with the conventional S-D class of symmetric distributions (γ4=0.5) the lower-bound of L-kurtosis τ4=-0.25 is obtained as γ3→0.Table1 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 L-correlation and subsequently develop the methodology for simulating S-D distributions with specified L-correlations.Table 1
Four S-D distributions and their associatedL-moments, conventional moments (C-moments), and S-D parameters used in the simulation. The values of L-skew (τ3) and L-kurtosis (τ4) are based on (3.4) and (3.5). The values of skew (α3) and kurtosis (α4) are based on the expressions in (A.2) of the Appendix.
Distribution
L-moments
C-moments
S-D parameters
L-moment
C-moment
(1)  Moderate positive skew and platykurtic
λ
1
=
0
α
1
=
0
γ
1
=
-
0.8429
γ
1
=
-
0.8427
λ
2
=
1
/
π
α
2
=
1
γ
2
=
4.6222
γ
2
=
4.6213
τ
3
=
0.1647
α
3
=
0.5524
γ
3
=
1.5
γ
3
=
1.5
τ
4
=
0.01606
α
4
=
-
0.9261
γ
4
=
0.25
γ
4
=
0.25
(2)  Symmetric and mesokurtic
λ
1
=
0
α
1
=
0
γ
1
=
0
γ
1
=
0
λ
2
=
1
/
π
α
2
=
1
γ
2
=
22.568
γ
2
=
21.166
τ
3
=
0
α
3
=
0
γ
3
=
3
γ
3
=
3
τ
4
=
0.2857
α
4
=
0.7692
γ
4
=
0.50
γ
4
=
0.50
(3)  Positive skew and leptokurtic
λ
1
=
0
α
1
=
0
γ
1
=
-
0.2344
γ
1
=
-
0.1373
λ
2
=
1
/
π
α
2
=
1
γ
2
=
59196.04
γ
2
=
34569.55
τ
3
=
0.3314
α
3
=
3
γ
3
=
12.89
γ
3
=
12.89
τ
4
=
0.7011
α
4
=
15
γ
4
=
0.4853
γ
4
=
0.4853
(4)  Negative skew  and leptokurtic
λ
1
=
0
α
1
=
0
γ
1
=
0.6148
γ
1
=
0.3043
λ
2
=
1
/
π
α
2
=
1
γ
2
=
2960.06
γ
2
=
1465.19
τ
3
=
-
0.8414
α
3
=
-
4.181
γ
3
=
15
γ
3
=
15
τ
4
=
0.6483
α
4
=
18.21
γ
4
=
0.70
γ
4
=
0.70


## 4.L-Correlations for the Schmeiser-Deutsch (S-D) Class of Distributions

TheL-correlation [17, 18] is introduced by considering two random variables Yj and Yk with distribution functions F(Yj) and F(Yk), respectively. The second L-moments of Yj and Yk can alternatively be expressed as
(4.1)λ2(Yj)=2Cov(Yj,F(Yj)),(4.2)λ2(Yk)=2Cov(Yk,F(Yk)).
The second L-comoments of Yj toward Yk and Yk toward Yj are
(4.3)λ2(Yj,Yk)=2Cov(Yj,F(Yk)),(4.4)λ2(Yk,Yj)=2Cov(Yk,F(Yj)).
As such, the L-correlations of Yj toward Yk and Yk toward Yj are expressed as
(4.5)ηjk=λ2(Yj,Yk)λ2(Yj),(4.6)ηkj=λ2(Yk,Yj)λ2(Yk).
The L-correlation in (4.5) (or (4.6)) is bounded such that -1≤ηjk≤1 where a value of ηjk=1 (ηjk=-1) indicates a strictly increasing (decreasing) monotone relationship between the two variables. In general, we would also note that ηjk≠ηkj.In the context of theL-moment-based S-D class of distributions, suppose it is desired to simulate a T-variate distribution based on quantile functions of the forms in (1.1) with a specified L-correlation matrix and where each distribution has its own specified values of τ3 and τ4. Let Z1,…,ZT denote standard normal variables where the distribution functions and bivariate density function associated with Zj and Zk are expressed as
(4.7)Φ(zj)=Pr{Zj≤zj}=∫-∞zj(2π)-1/2exp{-wj22}dwj,(4.8)Φ(zk)=Pr{Zk≤zk}=∫-∞zk(2π)-1/2exp{-wk22}dwk,(4.9)fjk=(2π(1-ρjk2)1/2)-1exp{-(2(1-ρjk2))-1(zj2+zk2-2ρjkzjzk)}.
Using (4.7), it follows that the jth S-D distribution associated with (1.1) can be expressed as qj(Φ(Zj)) since Φ(Zj) is zero-one uniformly distributed. As such, using (4.5), the L-correlation of qj(Φ(Zj)) toward qk(Φ(Zk)) can be evaluated using the solved values of the parameters for qj(Φ(Zj)), a specified intermediate correlation (IC) ρjk in (4.9), and the following integral generally expressed as
(4.10)ηjk=2π∬-∞+∞qj(Φ(zj),γ1j,γ2j,γ3j,γ4j)Φ(zk)fjkdzjdzk,
where it is required that the location and scale parameters (γ1j, γ2j) in (4.10) are to be solved such that qj(Φ(Zj)) will have the values of λ1=0 and λ2=1/π in (3.2) and (3.3), that is, set to the values of λ1 and λ2 associated with the unit normal distribution. Note that this requirement is not a limitation as the L-correlation is invariant to linear transformations [17]. Further, we would point out that the purpose of the IC (ρjk) in (4.9) and (4.10) is to adjust for the effect of the transformation qj(Φ(Zj)), which is induced by the parameters, such that qj(Φ(Zj)) has its specified L-correlation (ηjk) toward qk(Φ(Zk)). Analogously, the L-correlation of qk(Φ(Zk)) toward qj(Φ(Zj)) is expressed as
(4.11)ηkj=2π∬-∞+∞qk(Φ(zk),γ1k,γ2k,γ3k,γ4k)Φ(zj)fjkdzkdzj.
Note for the special case that if qj(Φ(Zj)) in (4.10) and qk(Φ(Zk)) in (4.11) have the same parameters, that is, γ1j=γ1k;γ2j=γ2k; γ3j=γ3k; γ4j=γ4k, then ηjk=ηkj. Provided in Algorithm 1 is source code written in Mathematica [20] that implements the computation of an IC (ρjk) based on (4.10). The details for simulating S-D distributions with specified values of L-skew, L-kurtosis, and L-correlations are described in the next section.Algorithm 1:Mathematica source code for computing intermediate correlations for specified L-correlations. The example is for distributions j=1 and k=2 (η12) in Table 1. See Tables 2 and 4.
(* Intermediate Correlation for Distributions 1 and 2. See Table4. *)
ρ
j
k
=
0.714645;
Needs[“MultivariateStatistics”]
f
j
k
=
PDF[MultinormalDistribution[{0,0}, {{1,ρjk}, {ρjk,1}}], {Zj,Zk}];
F
j
=
CDF[NormalDistribution[0, 1], Zj];
F
k
=
CDF[NormalDistribution[0, 1], Zk];
(* S-DL-moment Parameters for Distribution j=1 in Table 1*)
γ
1
=
-
0.8429;
γ
2
=
4.6222;
γ
3
=
1.5;
γ
4
=
0.25;
(* Quantile function for (1.1) *)
q
1
j
=
γ
1
j
-
γ
2
j
*
(
γ
4
j
-
F
j
)
γ
3
j;
q
2
j
=
γ
1
j
+
γ
2
j
*
(
F
j
-
γ
4
j
)
γ
3
j;
X
j
=Piecewise[{{q1j, Fj≤γ4j}, {q2j,Fj>γ4j}}]
(* Compute the specified correlation in Table2 for Distributions 1 and 2 *)
η
j
k
=
2
π
*
NIntegrate
[
X
j
*
F
k
*
f
j
k
,
{
Z
j
,
-
10
,
10
}
,
{
Z
k
,
-
10
,
10
}
,
Method
→
MultiDimensional
]
0.70


## 5. The Procedure for Simulation and Monte Carlo Study

To implement a method for simulating S-D distributions with specifiedL-moments and L-correlations we suggest the following six steps.(1)
Specify theL-moments for T transformations of the forms in (1.1), that is, q1(Φ(Z1)),…,qT(Φ(ZT)) and obtain the solutions for the parameters of γ1j, γ2j, γ3j, and γ4j by solving (3.2)–(3.5) with λ1=0 and λ2=1/π  using the specified values of L-skew (τ3) and L-kurtosis (τ4) for each distribution. Specify a T×T matrix of L-correlations (ηjk) for qj(Φ(Zj)) toward qk(Φ(Zk)), where j<k∈{1,2,…,T}.(2)
Compute the (Pearson) intermediate correlations (ICs)ρjk by substituting the solutions of the parameters from Step (1) into (4.10) and then numerically integrate to solve for ρjk (see Algorithm 1 for an example). Repeat this step separately for all T(T-1)/2 pairwise combinations of correlations.(3)
Assemble the ICs into aT×T 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 generateT standard normal variables (Z1,…,ZT) correlated at the intermediate levels as follows:
(5.1)Z1=a11V1,Z2=a12V1+a22V2,⋮Zj=a1jV1+a2jV2+⋯+aijVi+⋯+ajjVj,⋮ZT=a1TV1+a2TV2+⋯+aiTVi+⋯+ajTVj+⋯+aTTVT,
where V1,…,VT are independent standard normal random variables and where aij represents the element in the ith row and the jth column of the matrix associated with the Cholesky factorization performed in Step (3).(5)
SubstituteZ1,…,ZT from Step (4) into the following Taylor series-based expansion for the standard normal cdf [21]:
(5.2)Φ(Zj)=(12)+ϕ(Zj){Zj+Zj33+Zj5(3·5)+Zj7(3·5·7)+⋯},
where ϕ(Zj) denotes the standard normal pdf and where the absolute error associated with (5.2) is less than 8×10-16.(6)
Substitute the zero-one uniform deviates,Φ(Zj), generated from Step (5) into the T equations of the form of qj(Φ(Zj)) to generate the S-D distributions with the specified L-moments and L-correlations.To demonstrate the steps above, and evaluate the proposed method, a comparison between the proposed L-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 L-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 Z1,…,Z4 with the specified ICs by making use of the formulae given in (5.1) of Step (4) with T=4. The values of Z1,…,Z4 are subsequently transformed to Φ(Z1),…,Φ(Z4) using (5.2) and then substituted into equations of the forms in (1.1) to produce q1(Φ(Z1)),…,q4(Φ(Z4)) for both methods.Table 2
Specified correlation matrix for the distributions in Table1.
1
2
3
4
1
1
2
0.70
1
3
0.60
0.60
1
4
0.35
0.60
0.40
1Table 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
1Table 4
Intermediate correlations for theL-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
1Table 5
Cholesky decompositions for the conventional moment procedure.
a
11
=
1
a
12
=
0.761371
a
13
=
0.945345
a
14
=
0.810720
0
a
22
=
0.648316
a
23
=
0.144858
a
24
=
0.522892
0
0
a
33
=
0.292128
a
34
=
0.232914
0
0
0
a
44
=
0.122753Table 6
Cholesky decompositions for theL-moment procedure.
a
11
=
1
a
12
=
0
.
714645
a
13
=
0.615758
a
14
=
0.362634
0
a
22
=
0.699487
a
23
=
0.207094
a
24
=
0.465703
0
0
a
33
=
0.760233
a
34
=
0.016269
0
0
0
a
44
=
0.807064Algorithm 2:Mathematica source code for computing intermediate correlations for specified conventional Pearson correlations. The example is for distributions j=1 and k=2 (ρ12*) in Table 1. See also Tables 2 and 3.
(* Intermediate Correlation for Distributions 1 and 2. See Table3. *)
ρ
j
k
= 0.761371;
Needs[“MultivariateStatistics”]
f
j
k
=
PDF[MultinormalDistribution[{0,0},{{1,ρjk},{ρjk,1}}],{Zj,Zk}];
F
j
=
CDF[NormalDistribution[0,1],Zj];
F
k
=
CDF[NormalDistribution[0,1],Zk];
(* S-D C-Moment Parameters Distributionsj=1 and k=2 in Table 1*)
γ
1
j
=
-
0.8427;
γ
2
j
=
4.6213;
γ
3
j
=
1.5;
γ
4
j
=
0.25;
γ
1
k
=
0.0;
γ
2
k
=
21.166;
γ
3
k
=
3.0;
γ
4
k
=
0.5;
(* Quantile functions for (1.1) *)
q
1
j
=
γ
1
j
-
γ
2
j
*
(
γ
4
j
-
F
j
)
γ
3
j;
q
2
j
=
γ
1
j
+
γ
2
j
*
(
F
j
-
γ
4
j
)
γ
3
j;
X
j
= Piecewise[{{q1j,Fj≤γ4j},{q2j,Fj>γ4j}}]
q
1
k
=
γ
1
k
-
γ
2
k
*
(
γ
4
k
-
F
k
)
γ
3
k;
q
2
k
=
γ
1
k
+
γ
2
k
*
(
F
k
-
γ
4
k
)
γ
3
k;
X
k
= Piecewise[{{q1k,Fk≤γ4k},{q2k,Fk>γ4k}}]
(* Compute the specified correlation in Table2 for Distributions 1 and 2. *)
ρ
j
k
*
=NIntegrate[Xj*Xk*fjk, {Zj,-10,10}, {Zk,-10,10}, Method→MultiDimensional]
0.70In 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 (α3), kurtosis (α4), and Pearson correlation (ρjk*); (ii) L-skew (τ3), L-kurtosis (τ4), and L-correlation (ηjk). All estimates were based on sample sizes of n=25 and n=1000. The formulae used for computing estimates of α3,4 were based on Fisher’s k-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 α3,4=0 for the standard normal distribution). The formulae used for computing estimates of τ3,4 were Headrick’s equations (2.4) and (2.6) [22]. The estimate for ρjk* was based on the usual formula for the Pearson product-moment of correlation statistic and the estimate for ηjk was computed based on (4.5) using the empirical forms of the cdfs in (4.1) and (4.3). The estimates for ρjk* and ηjk were both transformed using Fisher’s z′ 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 (α3,4, τ3,4, zρjk*′, zηjk′) 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 zρjk*′ and zηjk′ were transformed back to their original metrics (i.e., estimates for ρjk* and ηjk). Further, if a parameter (P) was outside its associated bootstrap C.I., then an index of relative bias (RB) was computed for the estimate (E) as: RB=((E-P)/P)×100. 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 (α3) and kurtosis (α4) results for distributions 2 and 4 in Table 1.
Dist
Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %
n
=
25
2
α
3
=
0.0
0.0019
−0.0055,  0.0094
0.00381
—
α
4
=
0.7692
1.166
1.1496,  1.1821
0.00825
51.59
4
α
3
=
-
4.181
−3.581
−3.5917,  −3.5702
0.00549
−14.35
α
4
=
18.21
13.61
13.5278,  13.6973
0.04337
−25.26
n
=
1000
2
α
3
=
0.0
−0.0005
−0.0013,  0.0004
0.00042
—
α
4
=
0.7692
0.7798
0.7778,  0.7817
0.00100
1.38
4
α
3
=
-
4.181
−4.190
−4.1936,  −4.1866
0.00177
0.22
α
4
=
18.21
18.44
18.4046,  18.4738
0.01775
1.26Table 8
L-skew (τ3) and L-kurtosis (τ4) results for distributions 2 and 4 in Table 1.
Dist
Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %
n
=
25
2
τ
3
=
0.0
0.0002
−0.0014,  0.0018
0.00079
—
τ
4
=
0.2857
0.2981
0.2968,  0.2993
0.00063
4.34
4
τ
3
=
-
0.8414
−0.8638
−0.8647,  −0.8628
0.00049
2.66
τ
4
=
0.6483
0.6956
0.6935,  0.6975
0.00102
7.30
n
=
1000
2
τ
3
=
0.0
−0.0001
−0.0003,  0.0001
0.00010
—
τ
4
=
0.2857
0.2860
0.2859,  0.2862
0.00008
0.11
4
τ
3
=
-
0.8414
−0.8422
−0.8423,-0.8420
0.00007
0.095
τ
4
=
0.6483
0.6498
0.6495,  0.6501
0.00015
0.23Table 9
Correlation results for the conventional moment proceduren=25.
Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %
ρ
12
*
=
0.70
0.7083
0.7071,  0.7095
0.00121
1.19
ρ
13
*
=
0.60
0.6396
0.6380,  0.6414
0.00146
6.60
ρ
14
*
=
0.35
0.3837
0.3826,  0.3848
0.00065
9.63
ρ
23
*
=
0.60
0.6342
0.6327,  0.6356
0.00122
5.70
ρ
24
*
=
0.60
0.6342
0.6330,  0.6353
0.00097
5.70
ρ
34
*
=
0.40
0.4958
0.4919,  0.5000
0.00272
23.95Table 10
Correlation results for theL-moment procedure n=25.
Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %
η
12
=
0.70
0.7155
0.7140,  0.7170
0.00152
2.21
η
13
=
0.60
0.6155
0.6137,  0.6173
0.00148
2.58
η
14
=
0.35
0.3588
0.3562,  0.3611
0.00143
2.51
η
23
=
0.60
0.6116
0.6099,  0.6136
0.00149
1.93
η
24
=
0.60
0.6057
0.6038,  0.6075
0.00151
0.95
η
34
=
0.40
0.4217
0.4177,  0.4259
0.00257
5.43Table 11
Correlation results for the conventional moment proceduren=1000.
Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %
ρ
12
*
=
0.70
0.7001
0.6999,  0.7003
0.00018
—
ρ
13
*
=
0.60
0.6004
0.6001,  0.6006
0.00021
0.07
ρ
14
*
=
0.35
0.3502
0.3500,  0.3504
0.00011
—
ρ
23
*
=
0.60
0.6002
0.5999,  0.6004
0.00019
—
ρ
24
*
=
0.60
0.6004
0.6002,  0.6006
0.00015
0.07
ρ
34
*
=
0.40
0.4009
0.4004,  0.4013
0.00027
0.23Table 12
Correlation results for theL-moment procedure n=1000.
Parameter
Estimate
95% bootstrap C.I.
Standard error
Relative bias %
η
12
=
0.70
0.7002
0.7000,  0.7004
0.00022
—
η
13
=
0.60
0.6005
0.6003,  0.6008
0.00021
0.08
η
14
=
0.35
0.3501
0.3498,  0.3505
0.00021
—
η
23
=
0.60
0.6000
0.5997,  0.6002
0.00022
—
η
24
=
0.60
0.6002
0.5999,  0.6005
0.00021
—
η
34
=
0.40
0.3984
0.3979,  0.3988
0.00035
−0.40


## 6. Discussion and Conclusion

One of the primary advantages thatL-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 L-moment ratios (τ3,τ4) have over their corresponding conventional moment based counterparts (α3,α4) is obvious. For example, with samples of size n=25 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 L-skew and L-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 L-skew and L-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 (n=25) indicates RSE measures of RSE(α^3) = 0.1533% and RSE(α^4) =0.3187% compared with RSE(τ^3) =0.0567% and RSE(τ^4) =0.1466%. This demonstrates that L-skew and L-kurtosis have more precision because they have less variance around their estimates.Presented in Tables9–12 are the results associated with the conventional Pearson and L-correlations. Inspection of Tables 9 and 10 indicates that the L-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, n=25, ρ12*=0.40) the relative bias for Distributions 3 and 4  was 23.95% for the Pearson correlation compared to 5.43% for the L-correlation (Table 10, n=25, η12=0.40). For large sample sizes (Tables 11 and 12, n=1000), both procedures performed adequately as their estimates were in close proximity with their respective population parameters.In summary, the proposedL-moment-based S-D class of distributions is an attractive alternative to the conventional moment-based S-D system. In particular, the L-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 L-moment-based method.



---

*Source: 475781-2012-10-16.xml*