MLprofGB2 {GB2} | R Documentation |
Maximum Likelihood Estimation of the GB2 Based on the Profile Log-likelihood
Description
profml.gb2
performs maximum likelihood estimation based on the profile log-likelihood through the general-purpose optimisation function optim
from package stats
.
Usage
profml.gb2(z, w=rep(1, length(z)), method=1, hess = FALSE)
Arguments
z |
numeric; vector of data values. |
w |
numeric; vector of weights. Must have the same length as |
method |
numeric; the method to be used by |
hess |
logical; By default, |
Details
Two methods are considered: BFGS and L-BFGS-B (see optim
documentation for more details). Initial values of the parameters to be optimized over (a
and b
)
are given from the Fisk distribution. The function to be maximized by optim
is the negative of the profile log-likelihood (proflogl.gb2
)
and the gradient is equal to the negative of the scores (profscores.gb2
).
Value
A list with 1 argument: opt1
for the output of the BFGS fit or opt2
for the output of the L-BFGS fit. Further values are given by the values of optim
.
Author(s)
Monique Graf
References
Graf, M., Nedyalkova, D., Muennich, R., Seger, J. and Zins, S. (2011) AMELI Deliverable 2.1: Parametric Estimation of Income Distributions and Indicators of Poverty and Social Exclusion. Technical report, AMELI-Project.
See Also
optim
for the general-purpose optimization, link{ml.gb2}
for the fit using the full log-likelihood and fisk
for the Fisk distribution.
Examples
library(laeken)
data(eusilc)
# Income
inc <- as.vector(eusilc$eqIncome)
# Weights
w <- eusilc$rb050
# Data set
d <- data.frame(inc,w)
d <- d[!is.na(d$inc),]
# Truncate at 0
inc <- d$inc[d$inc > 0]
w <- d$w[d$inc > 0]
# Fit using the profile log-likelihood
fitp <- profml.gb2(inc, w)$opt1
# Fitted GB2 parameters
ap <- fitp$par[1]
bp <- fitp$par[2]
pp <- prof.gb2(inc, ap, bp, w)[3]
qp <- prof.gb2(inc, ap, bp, w)[4]
# Profile log-likelihood
proflik <- fitp$value
# If we want to compare the indicators
## Not run:
# GB2 indicators
indicp <- round(main.gb2(0.6,ap,bp,pp,qp), digits=3)
# Empirical indicators
indice <- round(main.emp(inc,w), digits=3)
## End(Not run)
# Plots
## Not run: plotsML.gb2(inc,ap,bp,pp,qp,w)