lmomco (version 1.3.3)

pargpaRC: Estimate the Parameters of the Generalized Pareto Distribution with Right-Tail Censoring

Description

This function estimates the parameters ($\xi$, $\alpha$, and $\kappa$) of the Generalized Pareto distribution given the B-type L-moments (through the B-type Probability-Weighted Moments) of the data under right censoring conditions (see pwmRC). The B-type L-moments in terms of the parameters are

$$\lambda^B_1 = \xi + \alpha m_1 \mbox{,}$$ $$\lambda^B_2 = \alpha (m_1 - m_2) \mbox{,}$$ $$\lambda^B_3 = \alpha (m_1 - 3m_2 + 2m_3)\mbox{,}$$ $$\lambda^B_4 = \alpha (m_1 - 6m_2 + 10m_3 - 5m_4)\mbox{, and}$$ $$\lambda^B_5 = \alpha (m_1 - 10m_2 + 30m_3 - 35m_4 + 14m_5)\mbox{,}$$

where $m_r = \lbrace 1-(1-\zeta)^{r+\kappa} \rbrace/(r+\kappa)$ and $\zeta$ is the right-tail censor fraction or the probability $\mathrm{Pr}\lbrace \rbrace$ that $x$ is less than the quantile at $\zeta$ nonexceedance probability: ($\mathrm{Pr}\lbrace x < X(\zeta) \rbrace$). Finally, the RC in the function name is to denote Right-tail Censoring.

Usage

pargpaRC(lmom,zeta=1,xi=NULL,lower=-1,upper=20,checklmom=TRUE)

Arguments

lmom

A B-type L-moment object created by a function such as pwm2lmom from B-type Probability-Weighted Moments from pwmRC.

zeta

The censoring fraction. The number of samples observed (noncensored) divided by the total number of samples.

xi

The lower limit of the distribution. If $\xi$ is known then alternative algorithms are used.

lower

The lower value for $\kappa$ for a call to the optimize function. For the L-moments of the distribution to be valid $\kappa > -1$.

upper

The upper value for $\kappa$ for a call to the optimize function. Hopefully, a large enough default is chosen for real-world data sets.

checklmom

Should the lmom be checked for validity using the are.lmom.valid function. Normally this should be left as the default and it is very unlikely that the L-moments will not be viable (particularly in the $\tau_4$ and $\tau_3$ inequ

Value

• An R list is returned.
• type
  The type of distribution: gpa.
• para
  The parameters of the distribution.
• zeta
  The right-tail censoring fraction.
• source
  The source of the parameters: pargpaRC.
• optim
  The list returned by the optimize function.

Details

The optimize function is used to numerically solve for the shape parameter $\kappa$. No test or evaluation is made on the quality of the minimization. Users should consult the contents of the optim portion of the returned list. Finally, this function should return the same parameters if $\zeta=1$ as the pargpa function.

References

Hosking, J.R.M., 1990, L-moments---Analysis and estimation of distributions using linear combinations of order statistics: Journal of the Royal Statistical Society, Series B, vol. 52, p. 105--124.

Hosking, J.R.M., 1995, The use of L-moments in the analysis of censored data, in Recent Advances in Life-Testing and Reliability, edited by N. Balakrishnan, chapter 29, CRC Press, Boca Raton, Fla., pp. 546--560.

See Also

lmomgpa, lmomgpaRC, pargpa, cdfgpa, quagpa

Examples

n         <- 60 # samplesize
para      <- vec2par(c(1500,160,.3),type="gpa") # build a GPA parameter set
fakedata  <- quagpa(runif(n),para) # generate n simulated values
threshold <- 1700 # a threshold to apply the simulated censoring
fakedata  <- sapply(fakedata,function(x) { if(x > threshold) 
                                           return(threshold) else return(x)
                                         })
lmr       <- lmoms(fakedata) # Ordinary L-moments without considering
                             # that the data is censored
estpara   <- pargpa(lmr) # Estimated parameters of parent 

pwm2     <- pwmRC(fakedata,threshold=threshold) # compute censored PWMs
typeBpwm <- pwm2$Bbetas # the B-type PWMs
zeta     <- pwm2$zeta # the censoring fraction

cenpara <- pargpaRC(pwm2lmom(typeBpwm),zeta=zeta) # Estimated parameters
F       <- nonexceeds() # nonexceedance probabilities for plotting purposes

# Visualize some data
plot(F,quagpa(F,para), type='l', lwd=3) # The true distribution
lines(F,quagpa(F,estpara), col=3) # Green estimated in the ordinary fashion
lines(F,quagpa(F,cenpara), col=2) # Red, consider that the data is censored
# now add in what the drawn sample looks like.
PP <- pp(fakedata) # plotting positions of the data
points(PP,sort(fakedata)) # sorting is needed!
# Interpretation. You should see that the red line more closely matches
# the heavy black line. The green line should be deflected to the right
# and pass through the values equal to the threshold, which reflects the
# much smaller L-skew of the ordinary L-moments compared to the type-B
# L-moments.


# Assertion, given some PWMs or L-moments, if zeta=1 then the parameter
# estimates must be identical. The following provides a demonstration.
para1 <- pargpaRC(pwm2lmom(typeBpwm),zeta=1)
para2 <- pargpa(pwm2lmom(typeBpwm))
str(para1)
str(para2)

# Assertion as previous assertion, let us trigger different optimizer
# algorithms with a non-NULL xi parameter and see if the two parameter
# lists are the same.
para1 <- pargpaRC(pwm2lmom(typeBpwm), zeta=zeta)
para2 <- pargpaRC(pwm2lmom(typeBpwm), xi=para1$para[1], zeta=zeta)
str(para1)
str(para2)