New distributions are still being suggested for better fitting of a distribution to data, as it is one of the most fundamental problems in terms of the parametric approach. One of such is weighted Lindley (WL) distribution (Ghitany et al. 2011). Even though WL distribution has become increasingly popular as a possible alternative to traditional distributions such as gamma and log normal distributions, fitting it to data has rarely been addressed in existing R packages. This is the reason we present the WLinfer package that implements overall statistical inference for WL distribution. In particular, WLinfer enables one to conduct the goodness of fit test, point estimation, bias correction, interval estimation, and the likelihood ratio test simply with the WL
function which is at the core of this package. To assist users who are unfamiliar with WL distribution, we present a brief review followed by an illustrative example with R codes.
Weighted Lindley (WL) distribution has recently received considerable
attention since it provides a more flexible fit to data from various
fields than traditional widely-used distributions such as exponential,
log normal, and gamma distributions (Ghitany et al. 2011; Mazucheli et al. 2013).
The probability density function (pdf) of WL distribution is given by
Despite this difficulty, there is no existing R package that implements this comprehensive process for WL distribution. To the best of our knowledge, mle.tools (Mazucheli et al. 2017) is the only R package that enables one to obtain maximum likelihood (ML) estimates for WL distribution with asymptotic variance and bias correction. However, they are just limited to ML estimates. An R package fitdistrplus (Delignette-Muller and Dutang 2015) which fits certain univariate distributions to data sets is not applicable to WL distribution, while LindleyR package (Mazucheli et al. 2016) is exclusively for the use of pdf, cumulative distribution function (cdf), quantile, and random number generation, not statistical inference itself.
Based on this motivation, we present an R package
WLinfer by providing
various estimation methods in addition to maximum likelihood estimator
(MLE), such as method of moment estimator (MME), modified method of
moment estimator (MME
The remainder of this paper is organized as follows. First, we briefly
review the theoretical results for each inferential step ranging from
the goodness of fit test to the likelihood ratio test, while introducing
relevant arguments of the WL
function. We then provide an example to
illustrate the whole output of WL
function and how to use WL
function which implements the whole statistical inference at once.
Finally, we conclude with summarizing remarks.
The WLinfer package is available from the Comprehensive R Archive Network (CRAN) at https://CRAN.R-project.org/package=WLinfer. R code for the examples demonstrated herein has been provided as supplementary material. The supplementary code has been tested with WLinfer version 1.0.0, and results presented herein have been produced with this version.
The first step in fitting a distribution to data is to check whether the
data can be assumed as generated from the distribution. For this
purpose, many goodness of fit tests have been developed. Among them,
WLinfer includes three popular tests: the Kolmogorov-Smirnov test, the
Anderson-Darling test, and the Cramér-von Mises test. Uncertainty that
could occur from using estimates instead of true values are not
considered here. The argument dist
test
of WL
function is
specified by one of either "ks"
, "cvm"
, "ad"
, or "all"
. The
default is dist
test="ks"
while dist
test=all"
returns the
results of all three tests.
WLinfer package considers four estimators: MLE, MLE
Each estimation method can be implemented by choosing est
method
of WL
function, among "MLE"
, "MLEc"
, "MME"
, and "MMEm"
. The
default refers to est
method = "MLEc"
. If point estimation is
solely of interest, one can separately obtain point estimates with the
following codes:
data(fail_fiber) # fail_fiber is included in WLinfer package.
MLEc_WL(fail_fiber) # Closed form MLE-like estimator
MME_WL(fail_fiber) # Method of moment estimator
MMEm_WL(fail_fiber) # Modified method of moment estimator
MLE_WL(fail_fiber,init =1) # Initial value for phi is required.
We use fail
fiber
data set (Bader and Priest 1982) which is included in the
WLinfer package. This data set consists of 65 observations of the
strength measurement of carbon fiber.
Unless sample size is sufficient for achieving consistency, bias
correction is highly recommended since all the aforementioned estimators
are upwardly biased in the case of small samples
(Mazucheli et al. 2013; Wang and Wang 2017; Kim and Jang 2020). To this end, the WLinfer
package provides three bias correction methods: Cox and Snell’s method
(Cox and Snell 1968; Cordeiro and Klein 1994) for MLE and MLE
Cox and Snell’s method
where
Firth’s method
The estimator corrected by Firth’s method is obtained by solving
modified likelihood equations given by
Bootstrap method
For bias correction, the ‘bias
cor
’ argument should be specified
since the default is bias
cor=NULL
. Note, unlike Cox and Snell’s
method that works with both MLE and MLEboot_iter
, while the default is 1000.
WL(fail_fiber, est_method = "MLE",bias_cor = "firth") # Firth's method
WL(fail_fiber, est_method = "MLEc",bias_cor = "coxsnell") # Cox and Snell's method
WL(fail_fiber, est_method = "MLEc",bias_cor = "boots") # The bootstrap method
MLE, MLE
One problem is that the lower bound of estimated intervals can be
negative. A confidence interval with negative lower bound may not be
useful since the parameters for WL distribution must be positive.
Therefore, in case of negative lower bounds, confidence intervals for
When asymptotic distribution is not available or the sample size is
insufficient, the bootstrap confidence interval can be a good
alternative. The basic bootstrap confidence interval is simply given by
The default arguments of WL
function for interval estimation are
CI
method="asymp"
, CI
scale="normal"
,
CI
side="two"
, and CI
alpha=0.05
. If any bias correction
method is used, CI
method="boot"
is automatically chosen.
CI
scale="exp"
is recommended in the case of a negative lower
confidence limit and one-sided intervals can be obtained by choosing
CI
side="one"
.
WL(fail_fiber,est_method="MLEc")$CI_list #asymptotic CI for MLEc
WL(fail_fiber,est_method="MLE")$CI_list #asymptotic CI for MLE
# Bootstrap CI for MLEc corrected by Cox and Snell's method
WL(fail_fiber,est_method = "MLEc", bias_cor = "coxsnell")$CI_list
Let WL
function unless
wilks
test="FALSE"
is selected. The significance level and types
of null hypotheses are set by wilks
alpha
and wilks
side
.
The default setting is wilks
alpha=0.05
and
wilks
side="two"
.
The LRT can be separately conducted with the wilks.test
function apart
from other inferences. This function is especially useful when testing
several possible values or regions. By specifying arguments estimator
and side
, both simple and composite null hypotheses can be tested.
H
wilks.test(fail_fiber,estimator = MME_WL(fail_fiber),side = "two")
H
wilks.test(fail_fiber,estimator=c(1,1),side="less")
WL
functionIn this section, we provide an example for illustrating how to conduct
the aforementioned steps of statistical inference from the beginning
using the WL
function. The WL
function returns an S3 object of class
‘WL’ which is assigned to fiber
in this example. For this object of
class ‘WL’, summary
and plot
functions are provided.
> data("fail_fiber")
> fiber = WL(fail_fiber,dist_test = "ks", est_method ="MLEc",
wilks_alpha=0.05, wilks_side="two")
> summary(fiber)
Data: fail_fiber
Data summary:
Mean: 2.244, Variance: 0.1728
Min 1st Qu Median 3rd Qu Max
1.339 1.931 2.272 2.558 3.174
Kolmogorov-Smirnov test for alpha=0.05
D = 0.0723, p-value: 0.8864
>> This data follows the weighted Lindley distribution with estimated parameters.
Estimation method (bias correction): MLEc(None)
lambda: 12.9318, phi: 28.3329
Variance of lambda & phi:
Var(lambda)= 5.1126, Var(phi)= 25.2938
Two-sided confidence interval for 95%
(Asymptotic method & Normal scaled )
CI for lambda: (8.5001, 17.3635)
CI for phi: (18.4757, 38.1902)
Two-sided Wilks' theorem test for estimated parameters
X = 0.0024, p-value: 0.9988
>>The null hypothesis cannot be rejected.
The simple descriptive statistic is provided by summary(fiber)
,
followed by the result of the goodness of fit test. This data set can be
assumed to originate from WL distribution since the null hypothesis of
the KS test is not rejected. We chose MLEVariance of lambda & phi
in the output means estimated asymptotic
variances obtained by replacing
To obtain some helpful plots, we use plot
function. Four plots are
returned by plot(fiber)
as per Figure 1: the
histogram with estimated density function, the boxplot for detecting
outliers, QQ plot, and contour plot with point estimates. The first
three plots provide insight on how much WL distribution is suitable for
the given data, along with the goodness of fit tests. The last contour
plot provides the surface of log likelihood near point estimates. Note
that the density function in the histogram and the quantiles for the QQ
plot are calculated based on estimated parameters.
We presented an R package WLinfer that implements a goodness of fit test, several types of point estimation, bias correction, interval estimation, and the likelihood ratio test, and provide some useful plots. We supply a set of simple codes and an illustrative example of how to apply this package in practice. This package could practically assist practitioners, removing the need to make codes from scratch.
The corresponding author’s research was supported by the Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education(2018R1D1A1B07045603) and a National Research Foundation of Korea (NRF) grant funded by the Korean government (MSIT) (2021R1A4A5032622).
WLinfer, mle.tools, fitdistrplus, LindleyR
ActuarialScience, Distributions, Survival
This article is converted from a Legacy LaTeX article using the texor package. The pdf version is the official version. To report a problem with the html, refer to CONTRIBUTE on the R Journal homepage.
Text and figures are licensed under Creative Commons Attribution CC BY 4.0. The figures that have been reused from other sources don't fall under this license and can be recognized by a note in their caption: "Figure from ...".
For attribution, please cite this work as
Jang, et al., "WLinfer: Statistical Inference for Weighted Lindley Distribution", The R Journal, 2023
BibTeX citation
@article{RJ-2022-042, author = {Jang, Yu-Hyeong and Kim, SungBum and Jung, Hyun-Ju and author), Hyoung-Moon Kim (Corresponding}, title = {WLinfer: Statistical Inference for Weighted Lindley Distribution}, journal = {The R Journal}, year = {2023}, note = {https://rjournal.github.io/}, volume = {14}, issue = {3}, issn = {2073-4859}, pages = {13-19} }