# A Method for Simulating Burr Type III and Type XII Distributions throughL-Moments and L-Correlations

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

---

## Abstract

This paper derives the Burr Type III and Type XII family of distributions in the contexts of univariateL-moments and the L-correlations. Included is the development of a procedure for specifying nonnormal distributions with controlled degrees of L-skew, L-kurtosis, and L-correlations. The procedure can be applied in a variety of settings such as statistical modeling (e.g., forestry, fracture roughness, life testing, operational risk, etc.) and Monte Carlo or simulation studies. Numerical examples are provided to demonstrate that L-moment-based Burr distributions are superior to their conventional moment-based analogs in terms of estimation and distribution fitting. Evaluation of the proposed procedure also demonstrates that the estimates of L-skew, L-kurtosis, and L-correlation are substantially superior to their conventional product moment-based counterparts of skew, kurtosis, and Pearson correlations in terms of relative bias and relative efficiency—most notably when heavy-tailed distributions are of concern.


---

## Body

## 1. Introduction

Burr [1] introduced twelve cumulative distribution functions (cdfs) with the primary purpose of fitting distributions to real-world data. Two popular Burr cdfs are referred to as the Burr Type III and Type XII distributions. The specific forms of these two cdfs are given as [1, Equations (11), (20)]
(1)F(x)(TypeIII)=(1+xc)-k,F(x)(TypeXII)=1-(1+xc)-k,
where x∈(0,∞), c,  and  k are real-valued parameters that determine the mean, variance, skew, and kurtosis of a distribution. The parameter c is negative for the Type III and positive for the Type XII distributions, whereas the parameter k is positive for both the distributions [2].The Burr Type III and Type XII distributions attract special attention because they include several families of nonnormal distributions (e.g., the Gamma distribution) with varying degrees of skew and kurtosis [2–5]. Further, these distributions have been used in a variety of applied mathematics contexts. Some examples include modeling events associated with forestry [6, 7], fracture roughness [8, 9], life testing [10, 11], operational risk [12], option market price distributions [13], meteorology [14], modeling crop prices [15], software reliability growth [16], reliability analysis [17], and in the context of Monte Carlo or simulation studies (e.g., [2]).The quantile functions associated with (1) are expressed as [2, Equations (5), (6)]
(2)q(u)(TypeIII)=F-1(x)=(u-1/k-1)1/c,(3)q(u)(TypeXII)=F-1(x)=((1-u)-1/k-1)1/c,
where u~iidU(0,1) with pdf and cdf as 1 and u, respectively. The shape of a Burr distribution associated with (2) or (3) is contingent on the values of the shape parameters (c and k), which can be determined by simultaneously solving equations (16) and (17) from [2, p. 2211] for given values of skew and kurtosis.In order for (2) or (3) to produce a valid Burr pdf, the quantile function q(u) is required to be a strictly increasing monotone function [2]. This requirement implies that an inverse function (q-1) exists. As such, the cdf associated with (2) or (3) can be expressed as F(q(u))=F(u)=u, and subsequently differentiating this cdf with respect to u will yield the parametric form of the pdf for q(u) as f(q(u))=1/q′(u). We would also note that the simple closed form expressions for the pdfs associated with (1) can be given as [1]
(4)f(x)(TypeIII)=-ckxc-1(1+xc)-(k+1),(5)f(x)(TypeXII)=ckxc-1(1+xc)-(k+1).Some of the problems associated with conventional moment-based estimates are that they can be (a) substantially biased, (b) highly dispersed, or (c) influenced by outliers [18, 19] and thus may not be good representatives of the true parameters. To demonstrate, Figure 1 gives the graphs of the pdf and cdf associated with Burr Type XII distribution with skew (γ3)=3 and kurtosis (γ4)=21. These values of skew and kurtosis have been used in a number of studies [18, 20–22]. Table 1 gives the parameters and sample estimates of skew and kurtosis for the distribution in Figure 1. Inspection of Table 1 indicates that the bootstrap estimates (g3 and g4) of skew and kurtosis (γ3 and γ4) are substantially attenuated below their corresponding parameter values with greater bias and variance as the order of the estimate increases. Specifically, for sample size of n=50, the values of the estimates are only 63.53%, and 22.81% of their corresponding parameters, respectively. The estimates (g3 and g4) of skew and kurtosis (γ3 and γ4) in Table 1 were calculated based on Fisher’s k-statistics formulae (see, e.g., [23, pages 299-300]), currently used by most commercial software packages such as SAS, SPSS, and Minitab, and so forth, for computing the values of skew and kurtosis (where γ3,4=0 for the standard normal distribution).Table 1
Conventional moment-based parameters (γ3,γ4) of skew and kurtosis and their estimates (g3,g4) for the pdf in Figure 1. Each bootstrapped estimate (Estimate), associated 95% bootstrap confidence interval (95% Bootstrap C.I.), and the standard error (St. error) were based on resampling 25,000 statistics. Each statistic was based on a sample size of n=50.
Parameter
Estimate
95% Bootstrap C.I.
St. error
γ
3
=
3.00
g
3
=
1.9060
(1.8955, 1.9163)
0.0053
γ
4
=
21.00
g
4
=
4.7900
(4.7199, 4.8525)
0.0338The pdf and cdf of the Burr Type XII distribution with skew(γ3)=3 and kurtosis (γ4)=21. The solved values of c and k used in (5) and (1) are c=1.187185 and k=5.199500 which are also associated with the parameters in Tables 1 and 2.
(a)
(b)Another unfavorable quality of conventional moment-based estimators of skew and kurtosis is that their values are algebraically bounded by the sample size(n) such that |g3|≤n and g4≤n [18]. This constraint implies that if a researcher wants to simulate non-normal data with kurtosis γ4=21 as in Table 1, and drawing a sample of size n=15 from this population, the largest possible value of the computed estimate (g4) of kurtosis (γ4) is only 15, which is only 71.43% of the parameter.The method ofL-moments introduced by Hosking [19] is an attractive alternative to conventional moments and can be used for describing theoretical probability distributions, fitting distributions to real-world data, estimating parameters, and testing of hypotheses [18, 19, 24, 25]. In these contexts, we note that the L-moment based estimators of L-skew, L-kurtosis, and L-correlation have been introduced to address the limitations associated with conventional moment-based estimators [18, 19, 24–31]. Some qualities of L-moments that make them superior to conventional moments are that they (a) exist for any distribution with finite mean, and their estimates are (b) nearly unbiased for any sample size and less affected from sampling variability, (c) more robust in the presence of outliers in the sample data, and (d) are not algebraically bounded by the sample size [18, 19, 24, 25]. For example, the estimates (𝓉3 and 𝓉4) of L-skew and L-kurtosis (τ3 and τ4) in Table 2 are relatively closer to their respective parameter values with much smaller variance compared to their conventional moment-based counterparts as in Table 1. Inspection of Table 2 shows that for the sample size of n=50, the values of the estimates are on average 97.71% and 97.28% of their corresponding parameters.Table 2
L-moment-based parameters (τ3,τ4) of L-skew and L-kurtosis and their estimates (𝓉3,𝓉4) for the pdf in Figure 1. Each bootstrapped estimate (Estimate), associated 95% bootstrap confidence interval (95% Bootstrap C.I.), and the standard error (St. error) were based on resampling 25,000 statistics. Each statistic was based on a sample size of n=50.
Parameter
Estimate
95% Bootstrap C.I.
St. error
τ
3
=
0.3536
𝓉
3
=
0.3455
(0.3445, 0.3465)
0.0005
τ
4
=
0.2020
𝓉
4
=
0.1965
(0.1955, 0.1975)
0.0005In view of the above, the main purpose of this study is to characterize the Burr Type III and Type XII distributions through the method ofL-moments in order to obviate the problems associated with conventional moment-based estimators. Another purpose of this study is to develop an L-correlation-based methodology to simulate correlated Burr Type III and Type XII distributions. Specifically, in Section 2, a brief introduction to univariate L-moments is provided. The systems of equations associated with the Type III and Type XII distributions are subsequently derived for determining the shape parameters (c and k) for user specified values of L-skew (τ3) and L-kurtosis (τ4). In Section 3, a comparison between conventional and L-moment-based Burr Type III and Type XII distributions is presented in the contexts of estimation and distribution fitting. Numerical examples based on Monte Carlo simulation are also provided to confirm the methodology and demonstrate the advantages that L-moments have over conventional moments. In Section 4, an introduction to the coefficient of L-correlation is provided, and the methodology for solving for intermediate correlations for specified L-correlation structure is subsequently presented. In Section 5, the steps for implementing the proposed L-moment procedure are described for simulating non-normal Burr Type III and Type XII distributions with controlled skew (L-skew), kurtosis (L-kurtosis), and Pearson correlations (L-correlations). Numerical examples and the results of simulation are also provided to confirm the derivations and compare the new procedure with the conventional moment-based procedure. In Section 6, the results of the simulation are discussed.


## 2. Methodology

### 2.1. Theoretical and Empirical Definitions ofL-Moments

L-moments can be expressed as certain linear combinations of probability weighted moments (PWMs). Let  X1,…,Xi,…,Xn be identically and independently distributed random variables each with pdf  f(x),  cdf F(x), and the quantile function F-1(x), then the PWMs are defined as [18, Equation (6)]
(6)βr=∫F-1(x){F(x)}rf(x)dx,
where r=0,1,2,3. The first four L-moments (λi=1,…,4) associated with X can be expressed in simplified forms as [25, pages 20–22]
(7)λ1=β0,λ2=2β1-β0,λ3=6β2-6β1+β0,λ4=20β3-30β2+12β1-β0,
where the coefficients associated with βr=0,…,3 in (7) are obtained from shifted orthogonal Legendre polynomials and are computed as in [25, pages 20–22] or in [18, pages 4–5].The notationsλ1 and λ2 denote the location and scale parameters. Specifically, in the literature of L-moments, λ1 is referred to as the L-location parameter which is equal to the arithmetic mean, and λ2 (>0) is referred to as the L-scale parameter and is one-half of Gini’s coefficient of mean difference [23, pages 47-48]. Dimensionless L-moment ratios are defined as the ratios of higher-order L-moments (i.e., λ3 and λ4) to λ2. Thus, τ3=λ3/λ2 and τ4=λ4/λ2 are, respectively, the indices of L-skew and L-kurtosis. In general, the indices of L-skew and L-kurtosis are bounded in the interval -1<τ3,4<1, and as in conventional moment theory, a symmetric distribution has L-skew equal to zero [18]. The boundary region for L-skew (τ3) and L-kurtosis (τ4) for a continuous distribution is given by the inequality (see, [32])
(8)5τ32-14<τ4<1.EmpiricalL-moments for a sample (of size n) of real-world data are expressed as linear combinations of the unbiased estimators of the PWMs based on sample order statistics X1:n≤X2:n≤⋯≤Xn:n. Specifically, the unbiased estimators of the PWMs are given as [19, pages 113-114]
(9)br=1n∑i=r+1n(i-1)(i-2)⋯(i-r)(n-1)(n-2)⋯(n-r)Xi:n,
where r=0,1,2,3 and b0 is the sample mean. The first four sample L-moments (ℓ1,ℓ2,ℓ3,ℓ4) are obtained by substituting br from (9) instead of βr in (7). The sample L-moment ratios (i.e., L-skew and L-kurtosis) are denoted by 𝓉3 and 𝓉4, where 𝓉3=ℓ3/ℓ2 and 𝓉4=ℓ4/ℓ2.


### 2.2.L-Moments for Burr Type III Distributions

SubstitutingF-1(x)=(u-1/k-1)1/c from (2), F(x)=u, and f(x)=1 in (6), the rth PWM for the Burr Type III distributions is given by
(10)βr=∫01(u-1/k-1)1/curdu.
Equation (10) can also be expressed as
(11)βr=∫01k(u1/k)(k+kr-1/c)-1(1-u1/k)1/c(k-1u(1/k)-1du).
Let u1/k=x. Then, k-1u(1/k)-1du=dx. Substituting in (11), the rth PWM can be expressed as
(12)βr=k∫01x(k+kr-1/c)-1(1-x)(1+1/c)-1dx.
The integral in (12) is a beta function, Beta(α,β), where α=k+kr-1/c and β=1+1/c such that Real(α),Real(β)>0.Integrating (12) for βr=0,1,2,3 and substituting these PWMs into (7) and simplifying gives the following system of equations for the Type III distributions:
(13)λ1=kΓ[k-1/c]Γ[1+1/c]Γ[1+k],(14)λ2=Γ[1+1c](Γ[2k-1/c]Γ[2k]-Γ[k-1/c]Γ[k]),(15)τ3={Γ[1+3k](6Γ[1+k]Γ[2k-1c]-Γ[1+2k]Γ[k-1c])-6Γ[1+k]Γ[1+2k]Γ[3k-1c]}×{Γ[1+3k]×(Γ[1+2k]Γ[k-1c]-2Γ[1+k]Γ[2k-1c])}-1,τ4={Γ[1+4k]×{Γ[1+3k]×(Γ[1+2k]Γ[k-1c]-12Γ[1+k]Γ[2k-1c])+30Γ[1+k]Γ[1+2k]Γ[3k-1c]}-20Γ[1+k]Γ[1+2k]Γ[1+3k]Γ[4k-1c]}×{[2k-1c]Γ[1+3k]×Γ[1+4k]×(Γ[1+2k]Γ[k-1c]-2Γ[1+k]Γ[2k-1c])}-1.


### 2.3.L-Moments for the Burr Type XII Distributions

SubstitutingF-1(x)=((1-u)-1/k-1)1/c from (3), F(x)=u, and f(x)=1 in (6), the rth PWM for the Burr Type XII distributions is given by
(16)βr=∫01((1-u)-1/k-1)1/curdu.
After some manipulations, (16) can be expressed as
(17)βr=∫01-kur((1-u)1/k)(k-1/c)-1(1-(1-u)1/k)1/c×(-k-1(1-u)(1/k)-1du).
Let (1-u)1/k=x. Then, -k-1(1-u)(1/k)-1du=dx. Substituting in (17), the rth PWM can be rewritten as
(18)βr=k∫01(1-xk)rx(k-1/c)-1(1-x)(1+1/c)-1dx.Integrating (18) for βr=0,1,2,3 and substituting these PWMs into (7) and simplifying gives the following system of equations for the Type XII distributions:
(19)λ1=kΓ[k-1/c]Γ[1+1/c]Γ[1+k],(20)λ2=-Γ[1+1c](Γ[2k-1/c]Γ[2k]-Γ[k-1/c]Γ[k]),(21)τ3={Γ[1+3k](-6Γ[1+k]Γ[2k-1c]+Γ[1+2k]Γ[k-1c])+6Γ[1+k]Γ[1+2k]Γ[3k-1c]}×{Γ[1+3k]×(Γ[1+2k]Γ[k-1c]-2Γ[1+k]Γ[2k-1c])}-1,τ4={[3k-1c]Γ[k]Γ[1+4k]×{[3k-1c]Γ[1+3k]×{Γ[1+2k]Γ[k-1c]-12Γ[1+k]Γ[2k-1c]}+30Γ[1+k]Γ[1+2k]Γ[3k-1c]}-40Γ[2k]Γ[1+k]2Γ[1+3k]Γ[4k-1c]}×{[2k-1c]Γ[k]Γ[1+3k]×Γ[1+4k]×{Γ[1+2k]Γ[k-1c]-2Γ[1+k]Γ[2k-1c]}}-1.Given specified values ofL-skew (τ3) and L-kurtosis (τ4) the systems of equations (15) and (21) can be simultaneously solved for real values of c and k. The solved values of c and k can be substituted in (2) and (3), respectively, for generating the Burr Type III and Type XII distributions. Further, the solved values of c and k can be substituted in (13)-(14) and (19)-(20), respectively, for computing the values of L-mean (λ1) and L-scale (λ2) associated with the Type III and Type XII distributions. In the next section, two examples are provided to demonstrate the aforementioned methodology and the advantages that L-moments have over conventional moments in the contexts of estimation and distribution fitting.


## 2.1. Theoretical and Empirical Definitions ofL-Moments

L-moments can be expressed as certain linear combinations of probability weighted moments (PWMs). Let  X1,…,Xi,…,Xn be identically and independently distributed random variables each with pdf  f(x),  cdf F(x), and the quantile function F-1(x), then the PWMs are defined as [18, Equation (6)]
(6)βr=∫F-1(x){F(x)}rf(x)dx,
where r=0,1,2,3. The first four L-moments (λi=1,…,4) associated with X can be expressed in simplified forms as [25, pages 20–22]
(7)λ1=β0,λ2=2β1-β0,λ3=6β2-6β1+β0,λ4=20β3-30β2+12β1-β0,
where the coefficients associated with βr=0,…,3 in (7) are obtained from shifted orthogonal Legendre polynomials and are computed as in [25, pages 20–22] or in [18, pages 4–5].The notationsλ1 and λ2 denote the location and scale parameters. Specifically, in the literature of L-moments, λ1 is referred to as the L-location parameter which is equal to the arithmetic mean, and λ2 (>0) is referred to as the L-scale parameter and is one-half of Gini’s coefficient of mean difference [23, pages 47-48]. Dimensionless L-moment ratios are defined as the ratios of higher-order L-moments (i.e., λ3 and λ4) to λ2. Thus, τ3=λ3/λ2 and τ4=λ4/λ2 are, respectively, the indices of L-skew and L-kurtosis. In general, the indices of L-skew and L-kurtosis are bounded in the interval -1<τ3,4<1, and as in conventional moment theory, a symmetric distribution has L-skew equal to zero [18]. The boundary region for L-skew (τ3) and L-kurtosis (τ4) for a continuous distribution is given by the inequality (see, [32])
(8)5τ32-14<τ4<1.EmpiricalL-moments for a sample (of size n) of real-world data are expressed as linear combinations of the unbiased estimators of the PWMs based on sample order statistics X1:n≤X2:n≤⋯≤Xn:n. Specifically, the unbiased estimators of the PWMs are given as [19, pages 113-114]
(9)br=1n∑i=r+1n(i-1)(i-2)⋯(i-r)(n-1)(n-2)⋯(n-r)Xi:n,
where r=0,1,2,3 and b0 is the sample mean. The first four sample L-moments (ℓ1,ℓ2,ℓ3,ℓ4) are obtained by substituting br from (9) instead of βr in (7). The sample L-moment ratios (i.e., L-skew and L-kurtosis) are denoted by 𝓉3 and 𝓉4, where 𝓉3=ℓ3/ℓ2 and 𝓉4=ℓ4/ℓ2.


## 2.2.L-Moments for Burr Type III Distributions

SubstitutingF-1(x)=(u-1/k-1)1/c from (2), F(x)=u, and f(x)=1 in (6), the rth PWM for the Burr Type III distributions is given by
(10)βr=∫01(u-1/k-1)1/curdu.
Equation (10) can also be expressed as
(11)βr=∫01k(u1/k)(k+kr-1/c)-1(1-u1/k)1/c(k-1u(1/k)-1du).
Let u1/k=x. Then, k-1u(1/k)-1du=dx. Substituting in (11), the rth PWM can be expressed as
(12)βr=k∫01x(k+kr-1/c)-1(1-x)(1+1/c)-1dx.
The integral in (12) is a beta function, Beta(α,β), where α=k+kr-1/c and β=1+1/c such that Real(α),Real(β)>0.Integrating (12) for βr=0,1,2,3 and substituting these PWMs into (7) and simplifying gives the following system of equations for the Type III distributions:
(13)λ1=kΓ[k-1/c]Γ[1+1/c]Γ[1+k],(14)λ2=Γ[1+1c](Γ[2k-1/c]Γ[2k]-Γ[k-1/c]Γ[k]),(15)τ3={Γ[1+3k](6Γ[1+k]Γ[2k-1c]-Γ[1+2k]Γ[k-1c])-6Γ[1+k]Γ[1+2k]Γ[3k-1c]}×{Γ[1+3k]×(Γ[1+2k]Γ[k-1c]-2Γ[1+k]Γ[2k-1c])}-1,τ4={Γ[1+4k]×{Γ[1+3k]×(Γ[1+2k]Γ[k-1c]-12Γ[1+k]Γ[2k-1c])+30Γ[1+k]Γ[1+2k]Γ[3k-1c]}-20Γ[1+k]Γ[1+2k]Γ[1+3k]Γ[4k-1c]}×{[2k-1c]Γ[1+3k]×Γ[1+4k]×(Γ[1+2k]Γ[k-1c]-2Γ[1+k]Γ[2k-1c])}-1.


## 2.3.L-Moments for the Burr Type XII Distributions

SubstitutingF-1(x)=((1-u)-1/k-1)1/c from (3), F(x)=u, and f(x)=1 in (6), the rth PWM for the Burr Type XII distributions is given by
(16)βr=∫01((1-u)-1/k-1)1/curdu.
After some manipulations, (16) can be expressed as
(17)βr=∫01-kur((1-u)1/k)(k-1/c)-1(1-(1-u)1/k)1/c×(-k-1(1-u)(1/k)-1du).
Let (1-u)1/k=x. Then, -k-1(1-u)(1/k)-1du=dx. Substituting in (17), the rth PWM can be rewritten as
(18)βr=k∫01(1-xk)rx(k-1/c)-1(1-x)(1+1/c)-1dx.Integrating (18) for βr=0,1,2,3 and substituting these PWMs into (7) and simplifying gives the following system of equations for the Type XII distributions:
(19)λ1=kΓ[k-1/c]Γ[1+1/c]Γ[1+k],(20)λ2=-Γ[1+1c](Γ[2k-1/c]Γ[2k]-Γ[k-1/c]Γ[k]),(21)τ3={Γ[1+3k](-6Γ[1+k]Γ[2k-1c]+Γ[1+2k]Γ[k-1c])+6Γ[1+k]Γ[1+2k]Γ[3k-1c]}×{Γ[1+3k]×(Γ[1+2k]Γ[k-1c]-2Γ[1+k]Γ[2k-1c])}-1,τ4={[3k-1c]Γ[k]Γ[1+4k]×{[3k-1c]Γ[1+3k]×{Γ[1+2k]Γ[k-1c]-12Γ[1+k]Γ[2k-1c]}+30Γ[1+k]Γ[1+2k]Γ[3k-1c]}-40Γ[2k]Γ[1+k]2Γ[1+3k]Γ[4k-1c]}×{[2k-1c]Γ[k]Γ[1+3k]×Γ[1+4k]×{Γ[1+2k]Γ[k-1c]-2Γ[1+k]Γ[2k-1c]}}-1.Given specified values ofL-skew (τ3) and L-kurtosis (τ4) the systems of equations (15) and (21) can be simultaneously solved for real values of c and k. The solved values of c and k can be substituted in (2) and (3), respectively, for generating the Burr Type III and Type XII distributions. Further, the solved values of c and k can be substituted in (13)-(14) and (19)-(20), respectively, for computing the values of L-mean (λ1) and L-scale (λ2) associated with the Type III and Type XII distributions. In the next section, two examples are provided to demonstrate the aforementioned methodology and the advantages that L-moments have over conventional moments in the contexts of estimation and distribution fitting.


## 3. Comparison ofL-Moments with Conventional Moments

### 3.1. Estimation

An example to demonstrate the advantages ofL-moment-based estimation over conventional moment-based estimation is provided in Figure 2 and Tables 4–7. Given in Figure 2 are the pdfs of the F distribution (3, 10), Chi-square (df=1), extreme value (0, 1), and logistic (0, 1) distributions superimposed, respectively, by the Burr Type XII, Type III, Type XII, and Type III pdfs (dashed curves) in both (a) conventional moment- and (b) L-moment-based systems. The conventional moment-based parameters of skew (γ3) and kurtosis (γ4) associated with these four distributions, given in Table 4, were computed by using equations (11)−(13) from [2, page 2211]. The values of shape parameters (c and k) given in Table 4 were determined by solving equations (16) and (17) from [2, page 2211]. The values of c and k were used in (4) and (5) to superimpose the conventional moment-based Burr Type XII, Type III, Type XII, and Type III distributions as shown in  Figure 2(a).The pdfs (dashed curves) of the four distributions: Distribution 1∶= Burr Type XII ≈F distribution (3, 10), Distribution 2 ∶= Burr Type III ≈ Chi square (df= 1), Distribution 3 ∶=   Burr Type XII ≈ Extreme Value (0, 1), and Distribution 4 ∶= Burr Type III ≈ Logistic (0, 1) superimposed by the (a) conventional moment- and (b) L-moment-based Burr distributions.
(a)
(b)TheL-moment-based parameters of L-skew (τ3) and L-kurtosis (τ4) associated with the four distributions in Figure 2, given in Table 5, were obtained in three steps as follows: (a) compute the values of PWMs (βr=0,1,2,3) using (6), (b) substitute these PWMs into (7) to obtain the values of the first four L-moments, and (c) compute the values of τ3 and τ4 using τ3=λ3/λ2 and τ4=λ4/λ2. The values of shape parameters (c and k) given in Table 5 were determined by solving the systems of equations (15) and (21), respectively. These values of c and k were used in (4) and (5) to superimpose the L-moment-based Burr Type XII, Type III, Type XII, and Type III distributions as shown in Figure 2(b).To superimpose the Type III or Type XII distribution, the quantile functionq(u) in (2) or (3) was transformed into (a) (X-σ-μS+Sq(u))/σ, and (b) (ℓ1λ2-λ1ℓ2+ℓ2q(u))/λ2, respectively, where (X-, S) and (μ, σ) are the values of mean, standard deviation, whereas (ℓ1, ℓ2) and (λ1, λ2) are the values of L-mean, L-scale obtained from the original distribution and the respective Burr Type III or Type XII  approximation, respectively.The advantages ofL-moment-based estimators over those based on conventional moments can also be demonstrated in the context of Burr Type III and Type XII distributions by considering the Monte Carlo simulation results associated with the indices for the percentage of relative bias (RB%) and standard error (St. error) reported in Tables 6 and 7.Specifically, a Fortran [33] algorithm was written to simulate 25,000 independent samples of sizes n=25 and n=1000, and the conventional moment-based estimates (g3 and g4) of skew and kurtosis (γ3 and   γ4) and the L-moment based estimates (𝓉3 and  𝓉4) of L-skew and L-kurtosis (τ3 and τ4) were computed for each of the (2×25,000) samples based on the parameters and the values of c and k listed in Tables 4 and 5. The estimates (g3 and g4) of γ3 and γ4 were computed based on Fisher’s k-statistics formulae [23, pages 47-48], whereas the estimates (𝓉3 and 𝓉4) of τ3 and τ4 were computed using (7) and (9). Bias-corrected accelerated bootstrapped average estimates (Estimate), associated 95% confidence intervals (95% Bootstrap C.I.), and standard errors (St. error) were obtained for each type of estimates using 10,000 resamples via the commercial software package Spotfire S+ [34]. Further, if a parameter was outside its associated 95% bootstrap C.I., then the percentage of relative bias (RB%) was computed for the estimate as
(22)RB%=100×(Estimate-Parameter)Parameter.The results in Tables6 and 7 illustrate that the L-moment-based estimators are superior to their conventional moment-based counterparts in terms of both smaller relative bias and error. These advantages are most pronounced in the context of smaller sample sizes and higher order moments. For example for the Distribution 1, given a sample of size n=25, the conventional moment-based estimates (g3 and g4) generated in the simulation were, on average, 41.46% and 5.62% of their corresponding parameters (γ3 and γ4). On the other hand, for the same Distribution 1, the L-moment-based estimates (𝓉3 and 𝓉4) generated in the simulation study were, on average, 93.88% and 92.03% of their corresponding parameters (τ3 and τ4). Thus, the relative biases of estimators based on L-moments are essentially negligible compared to those associated with the estimators based on conventional moments. Also, it can be verified that the standard errors associated with the estimates 𝓉3 and 𝓉4 are relatively much smaller and more stable than the standard errors associated with the estimates g3 and g4.Inspection of the graphs in Figures2(a) and 2(b)  and the Monte Carlo simulation results in Tables 6 and 7 indicate that the L-moment-based Type III and Type XII pdfs provide a more accurate approximation of the four distributions than those based on conventional moment theory.


### 3.2. Distribution Fitting

Figure3 shows the conventional moment- and the L-moment-based Burr Type XII pdfs superimposed on the histogram of ankle circumference data obtained from 252 adult males (http://lib.stat.cmu.edu/datasets/bodyfat) as cited in Headrick [35, page 48].Histograms of the ankle circumference data superimposed by (a) conventional moment- and (b)L-moment-based Burr Type XII distributions. To superimpose the Burr Type XII distribution (dashed curves), the quantile function q(u) from (3) was transformed as (a) (X-σ-μS+Sq(u))/σ, and (b) (ℓ1λ2-λ1ℓ2+ℓ2q(u))/λ2, respectively, where (X-, S) and (μ, σ) are the values of mean and standard deviation, whereas (ℓ1, ℓ2) and (λ1, λ2) are the values of L-mean and L-scale obtained from the actual data and the Burr Type XII distribution, respectively.
(a)
(b)The conventional moment-based estimates (g3 and g4) of skew and kurtosis (γ3 and γ4)  and the L-moment-based estimates (𝓉3 and 𝓉4) of  L-skew and L-kurtosis (τ3 and τ4) were computed for the sample of size n=252 participants. The estimates of γ3 and γ4 were computed based on Fisher’s k-statistics formulae [23, pages 47-48], whereas the estimates of τ3 and τ4 were computed using (7) and (9), respectively. These sample estimates were then used to solve for the values of shape parameters (candk) using (a) equations (16)−(17) from [2, page 2211] and (b) (21). The solved values of candk were subsequently used in (5) to superimpose the parametric plots of the Burr Type XII pdfs shown in Figure 3.Inspection of the two panels in Figure3 illustrates that the L-moment-based Type XII pdf provides a better fit to the sample data. The Chi-square goodness of fit statistics along with their corresponding p values given in Table 3 provide evidence that the conventional moment-based Type XII distribution does not provide a good fit to the actual data, whereas the L-moment-based Type XII distribution fits very well. The degrees of freedom for the Chi-square goodness of fit tests were computed as df=5=10 (class intervals) − 4 (estimates of the parameters) − 1 (sample size).Table 3
Chi-square goodness of fit statistics for the conventional(C) moment- and the L-moment-(L) based Burr Type XII approximations for the ankle circumference (n=252) data in Figure 3.
Percent
Expected
Obs. (C)
Obs. (L)
Ankle  circumference (C)
Ankle  circumference (L)
10
25.2
22
18
<21.4866
<21.3760
20
25.2
15
29
21.4866–21.7890
21.3760–21.8701
30
25.2
29
30
21.7890–22.0747
21.8701–22.2360
40
25.2
18
26
22.0747–22.3687
22.2360–22.5606
50
25.2
30
25
22.3687–22.6888
22.5606–22.8795
60
25.2
27
27
22.6888–23.0571
22.8795–23.2201
70
25.2
30
20
23.0571–23.5107
23.2201–23.6175
80
25.2
29
25
23.5107–24.1310
23.6175–24.1426
90
25.2
33
32
24.1310–25.1858
24.1426–25.0244
100
25.2
19
20
25.1858 or more
25.0244 or more
χ
2
=
13.6349
χ
2
=
7.6825
p
=
0.0181
p
=
0.1746Table 4
Conventional moment based parameters of the mean(μ), standard deviation (σ), skew (γ3), and kurtosis (γ4) along with their corresponding values of shape parameters (c and k) for the four distributions (dashed curves) in Figure 2 (Panel a).
Distribution
μ
σ
γ
3
γ
4
c
k
1
0.306182
0.362801
4.221159
59.454545
1.097740
4.481719
2
0.139427
0.263634
2.828427
12.0
−5.542062
0.027237
3
0.415523
0.242125
1.139547
2.4
2.025127
5.374401
4
0.937510
0.139248
0.0
1.2
−15.776368
0.540152Table 5
L-moment-based parameters of L-mean (λ1), L-scale (λ2), L-skew (τ3), and L-kurtosis (τ4) along with their corresponding values of shape parameters (c and k) for the four distributions (dashed curves) in Figure 2 (Panel b).
Distribution
λ
1
λ
2
τ
3
τ
4
c
k
1
0.441134
0.221816
0.387355
0.233548
1.227122
3.417227
2
0.349660
0.226233
0.464102
0.226497
−2.944030
0.126960
3
0.643787
0.182777
0.169925
0.150375
2.498521
2.913452
4
0.936479
0.085223
0.0
0.166667
−13.875566
0.562179Table 6
Skew (γ3) and kurtosis (γ4) results for the conventional moment procedure.
Distribution
Parameter
Estimate
95% Bootstrap C.I.
St. error
RB%
n
=
25
1
γ
3
=
4.2212
g
3
=
1.75
1.7397, 1.7597
0.00510
-
58.54
γ
4
=
59.4545
g
4
=
3.34
3.2911, 3.3871
0.02430
-
94.38
2
γ
3
=
2.8284
g
3
=
2.151
2.1428, 2.1596
0.00429
-
23.95
γ
4
=
12.0
g
4
=
4.215
4.1715, 4.2614
0.02306
-
64.88
3
γ
3
=
1.1395
g
3
=
0.7929
0.7859, 0.7994
0.00346
-
30.42
γ
4
=
2.4
g
4
=
0.5014
0.4786, 0.5220
0.01111
-
79.11
4
γ
3
=
0.0
g
3
=
-
0.051
−0.0585, −0.0430
0.00391
—
γ
4
=
1.2
g
4
=
0.3396
0.3243, 0.3552
0.00791
-
71.70
n
=
1000
1
γ
3
=
4.2212
g
3
=
3.606
3.5888, 3.6275
0.00972
-
14.57
γ
4
=
59.4545
g
4
=
26.29
25.8657, 26.7366
0.22370
-
55.78
2
γ
3
=
2.8284
g
3
=
2.761
2.7548, 2.7681
0.00335
-
2.38
γ
4
=
12.0
g
4
=
10.17
10.0615, 10.3107
0.06265
-
15.25
3
γ
3
=
1.1395
g
3
=
1.124
1.1219, 1.1266
0.00122
-
1.36
γ
4
=
2.4
g
4
=
2.266
2.2451, 2.2861
0.01048
-
5.58
4
γ
3
=
0.0
g
3
=
-
0.0045
−0.0064, −0.0025
0.00101
—
γ
4
=
1.2
g
4
=
1.158
1.1487, 1.1683
0.00498
-
3.5Table 7
L-skew (τ3) and L-kurtosis (τ4) results for the L-moment procedure.
Distribution
Parameter
Estimate
95% Bootstrap C.I.
St. error
RB%
n
=
25
1
τ
3
=
0.3874
𝓉
3
=
0.3637
0.3621, 0.3652
0.00078
-
6.12
τ
4
=
0.2335
𝓉
4
=
0.2149
0.2134, 0.2164
0.00078
-
7.97
2
τ
3
=
0.4641
𝓉
3
=
0.4446
0.4431, 0.4461
0.00076
-
4.20
τ
4
=
0.2265
𝓉
4
=
0.2024
0.2007, 0.2042
0.00087
-
10.64
3
τ
3
=
0.1699
𝓉
3
=
0.1614
0.1601, 0.1627
0.00066
-
5.00
τ
4
=
0.1504
𝓉
4
=
0.1472
0.1462, 0.1483
0.00054
-
2.13
4
τ
3
=
0.0
𝓉
3
=
-
0.0019
-
0.0033
,
-
0.0005
0.00071
—
τ
4
=
0.1667
𝓉
4
=
0.1653
0.1643, 0.1663
0.00051
-
0.84
n
=
1000
1
τ
3
=
0.3874
𝓉
3
=
0.3865
0.3863, 0.3868
0.00014
-
0.23
τ
4
=
0.2335
𝓉
4
=
0.2328
0.2325, 0.2331
0.00015
-
0.30
2
τ
3
=
0.4641
𝓉
3
=
0.4633
0.4630, 0.4636
0.00015
-
0.17
τ
4
=
0.2265
𝓉
4
=
0.2255
0.2251, 0.2258
0.00019
-
0.44
3
τ
3
=
0.1699
𝓉
3
=
0.1696
0.1694, 0.1698
0.00010
-
0.18
τ
4
=
0.1504
𝓉
4
=
0.1502
0.1501, 0.1504
0.00008
—
4
τ
3
=
0.0
𝓉
3
=
-
0.0002
-
0.0004
,
0.00002
0.00011
—
τ
4
=
0.1667
𝓉
4
=
0.1667
0.1665, 0.1668
0.00007
—


## 3.1. Estimation

An example to demonstrate the advantages ofL-moment-based estimation over conventional moment-based estimation is provided in Figure 2 and Tables 4–7. Given in Figure 2 are the pdfs of the F distribution (3, 10), Chi-square (df=1), extreme value (0, 1), and logistic (0, 1) distributions superimposed, respectively, by the Burr Type XII, Type III, Type XII, and Type III pdfs (dashed curves) in both (a) conventional moment- and (b) L-moment-based systems. The conventional moment-based parameters of skew (γ3) and kurtosis (γ4) associated with these four distributions, given in Table 4, were computed by using equations (11)−(13) from [2, page 2211]. The values of shape parameters (c and k) given in Table 4 were determined by solving equations (16) and (17) from [2, page 2211]. The values of c and k were used in (4) and (5) to superimpose the conventional moment-based Burr Type XII, Type III, Type XII, and Type III distributions as shown in  Figure 2(a).The pdfs (dashed curves) of the four distributions: Distribution 1∶= Burr Type XII ≈F distribution (3, 10), Distribution 2 ∶= Burr Type III ≈ Chi square (df= 1), Distribution 3 ∶=   Burr Type XII ≈ Extreme Value (0, 1), and Distribution 4 ∶= Burr Type III ≈ Logistic (0, 1) superimposed by the (a) conventional moment- and (b) L-moment-based Burr distributions.
(a)
(b)TheL-moment-based parameters of L-skew (τ3) and L-kurtosis (τ4) associated with the four distributions in Figure 2, given in Table 5, were obtained in three steps as follows: (a) compute the values of PWMs (βr=0,1,2,3) using (6), (b) substitute these PWMs into (7) to obtain the values of the first four L-moments, and (c) compute the values of τ3 and τ4 using τ3=λ3/λ2 and τ4=λ4/λ2. The values of shape parameters (c and k) given in Table 5 were determined by solving the systems of equations (15) and (21), respectively. These values of c and k were used in (4) and (5) to superimpose the L-moment-based Burr Type XII, Type III, Type XII, and Type III distributions as shown in Figure 2(b).To superimpose the Type III or Type XII distribution, the quantile functionq(u) in (2) or (3) was transformed into (a) (X-σ-μS+Sq(u))/σ, and (b) (ℓ1λ2-λ1ℓ2+ℓ2q(u))/λ2, respectively, where (X-, S) and (μ, σ) are the values of mean, standard deviation, whereas (ℓ1, ℓ2) and (λ1, λ2) are the values of L-mean, L-scale obtained from the original distribution and the respective Burr Type III or Type XII  approximation, respectively.The advantages ofL-moment-based estimators over those based on conventional moments can also be demonstrated in the context of Burr Type III and Type XII distributions by considering the Monte Carlo simulation results associated with the indices for the percentage of relative bias (RB%) and standard error (St. error) reported in Tables 6 and 7.Specifically, a Fortran [33] algorithm was written to simulate 25,000 independent samples of sizes n=25 and n=1000, and the conventional moment-based estimates (g3 and g4) of skew and kurtosis (γ3 and   γ4) and the L-moment based estimates (𝓉3 and  𝓉4) of L-skew and L-kurtosis (τ3 and τ4) were computed for each of the (2×25,000) samples based on the parameters and the values of c and k listed in Tables 4 and 5. The estimates (g3 and g4) of γ3 and γ4 were computed based on Fisher’s k-statistics formulae [23, pages 47-48], whereas the estimates (𝓉3 and 𝓉4) of τ3 and τ4 were computed using (7) and (9). Bias-corrected accelerated bootstrapped average estimates (Estimate), associated 95% confidence intervals (95% Bootstrap C.I.), and standard errors (St. error) were obtained for each type of estimates using 10,000 resamples via the commercial software package Spotfire S+ [34]. Further, if a parameter was outside its associated 95% bootstrap C.I., then the percentage of relative bias (RB%) was computed for the estimate as
(22)RB%=100×(Estimate-Parameter)Parameter.The results in Tables6 and 7 illustrate that the L-moment-based estimators are superior to their conventional moment-based counterparts in terms of both smaller relative bias and error. These advantages are most pronounced in the context of smaller sample sizes and higher order moments. For example for the Distribution 1, given a sample of size n=25, the conventional moment-based estimates (g3 and g4) generated in the simulation were, on average, 41.46% and 5.62% of their corresponding parameters (γ3 and γ4). On the other hand, for the same Distribution 1, the L-moment-based estimates (𝓉3 and 𝓉4) generated in the simulation study were, on average, 93.88% and 92.03% of their corresponding parameters (τ3 and τ4). Thus, the relative biases of estimators based on L-moments are essentially negligible compared to those associated with the estimators based on conventional moments. Also, it can be verified that the standard errors associated with the estimates 𝓉3 and 𝓉4 are relatively much smaller and more stable than the standard errors associated with the estimates g3 and g4.Inspection of the graphs in Figures2(a) and 2(b)  and the Monte Carlo simulation results in Tables 6 and 7 indicate that the L-moment-based Type III and Type XII pdfs provide a more accurate approximation of the four distributions than those based on conventional moment theory.


## 3.2. Distribution Fitting

Figure3 shows the conventional moment- and the L-moment-based Burr Type XII pdfs superimposed on the histogram of ankle circumference data obtained from 252 adult males (http://lib.stat.cmu.edu/datasets/bodyfat) as cited in Headrick [35, page 48].Histograms of the ankle circumference data superimposed by (a) conventional moment- and (b)L-moment-based Burr Type XII distributions. To superimpose the Burr Type XII distribution (dashed curves), the quantile function q(u) from (3) was transformed as (a) (X-σ-μS+Sq(u))/σ, and (b) (ℓ1λ2-λ1ℓ2+ℓ2q(u))/λ2, respectively, where (X-, S) and (μ, σ) are the values of mean and standard deviation, whereas (ℓ1, ℓ2) and (λ1, λ2) are the values of L-mean and L-scale obtained from the actual data and the Burr Type XII distribution, respectively.
(a)
(b)The conventional moment-based estimates (g3 and g4) of skew and kurtosis (γ3 and γ4)  and the L-moment-based estimates (𝓉3 and 𝓉4) of  L-skew and L-kurtosis (τ3 and τ4) were computed for the sample of size n=252 participants. The estimates of γ3 and γ4 were computed based on Fisher’s k-statistics formulae [23, pages 47-48], whereas the estimates of τ3 and τ4 were computed using (7) and (9), respectively. These sample estimates were then used to solve for the values of shape parameters (candk) using (a) equations (16)−(17) from [2, page 2211] and (b) (21). The solved values of candk were subsequently used in (5) to superimpose the parametric plots of the Burr Type XII pdfs shown in Figure 3.Inspection of the two panels in Figure3 illustrates that the L-moment-based Type XII pdf provides a better fit to the sample data. The Chi-square goodness of fit statistics along with their corresponding p values given in Table 3 provide evidence that the conventional moment-based Type XII distribution does not provide a good fit to the actual data, whereas the L-moment-based Type XII distribution fits very well. The degrees of freedom for the Chi-square goodness of fit tests were computed as df=5=10 (class intervals) − 4 (estimates of the parameters) − 1 (sample size).Table 3
Chi-square goodness of fit statistics for the conventional(C) moment- and the L-moment-(L) based Burr Type XII approximations for the ankle circumference (n=252) data in Figure 3.
Percent
Expected
Obs. (C)
Obs. (L)
Ankle  circumference (C)
Ankle  circumference (L)
10
25.2
22
18
<21.4866
<21.3760
20
25.2
15
29
21.4866–21.7890
21.3760–21.8701
30
25.2
29
30
21.7890–22.0747
21.8701–22.2360
40
25.2
18
26
22.0747–22.3687
22.2360–22.5606
50
25.2
30
25
22.3687–22.6888
22.5606–22.8795
60
25.2
27
27
22.6888–23.0571
22.8795–23.2201
70
25.2
30
20
23.0571–23.5107
23.2201–23.6175
80
25.2
29
25
23.5107–24.1310
23.6175–24.1426
90
25.2
33
32
24.1310–25.1858
24.1426–25.0244
100
25.2
19
20
25.1858 or more
25.0244 or more
χ
2
=
13.6349
χ
2
=
7.6825
p
=
0.0181
p
=
0.1746Table 4
Conventional moment based parameters of the mean(μ), standard deviation (σ), skew (γ3), and kurtosis (γ4) along with their corresponding values of shape parameters (c and k) for the four distributions (dashed curves) in Figure 2 (Panel a).
Distribution
μ
σ
γ
3
γ
4
c
k
1
0.306182
0.362801
4.221159
59.454545
1.097740
4.481719
2
0.139427
0.263634
2.828427
12.0
−5.542062
0.027237
3
0.415523
0.242125
1.139547
2.4
2.025127
5.374401
4
0.937510
0.139248
0.0
1.2
−15.776368
0.540152Table 5
L-moment-based parameters of L-mean (λ1), L-scale (λ2), L-skew (τ3), and L-kurtosis (τ4) along with their corresponding values of shape parameters (c and k) for the four distributions (dashed curves) in Figure 2 (Panel b).
Distribution
λ
1
λ
2
τ
3
τ
4
c
k
1
0.441134
0.221816
0.387355
0.233548
1.227122
3.417227
2
0.349660
0.226233
0.464102
0.226497
−2.944030
0.126960
3
0.643787
0.182777
0.169925
0.150375
2.498521
2.913452
4
0.936479
0.085223
0.0
0.166667
−13.875566
0.562179Table 6
Skew (γ3) and kurtosis (γ4) results for the conventional moment procedure.
Distribution
Parameter
Estimate
95% Bootstrap C.I.
St. error
RB%
n
=
25
1
γ
3
=
4.2212
g
3
=
1.75
1.7397, 1.7597
0.00510
-
58.54
γ
4
=
59.4545
g
4
=
3.34
3.2911, 3.3871
0.02430
-
94.38
2
γ
3
=
2.8284
g
3
=
2.151
2.1428, 2.1596
0.00429
-
23.95
γ
4
=
12.0
g
4
=
4.215
4.1715, 4.2614
0.02306
-
64.88
3
γ
3
=
1.1395
g
3
=
0.7929
0.7859, 0.7994
0.00346
-
30.42
γ
4
=
2.4
g
4
=
0.5014
0.4786, 0.5220
0.01111
-
79.11
4
γ
3
=
0.0
g
3
=
-
0.051
−0.0585, −0.0430
0.00391
—
γ
4
=
1.2
g
4
=
0.3396
0.3243, 0.3552
0.00791
-
71.70
n
=
1000
1
γ
3
=
4.2212
g
3
=
3.606
3.5888, 3.6275
0.00972
-
14.57
γ
4
=
59.4545
g
4
=
26.29
25.8657, 26.7366
0.22370
-
55.78
2
γ
3
=
2.8284
g
3
=
2.761
2.7548, 2.7681
0.00335
-
2.38
γ
4
=
12.0
g
4
=
10.17
10.0615, 10.3107
0.06265
-
15.25
3
γ
3
=
1.1395
g
3
=
1.124
1.1219, 1.1266
0.00122
-
1.36
γ
4
=
2.4
g
4
=
2.266
2.2451, 2.2861
0.01048
-
5.58
4
γ
3
=
0.0
g
3
=
-
0.0045
−0.0064, −0.0025
0.00101
—
γ
4
=
1.2
g
4
=
1.158
1.1487, 1.1683
0.00498
-
3.5Table 7
L-skew (τ3) and L-kurtosis (τ4) results for the L-moment procedure.
Distribution
Parameter
Estimate
95% Bootstrap C.I.
St. error
RB%
n
=
25
1
τ
3
=
0.3874
𝓉
3
=
0.3637
0.3621, 0.3652
0.00078
-
6.12
τ
4
=
0.2335
𝓉
4
=
0.2149
0.2134, 0.2164
0.00078
-
7.97
2
τ
3
=
0.4641
𝓉
3
=
0.4446
0.4431, 0.4461
0.00076
-
4.20
τ
4
=
0.2265
𝓉
4
=
0.2024
0.2007, 0.2042
0.00087
-
10.64
3
τ
3
=
0.1699
𝓉
3
=
0.1614
0.1601, 0.1627
0.00066
-
5.00
τ
4
=
0.1504
𝓉
4
=
0.1472
0.1462, 0.1483
0.00054
-
2.13
4
τ
3
=
0.0
𝓉
3
=
-
0.0019
-
0.0033
,
-
0.0005
0.00071
—
τ
4
=
0.1667
𝓉
4
=
0.1653
0.1643, 0.1663
0.00051
-
0.84
n
=
1000
1
τ
3
=
0.3874
𝓉
3
=
0.3865
0.3863, 0.3868
0.00014
-
0.23
τ
4
=
0.2335
𝓉
4
=
0.2328
0.2325, 0.2331
0.00015
-
0.30
2
τ
3
=
0.4641
𝓉
3
=
0.4633
0.4630, 0.4636
0.00015
-
0.17
τ
4
=
0.2265
𝓉
4
=
0.2255
0.2251, 0.2258
0.00019
-
0.44
3
τ
3
=
0.1699
𝓉
3
=
0.1696
0.1694, 0.1698
0.00010
-
0.18
τ
4
=
0.1504
𝓉
4
=
0.1502
0.1501, 0.1504
0.00008
—
4
τ
3
=
0.0
𝓉
3
=
-
0.0002
-
0.0004
,
0.00002
0.00011
—
τ
4
=
0.1667
𝓉
4
=
0.1667
0.1665, 0.1668
0.00007
—


## 4.L-Correlations for the Burr Type III and Type XII Distributions

LetYj and Yk be two random variables with cdfs F(Yj) and F(Yk), respectively. The second L-moments of Yj and Yk can be defined as [31]
(23)λ2(Yj)=2Cov(Yj,F(Yj)),λ2(Yk)=2Cov(Yk,F(Yk)).
The second L-comoments of Yj toward Yk and Yk toward Yj are given as
(24)λ2(Yj,Yk)=2Cov(Yj,F(Yk)),λ2(Yk,Yj)=2Cov(Yk,F(Yj)).
The L-correlations of Yj toward Yk and Yk toward Yj are subsequently defined as
(25)ηjk=λ2(Yj,Yk)λ2(Yj),(26)ηkj=λ2(Yk,Yj)λ2(Yk).
The L-correlation given in (25) (or, (26)) is bounded in the interval -1≤ηjk≤1. A value of ηjk=1(orηjk=-1) implies that Yj and Yk have a strictly and monotonically increasing (or decreasing) relationship. See Serfling and Xiao [31] for further details on the topics related to the L-correlation.The extension of the Burr Type III and Type XII distributions to multivariate level can be obtained by specifyingT  quantile functions as given in (2) and (3) with a specified L-correlation structure. Specifically, let Z1,…,ZT denote standard normal variables with cdfs and the joint pdf associated with Zj and Zk given by the following expressions:(27)Φ(Zj)=∫-∞zj(2π)-1/2exp{-vj22}dvj,(28)Φ(Zk)=∫-∞zk(2π)-1/2exp{-vk22}dvk,(29)fjk=(2π(1-rjk2)1/2)-1exp{-(2(1-rjk2))-1×(zj2+zk2-2rjkzjzk)(2(1-rjk2))-1},
where rjk in (29) is the intermediate correlation (IC) between Zj and Zk. Using the cdfs in (27) and (28) as zero-one uniform deviates, that is, Φ(Zj),   Φ(Zk)~U(0,1), the quantile function defined in either (2) or (3) can be expressed as a function of Φ(Zj),  or Φ(Zk) (e.g.,  qj(Φ(Zj)) or  qk(Φ(Zk))). Thus, the L-correlation of Yj=qj(Φ(Zj)) toward Yk=qk(Φ(Zk)) can be determined using (25) with the denominator standardized to λ2(Yj)=1/π for the standard normal distribution as
(30)ηjk=2π∬-∞∞xj(qj(Φ(Zj)))Φ(zk)fjkdzjdzk.
The variable xj(qj(Φ(Zj))) in (30) is the standardized quantile function of (2) or (3) such that it has an L-mean (or arithmetic mean) of zero and L-scale equal to that of the standard normal distribution. That is, the quantile function Yj=qj(Φ(Zj)) is standardized by a linear transformation as
(31)xj(qj(Φ(Zj)))=δ(qj(Φ(Zj))-λ1),
where λ1 is the mean from (13) or (19) and δ is a constant that scales λ2 in (14) or (20) and in the denominator of (25) to 1/π. In particular, δ for the Type III and Type XII distributions can be expressed as
(32)δ(TypeIII)=-Γ[k]Γ[2k]π{Γ[1+1/c](Γ[2k]Γ[k-1/c]-Γ[k]Γ[2k-1/c])},(33)δ(TypeXII)=Γ[k]Γ[2k]π{Γ[1+1/c](Γ[2k]Γ[k-1/c]-Γ[k]Γ[2k-1/c])}.
The next step is to use (30) to solve for the values of the T(T-1)/2  ICs (rjk) such that the T specified Type III and Type XII distributions have their specified L-correlation structure.Analogously, theL-correlation of Yk=qk(Φ(Zk)) toward Yj=qj(Φ(Zj)) is given as
(34)ηkj=2π∬-∞∞xk(qk(Φ(Zk)))Φ(zj)fjkdzkdzj.
Note that in general, the L-correlation of Yj=qj(Φ(Zj)) toward Yk=qk(Φ(Zk)) in (30) is not equal to the L-correlation of Yk=qk(Φ(Zk)) toward Yj=qj(Φ(Zj)) in (34). These L-correlations are equal only when the values of shape parameters c and k associated with qj(Φ(Zj)) and qk(Φ(Zk)) are equal (i.e., when the two distributions are the same). Provided in Algorithm 1 is a source code written in Mathematica [36, 37], which shows an example for computing ICs(rjk) for the L-correlation procedure. The steps for simulating correlated Burr Type III and Type XII distributions with specified values of L-skew (τ3), L-kurtosis (τ4), and with specified L-correlation structure are given in Section 5.Algorithm 1:Mathematica source code for computing intermediate correlations for specified L-correlations. The example is for distribution j=1 toward distribution k=2 (η12). See dashed curves of Distribution 1 and Distribution 2 in Figure 2 (Panel b), specified correlation in Table 8, and intermediate correlation in Table 10.
(* Intermediate Correlation *)
r
12
= 0.784047;
Needs[“MultivariateStatistics̀ ”]
f
12
= PDF[MultinormalDistribution[{0, 0}, {{1, r12}, {r12, 1}}], {Z1, Z2}];
Φ
1
= CDF[NormalDistribution[0, 1], Z1];
Φ
2
= CDF[NormalDistribution[0, 1], Z2];
(* Parameters for dashed curves of Distribution 1 and Distribution 2 in Figure2 (Panel b) *)
c
1 = 1.227122;
k
1 = 3.417227;
λ
1
=
0.441134;
δ
1
=
2.543505;
c
2 = −2.944030;
k
2 = 0.126960;
(* Quantile function from (3) *)
y
1
=
(
(
1
-
Φ
1
)
(
-
1
/
k
1
)
-
1
)
(
1
/
c
1
);
(* Standardizing constantsλ1 and δ1 were obtained, respectively, from  (19) and (33) *)
x
1
=
δ
1
*
(
y
1
-
λ
1
);
(* Compute the specifiedL-correlation *)
η
12
=
2
π
*
NIntegrate[x1*Φ2*f12, {Z1, −8, 8},  {Z2, −8, 8}, Method →{“MultiDimensionalRule”,
“Generators”→ 9}]
0.80


## 5. The Procedure for Monte Carlo Simulation with an Example

The procedure for simulating Burr Type III and Type XII distributions with specifiedL-moments and L-correlations can be summarized in the following six steps.(1)
Specify theL-moments for T transformations of the form in (2) and (3), that is, q1(Φ(z1)),…,qT(Φ(zT)), and obtain the solutions for the shape parameters c and k by simultaneously solving the systems of equations (15) and (21) for 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 values of intermediate (Pearson) correlations (ICs),rjk, by substituting the value of specified L-correlation (ηjk) and the solved values of c and k from Step (1) into the left- and the right-hand sides of (30), respectively, and then numerically integrating (30) to solve for rjk. See Algorithm 1 for an example. Repeat this step separately for all T(T-1)/2 pairwise combinations of ICs.(3)
Assemble the ICs computed in Step (2) into aT×T matrix and then decompose this matrix using Cholesky factorization. Note that this step requires the IC matrix to be positive definite.(4)
Use elements of the matrix resulting from Cholesky factorization of Step (3) to generateT standard normal variables (Z1,…,ZT) correlated at the IC levels as follows:
(35)Z1=a11V1Z2=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 is the element in the ith row and jth column of the matrix resulting from Cholesky factorization of Step (3).(5)
SubstituteZ1,…,ZT from Step (4) into the following Taylor series-based expansion for computing the cdf, Φ(Zj), of standard normal distribution [38]
(36)Φ(Zj)=(12)+ϕ(Zj){Zj+Zj33+Zj5(3·5)+Zj7(3·5·7)+⋯},
where ϕ(Zj) is the pdf of standard normal distribution and the absolute error associated with (36) is less than 8×10-16.(6)
Substitute the uniform(0,1) variables, Φ(Zj), generated in Step (5) into the T equations of the form qj(Φ(zj)) in (2) and (3) to generate the Burr Type III and Type XII distributions with specified values of L-skew (τ3), L-kurtosis (τ4), and with specified L-correlation structure.For the purpose of evaluating the proposed methodology and demonstrating the steps above, an example is subsequently provided to compare theL-correlation-based procedure with the conventional product moment-based (Pearson) correlation procedure. Specifically, the distributions in Figure 2 (dashed curves) are used as a basis for a comparison using the specified correlation matrix in Table 8 where both strong and moderate correlations are considered in a single matrix. Let the four distributions in Figure 2 be Y1=q1(Φ(Z1)), Y2=q2(Φ(Z2)), Y3=q3(Φ(Z3)), and Y4=q4(Φ(Z4)), where Y2 and Y4 are the quantile functions from (2), and Y1 and Y3 are the quantile functions from (3). The specified values of conventional and L-moments associated with these four distributions are given in Tables 4 and 5, respectively. Presented in Tables 9 and 10 are the intermediate correlations (ICs) obtained for the conventional product moment-based (Pearson) correlation and L-moment-based L-correlation procedures, respectively, for the distributions in Figure 2. Provided in Algorithm 2 is a source code written in Mathematica [36, 37], which shows an example for computing ICs (rjk) for the conventional product moment-based (Pearson) correlation procedure. See also Headrick et al. [2, pages 2217–2221] for a detailed methodology for simulating correlated Burr distributions through the method of Pearson correlation.Table 8
Specified correlations for the distributions in Figure2.
1
2
3
4
1
1.0
2
0.80
1.0
3
0.75
0.70
1.0
4
0.60
0.50
0.40
1.0Table 9
Intermediate correlations for the conventional product moment-based (Pearson) correlation procedure.
1
2
3
4
1
1.0
2
0.868005
1.0
3
0.823239
0.804802
1.0
4
0.718908
0.646685
0.415208
1.0Table 10
Intermediate correlations for theL-moment-based L-correlation procedure.
1
2
3
4
1
1.0
2
0.784047
1.0
3
0.731921
0.679639
1.0
4
0.579134
0.478463
0.395474
1.0Algorithm 2:Mathematica source code for computing intermediate correlations for specified conventional product moment-based (Pearson) correlations. The example is for distribution j=1 and distribution k=2 (ρ12). See dashed curves of Distribution 1 and Distribution 2 in Figure 2 (Panel a), specified correlation in Table 8 and intermediate correlation in Table 9.
(* Intermediate Correlation *)
r
12
= 0.868005;
Needs[“MultivariateStatistics̀ ”]
f
12
= PDF[MultinormalDistribution[{0, 0}, {{1, r12}, {r12, 1}}], {Z1, Z2}];
Φ
1
= CDF[NormalDistribution[0, 1], Z1];
Φ
2
= CDF[NormalDistribution[0, 1], Z2];
(* Parameters for dashed curves of Distribution 1 and Distribution 2 in Figure2 (Panel a) *)
c
1 = 1.097740;
k
1 = 4.481719;
c
2 = -5.542062;
k
2 = 0.027237;
(* Quantile function from (3) and (2) *)
y
1
=
(
(
1
-
Φ
1
)
(
-
1
/
k
1
)
-
1
)
(
1
/
c
1
)
y
2
=
(
Φ
2
(
-
1
/
k
2
)
-
1
)
(
1
/
c
2
);
(* Standardizing constantsμ1, μ2 and σ1, σ2 are obtained, respectively, using (14) and (15)
from Headrick et al. [2, page 2211] *)
x
1
=
(
y
1
-
μ
1
)
/
σ
1;
x
2
=
(
y
2
-
μ
2
)
/
σ
2;
(* Compute the specified conventional product moment-based (Pearson) correlation *)
ρ
12
= NIntegrate[x1*x2*f12, {Z1, −8, 8}, {Z2, −8, 8}, Method →{“MultiDimensionalRule”,
“Generators”→ 9}]
0.80Provided in Tables11 and 12 are the results of Cholesky factorization on the IC matrices in Tables 9 and 10, respectively. The elements of matrices in Tables 11 and 12 are used to generate Z1,…,Z4 correlated at the IC levels by making use of the formulae (35) in Step 4 with T=4. The values of  Z1,…,Z4 are then used in (36) to obtain the Taylor series-based approximations of the cdfs Φ(Z1), Φ(Z2),   Φ(Z3), and Φ(Z4), which are treated as uniform (0,1) variables. These uniform variables are used in (2) and (3) to obtain the quantile functions q1(Φ(Z1)),   q2(Φ(Z2)),   q3(Φ(Z3)),   and q4(Φ(Z4)) to generate the four distributions in Figure 2 that are correlated at the specified correlation level of Table 8.Table 11
Cholesky decomposition on the intermediate correlations from Table9.
a
11
=
1.0
a
12
=
0.868005
a
13
=
0.823239
a
14
=
0.718908
a
21
=
0.0
a
22
=
0.496556
a
23
=
0.181704
a
24
=
0.045655
a
31
=
0.0
a
32
=
0.0
a
33
=
0.537830
a
34
=
-
0.343828
a
41
=
0.0
a
42
=
0.0
a
43
=
0.0
a
44
=
0.602387Table 12
Cholesky decomposition on the intermediate correlations from Table10.
a
11
=
1.0
a
12
=
0.784047
a
13
=
0.731921
a
14
=
0.579134
a
21
=
0.0
a
22
=
0.620701
a
23
=
0.170417
a
24
=
0.039301
a
31
=
0.0
a
32
=
0.0
a
33
=
0.659734
a
34
=
-
0.053210
a
41
=
0.0
a
42
=
0.0
a
43
=
0.0
a
44
=
0.812544For the Monte Carlo simulation, Fortran [33] algorithm was written for both procedures to generate 25,000 independent sample estimates for the specified parameters of (a) conventional product moment-based (Pearson) correlation (ρjk) and (b) L-moment-based L-correlation (ηjk) based on samples of sizes n=25 and n=1000. The estimate for ρjk was based on the usual formula for the Pearson correlation statistic. The estimate of  ηjk was computed by substituting (23) and (24) into (25), where the empirical forms of the cdfs were used in (23) and (24). The sample estimates ρjk and ηjk were both transformed using Fisher’s  z′ transformations. Bias-corrected accelerated bootstrapped average estimates (Estimate), 95% bootstrap confidence intervals (95% Bootstrap C.I.), and standard errors (St. error) were obtained for the estimates associated with the parameters (z(ρjk)′ and z(ηjk)′) using 10,000 resamples via the commercial software package Spotfire S+ [34]. The bootstrap results associated with the estimates of  z(ρjk)′ and  z(ηjk)′ were transformed back to their original metrics. Further, if a parameter was outside its associated 95% bootstrap C.I., then the percentage of relative bias (RB%) was computed for the estimate as in (22). The results of this simulation are presented in Tables 13 and 14 and are discussed in Section 6.Table 13
Correlation results for the conventional product moment-based (Pearson) correlations.
Parameter
Estimate
95% Bootstrap C.I.
St. error
RB%
n
=
25
ρ
12
=
0.80
0.8252
(0.8239, 0.8266)
0.00220
3.15
ρ
13
=
0.75
0.7910
(0.7899, 0.7922)
0.00155
5.47
ρ
14
=
0.60
0.6467
(0.6452, 0.6481)
0.00124
7.78
ρ
23
=
0.70
0.7125
(0.7110, 0.7141)
0.00161
1.79
ρ
24
=
0.50
0.5167
(0.5148, 0.5185)
0.00132
3.34
ρ
34
=
0.40
0.4117
(0.4094, 0.4137)
0.00131
2.93
n
=
1000
ρ
12
=
0.80
0.8026
(0.8025, 0.8031)
0.00038
0.33
ρ
13
=
0.75
0.7547
(0.7545, 0.7551)
0.00034
0.63
ρ
14
=
0.60
0.6044
(0.6040, 0.6047)
0.00026
0.73
ρ
23
=
0.70
0.7003
(0.7001, 0.7006)
0.00026
0.04
ρ
24
=
0.50
0.5007
(0.5004, 0.5010)
0.00021
0.14
ρ
34
=
0.40
0.4003
(0.4000, 0.4006)
0.00020
—Table 14
Correlation results for theL-moment-based L-correlation procedure.
Parameter
Estimate
95% Bootstrap C.I.
St. error
RB%
n
=
25
η
12
=
0.80
0.8062
(0.8048, 0.8072)
0.00172
0.78
η
13
=
0.75
0.7571
(0.7556, 0.7585)
0.00172
0.95
η
14
=
0.60
0.6047
(0.6026, 0.6068)
0.00167
0.78
η
23
=
0.70
0.7086
(0.7070, 0.7105)
0.00180
1.23
η
24
=
0.50
0.5057
(0.5031, 0.5082)
0.00176
1.14
η
34
=
0.40
0.4054
(0.4030, 0.4078)
0.00145
1.35
n
=
1000
η
12
=
0.80
0.8001
(0.8000, 0.8003)
0.00026
—
η
13
=
0.75
0.7502
(0.7500, 0.7504)
0.00026
—
η
14
=
0.60
0.6001
(0.5997, 0.6004)
0.00026
—
η
23
=
0.70
0.7002
(0.6999, 0.7004)
0.00028
—
η
24
=
0.50
0.5002
(0.4998, 0.5006)
0.00028
—
η
34
=
0.40
0.4003
(0.3999, 0.4006)
0.00022
—


## 6. Discussion and Conclusion

One of the advantages thatL-moments have over conventional moments can be expressed in the context of estimation. The L-moment-based estimators of L-skew and L-kurtosis can be far less biased than the conventional moment-based estimators of skew and kurtosis when samples are drawn from the distributions with more severe departures from normality [18, 25–29, 31]. Inspection of the simulation results in Tables 6 and 7 clearly indicates that this is the case for the Burr Type III and Type XII distributions. That is, the superiority that estimates of L-moment ratios (τ3 and τ4) have over their corresponding conventional moment-based estimates of skew and kurtosis (γ3 and γ4) is obvious. For example, for samples of size n=25, the estimates of γ3 and γ4 for Distribution 1 were, on average, 41.46% and 5.62% of their associated parameters, whereas the estimates of τ3 and τ4 were 93.88% and 92.03% of their associated parameters. This advantage of L-moment-based estimates can also be expressed by comparing their relative standard errors (RSEs), where RSE={(St.error/Estimate)×100}. Comparing Tables 6 and 7, it is evident that the estimates of τ3 and τ4 are more efficient as their RSEs are considerably smaller than the RSEs associated with the conventional moment-based estimates of γ3 and γ4. For example, in terms of Distribution 1 in Figure 2, inspection of Tables 6 and 7 (for n=1000) indicates that RSE measures of RSE(𝓉3)=0.036% and RSE(𝓉4)=0.064% are considerably smaller than the RSE measures of RSE(g3)=0.270% and RSE(g4)=0.851%. This demonstrates that the estimates of L-skew and L-kurtosis have more precision because they have less variance around their bootstrapped estimates.Another advantage ofL-moments can be highlighted in the context of distribution fitting. Comparison of the four distributions in Figures 2(a) and 2(b) clearly indicates that L-moment-based Burr Type III and Type XII distributions provide a better fit to the theoretical distributions compared with their conventional moment-based counterparts. This advantage is most pronounced in the context of the first two distributions: Distribution 1 and Distribution 2, where L-moment-based Burr Type XII and Type III  in Figure 2(b) provide a better fit to the F distribution (3, 10) and Chi-square (df=1) distributions than their conventional moment-based counterparts in  Figure 2(a). In the context of fitting real-world data, the L-moment-based Burr Type XII in Figure 3(b)  provides a better fit to the ankle circumference data than the conventional moment-based Burr Type XII in Figure 3(a).Presented in Tables13 and 14 are the simulation results of conventional product moment-based (Pearson) correlations and L-moment-based L-correlations, respectively. Overall inspection of these tables indicates that the L-correlation is superior to Pearson correlation in terms of relative bias. For example, for n=25, the percentage of relative bias for the two distributions, Distribution 1 and Distribution 4, in Figure 2 was  7.78% for the Pearson correlation compared with only 0.78% for the L-correlation. It is also noted that the variability associated with bootstrapped estimates of L-correlation appears to be more stable than that of the bootstrapped estimates of Pearson correlation both within and across different conditions.In summary, the newL-moment-based procedure is an attractive alternative to the more traditional conventional moment-based procedure in the context of Burr Type III and Type XII distributions. In particular, the L-moment-based procedure has distinct advantages when distributions with large departures from normality are used. Finally, we note that Mathematica [36, 37] source codes are available from the authors for implementing both the conventional moment- and L-moment-based procedures.



---

*Source: 191604-2013-05-22.xml*