distrMod: MLEstimator – R documentation

Become an expert in R — Interactive courses, Cheat Sheets, certificates and more!

MLEstimator

Function to compute maximum likelihood estimates

Description

The function MLEstimator provides a general way to compute maximum likelihood estimates for a given parametric family of probability measures. This is done by calling the function MCEstimator which minimizes the negative log-Likelihood.

Usage

MLEstimator(x, ParamFamily, startPar = NULL, 
            Infos, trafo = NULL, penalty = 1e20,
            validity.check = TRUE, na.rm = TRUE, ...,
            .withEvalAsVar = TRUE, dropZeroDensity = TRUE, nmsffx = "",
            .with.checkEstClassForParamFamily = TRUE)

Arguments

`x`	(empirical) data
`ParamFamily`	object of class `"ParamFamily"`
`startPar`	initial information used by `optimize` resp. `optim`; i.e; if (total) parameter is of length 1, `startPar` is a search interval, else it is an initial parameter value; if `NULL` slot `startPar` of `ParamFamily` is used to produce it; in the multivariate case, `startPar` may also be of class `Estimate`, in which case slot `untransformed.estimate` is used.
`Infos`	character: optional informations about estimator
`trafo`	an object of class `MatrixorFunction` – a transformation for the main parameter
`penalty`	(non-negative) numeric: penalizes non valid parameter-values
`validity.check`	logical: shall return parameter value be checked for validity? Defaults to yes (`TRUE`)
`na.rm`	logical: if `TRUE`, the estimator is evaluated at `complete.cases(x)`.
`...`	further arguments to `criterion` or `optimize` or `optim`, respectively.
`.withEvalAsVar`	logical: shall slot `asVar` be evaluated (if `asvar.fct` is given) or just the call be returned?
`dropZeroDensity`	logical of length 1; shall observations with density zero be dropped? Optimizers like `optim` require finite values, so get problems when negative loglikelihood is evaluated.
`nmsffx`	character: a potential suffix to be appended to the estimator name.
`.with.checkEstClassForParamFamily`	logical: Should a the end of the function `.checkEstClassForParamFamily`; defaults to `TRUE`; can be switched off for computational time or because this is already checked in a calling wrapper function.

Details

The function uses mleCalc for method dispatch; this method by default calls mceCalc using the negative log-likelihood as criterion which should be minimized.

Value

An object of S4-class "MCEstimate" which inherits from class "Estimate".

Author(s)

Matthias Kohl Matthias.Kohl@stamats.de,
Peter Ruckdeschel peter.ruckdeschel@uni-oldenburg.de

Examples

#############################
## 1. Binomial data
#############################
## (empirical) data
# seed for reproducibility:
set.seed(20200306)
x <- rbinom(100, size=25, prob=.25)

## ML-estimate
MLEstimator(x, BinomFamily(size = 25))


#############################
## 2. Poisson data
#############################
## Example: Rutherford-Geiger (1910); cf. Feller~(1968), Section VI.7 (a)
x <- c(rep(0, 57), rep(1, 203), rep(2, 383), rep(3, 525), rep(4, 532), 
       rep(5, 408), rep(6, 273), rep(7, 139), rep(8, 45), rep(9, 27), 
       rep(10, 10), rep(11, 4), rep(12, 0), rep(13, 1), rep(14, 1))

## ML-estimate
MLEstimator(x, PoisFamily())


#############################
## 3. Normal (Gaussian) location and scale
#############################
## (empirical) data
# seed for reproducibility:
set.seed(20200306)
x <- rnorm(100)

## ML-estimate
MLEstimator(x, NormLocationScaleFamily())
## compare:
c(mean(x),sd(x))


#############################
## 4. Gamma model
#############################
## (empirical) data
# seed for reproducibility:
set.seed(20200306)
x <- rgamma(50, scale = 0.5, shape = 3)

## parametric family of probability measures
G <- GammaFamily(scale = 1, shape = 2)

## Maximum likelihood estimator
(res <- MLEstimator(x = x, ParamFamily = G))

## Asymptotic (CLT-based) confidence interval
confint(res)

## some profiling
par(mfrow=c(1,2))
plot(profile(res))
par(mfrow=c(1,1))

## implementation of ML-estimator of package MASS
require(MASS)
(res1 <- fitdistr(x, "gamma"))

## comparison
## shape
estimate(res)[2]
## rate
1/estimate(res)[1]

## minor differences due to the fact that by default, fitdistr uses
## BFGS, while we use Nelder-Mead instead

## log-likelihood
res1$loglik
## negative log-likelihood
criterion(res)


## explicitely transforming to
## MASS parametrization:
mtrafo <- function(x){
     nms0 <- names(c(main(param(G)),nuisance(param(G))))
     nms <- c("shape","rate")
     fval0 <- c(x[2], 1/x[1])
     names(fval0) <- nms
     mat0 <- matrix( c(0, -1/x[1]^2, 1, 0), nrow = 2, ncol = 2,
                     dimnames = list(nms,nms0))                          
     list(fval = fval0, mat = mat0)}

G2 <- G
trafo(G2) <- mtrafo
res2 <- MLEstimator(x = x, ParamFamily = G2)

old <- getdistrModOption("show.details")
distrModoptions("show.details" = "minimal")
res1
res2

## some profiling
par(mfrow=c(1,2))
plot(profile(res2))
par(mfrow=c(1,1))

#############################
## 5. Cauchy Location Scale model
#############################
(C <- CauchyLocationScaleFamily())
loc.true <- 1
scl.true <- 2

## (empirical) data
# seed for reproducibility:
set.seed(20200306)
x <- rcauchy(50, location = loc.true, scale = scl.true)

## Maximum likelihood estimator
(res <- MLEstimator(x = x, ParamFamily = C))
## Asymptotic (CLT-based) confidence interval
confint(res)

distrMod

Object Oriented Implementation of Probability Models

v2.8.4

LGPL-3

Authors

Matthias Kohl [aut, cph], Peter Ruckdeschel [cre, cph], R Core Team [ctb, cph] (for source file 'format.perc')

Initial release

2020-03-06