GB2 (version 2.1.1)

MLprofGB2: 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 z. By default w is a vector of 1.

method

numeric; the method to be used by optim. By default, codemethod = 1 and the used method is BFGS. If method = 2, method L-BFGS-B is used.

hess

logical; By default, hess = FALSE, the hessian matrix is not calculated.

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.

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).

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

# NOT RUN {
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)
# }
# NOT RUN {
# Plots
# }
# NOT RUN {
plotsML.gb2(inc,ap,bp,pp,qp,w)
# }