maxp: Maximum likelihood estimation

Description

Find the maximum likelihood estimate for p, also equal probabilities

Usage

maxp(H, startp=NULL, give=FALSE, fcm=NULL, fcv=NULL, SMALL=1e-6, n=1,
   show=FALSE, justlikes=FALSE, ...)
maxplist(Hlist, startp=NULL, give=FALSE, fcm=NULL, fcv=NULL, SMALL=1e-6, ...)
maxp_single(H, startp=NULL, give=FALSE, fcm=NULL, fcv=NULL, SMALL=1e-6,
   maxtry=100, ...)
maxp_single2(H, startp=NULL, give=FALSE, fcm=NULL, fcv=NULL, SMALL=1e-6,
   maxtry=100, ...)
maxp_simplex(H, n=100, show=FALSE, give=FALSE, ...)
maxp_lsl(HLSL, startp = NULL, give = FALSE, fcm = NULL, fcv = NULL, SMALL=1e-6, ...)
equalp(H)

Arguments

H: A hyper2 or hyper3 object
Hlist: A list with elements all hyper2 objects
HLSL: An lsl object
startp: A vector of probabilities specifying the start-point for optimization; if a full unit-sum vector, then the fill-up value will be removed by indep() (except for maxp_lsl())
give: Boolean, with default FALSE meaning to return just the evaluate (including fillup), and TRUE meaning to return the entire formal output of the optimization routine. In function maxp(), a non-logical value [such as 0] gives only the two most important pieces of information, viz the evaluate and a log-likelihood
fcm,fcv: Further problem-specific constraints
n: Number of start points to use
show: Boolean, with TRUE meaning to show successive estimates
justlikes: Boolean, with TRUE meaning to return just a vector of estimated likelihoods
SMALL: Numerical minimum for probabilities
maxtry: Integer specifying maximum number of times to try constrOptim() with slightly differing start points, to avoid a known R bug which reports wmmin is not finite, bugzilla id 17703
...: Further arguments which maxp() passes to constrOptim()

Author

Robin K. S. Hankin

Details

Function maxp() returns the maximum likelihood estimate for p, which has the unit sum constraint implemented.

Function maxplist() does the same but takes a list of hyper2 objects (for example, the output of ggrl()). Note that maxplist() does not have access to the gradient of the objective function, which makes it slow.

If function maxp() is given a suplist object it dispatches to maxplist().

Functions maxp_single() and maxp_single2() are helper functions which perform a single constrained optimization using base::constrOptim() or alabama::constrOptim.nl() respectively. The functions should produce identical (or at least very similar) results. They are used by maxp() and maxp_simplex() which dispatch to either maxp_single() or maxp_single2() depending on the value of option use_alabama. If TRUE, they will use (experimental) maxp_single2(), otherwise (default) maxp_single(). Function maxp_single() is prone to the “wmmin not finite” bug [bugzilla id 17703] but on the other hand is a bit slower. I am not sure which one is better at this time.

Function maxp_simplex() is intended for complicated or flat likelihood functions where finding local maxima might be a problem. It repeatedly calls maxp_single(), starting from a different randomly chosen point in the simplex each time. This function does not take fcm or fcv arguments, it operates over the whole simplex (hence the name). Further arguments, ..., are passed to maxp_single().

The functions do not work for the masterchef_series6 likelihood function. These require a bespoke optimization as shown in the vignette.

Function equalp() returns the value of \(p\) for which all elements are the same.

In functions maxp() etc, arguments fcm and fcv implement linear constraints to be passed to constrOptim(). These constraints are in addition to the usual nonnegativity constraints and unit-sum constraint, and are added to the ui and ci arguments of constrOptim() with rbind() and c() respectively. The operative lines are in maxp_single():


    UI <- rbind(diag(nrow = n - 1), -1, fcm)
    CI <- c(rep(SMALL, n - 1), -1 + SMALL, fcv)

where in UI, the first \(n-1\) rows enforce nonnegativity of \(p_i\), \(1\leq p < n\); row \(n\) enforces nonnegativity of the fillup value \(p_n\); and the remaining (optional) rows enforce additional linear constraints. Argument CI is a vector with corresponding elements.

Examples of their use are given in the “icons” vignette.

Examples

Run this code


maxp(icons)

W <- hyper2(pnames=letters[1:5])
W1 <- ggrl(W, 'a', letters[2:3],'d')  # W1 is a suplist object
if (FALSE) maxp(W1)  # takes a long time to maximize a suplist