MLE of the Dirichlet distribution via Newton-Rapshon {Compositional}

MLE of the Dirichlet distribution via Newton-Rapshon

Description

MLE of the Dirichlet distribution via Newton-Rapshon.

Usage

diri.nr(x, type = 1, tol = 1e-07)

Arguments

x

A matrix containing compositional data. Zeros are not allowed.

type

Type can either be 1, so that the Newton-Rapshon is used for the maximisation of the log-likelihood, as Minka (2012) suggested or it can be 1. In the latter case the Newton-Raphson algorithm is implemented involving matrix inversions. In addition an even faster implementation has been implemented (in C++) in the package Rfast and is used here.

tol

The tolerance level indicating no further increase in the log-likelihood.

Details

Maximum likelihood estimation of the parameters of a Dirichlet distribution is performed via Newton-Raphson. Initial values suggested by Minka (2003) are used. The estimation is super faster than "diri.est" and the difference becomes really apparent when the sample size and or the dimensions increase. In fact this will work with millions of observations. So in general, I trust this one more than "diri.est".

The only problem I have seen with this method is that if the data are concentrated around a point, say the center of the simplex, it will be hard for this and the previous methods to give estimates of the parameters. In this extremely difficult scenario I would suggest the use of the previous function with the precision parametrization "diri.est(x, type = "prec")". It will be extremely fast and accurate.

Value

A list including:

iter

The number of iterations required. If the argument "type" is set to 2 this is not returned.

loglik

The value of the log-likelihood.

param

The estimated parameters.

Author(s)

Michail Tsagris.

R implementation and documentation: Michail Tsagris [e-mail].

References

Thomas P. Minka (2003). Estimating a Dirichlet distribution. http://research.microsoft.com/en-us/um/people/minka/papers/dirichlet/minka-dirichlet.pdf

See Also

diri.est, diri.contour rdiri, ddiri, dda

Examples

x <- rdiri( 100, c(5, 7, 5, 8, 10, 6, 4) )
diri.nr(x)
diri.nr(x, type = 2)
diri.est(x)