MLE
Calculates maximum likelihood estimates for the parameters of one of several univariate probability distributions.
Required Arguments
X— Array containing the data. (Input)
IPDF — Specifies the probability density function. (Input)
|
Distribution |
IPDF |
size(PARAM) |
i |
PARAM(i) |
|
Discrete uniform |
0 |
1 |
1 |
scale - upper limit |
|
Bernoulli |
1 |
1 |
1 |
probability of success (mean) |
|
Binomial |
2 |
1 |
1 |
probability of success |
|
Negative binomial |
3 |
1 |
1 |
probability of success |
|
Poisson |
4 |
1 |
1 |
location (mean) - θ |
|
Geometric |
5 |
1 |
1 |
probability of success |
|
Continuous uniform |
6 |
2 |
1 |
scale - lower boundary |
|
Beta |
7 |
2 |
1 |
shape - p |
|
Exponential |
8 |
1 |
1 |
scale - b |
|
Gamma |
9 |
2 |
1 |
shape - k |
|
Weibull |
10 |
2 |
1 |
scale - λ |
|
Rayleigh |
11 |
1 |
1 |
scale - α |
|
Extreme value |
12 |
2 |
1 |
location - μ |
|
Generalized extreme value |
13 |
3 |
1 |
location - μ |
|
Pareto |
14 |
2 |
1 |
scale (lower boundary) xm |
|
Generalized Pareto |
15 |
2 |
1 |
scale - σ |
|
Normal |
16 |
2 |
1 |
location (mean) - μ |
|
Log-normal |
17 |
2 |
1 |
location (mean of log(x)) - μ |
|
Logistic |
18 |
2 |
1 |
location (mean) - μ |
|
Log-logistic |
19 |
2 |
1 |
scale (exp(mean)) - eμ |
|
Inverse Gaussian |
20 |
2 |
1 |
location (mean) - μ |
PARAM — Array of length p containing the parameter values, where p = size(PARAM) denotes the number of parameters (see the IPDF table above). On input, the values of PARAM are used as starting values, when USEMM = .FALSE.. On output, final parameter estimates replace the starting values. (Input/Output)
Optional Arguments
NOBS — Number of observations to use in the analysis. (Input)
Default: NOBS = size(X).
PARAMLB — Array of length p containing the lower bounds of the parameters. (Input)
Default: The default lower bound depends on the range of the parameter. That is, if PARAM(i) is positive, PARAMLB(i) = 0.01. If PARAM(i) is non‑negative (≥0), then PARAMLB(i) = 0.0. If PARAM(i) can be any real value, then PARAMLB(i) = ‑10000.00. Exceptions are PARAMLB(i) = 0.25 for the scale parameter of the extreme value distribution, PARAMLB(i) = ‑5.0 for the shape parameter of the generalized Pareto distribution, and PARAMLB(i) = ‑10.0 for the shape parameter of the generalized extreme value distribution.
PARAMUB — Array of length p containing the upper bounds of the parameters. (Input)
Default: PARAMUB(i) = 10000.00. Exceptions are PARAMUB(i) = 5.0 for the shape parameter of the generalized Pareto distribution and PARAMUB(i) = 10.0 for the shape parameter of the generalized extreme value distribution
USEMM — Logical. If .true., starting values are set to the method of moments estimates. (Input)
If USEMM = .FALSE., PARAM values are used.
Default: USEMM = .TRUE..
XSCALE — Array of length p containing the scaling factors for the parameters. XSCALE is used in the routine BCONF mainly in scaling the gradient and the distance between two points. See BCONF in the Math Libray, Chapter 8, “Optimization” for details.
Default: XSCALE =1.0.
MAXIT — Maximum number of iterations. (Input)
Default: MAXIT = 100.
MAXFUN — Maximum number of function evaluations. (Input)
Default: MAXFUN = 400.
MAXGRAD — Maximum number of gradient evaluations. (Input)
Default: MAXGRAD = 400.
IPRINT — Printing option (Input)
|
IPRINT |
Action |
|
0 |
No printing |
|
1 |
Print final results only |
|
2 |
Print intermediate and final results |
Default: IPRINT = 0.
MLOGLIKE — Minus log‑likelihood evaluated at the parameter estimates. (Output)
SE — Array of length p containing the standard errors of the parameter estimates. (Output)
HESS — Array of size p by p containing the Hessian matrix. (Output)
FORTRAN 90 Interface
Generic: CALL MLE (X,IPDF,PARAM [, …])
Specific: The specific interface names are S_MLE and D_MLE.
Description
Routine MLE calculates maximum likelihood estimates for the parameters of a univariate probability distribution, where the distribution is one specified by IPDF and where the input data X is (assumed to be) a random sample from that distribution.
Let {xi, i = 1, …, N} represent a random sample from a probability distribution with density function f(x∣θ), which depends on a vector θ ∈ ℜp containing the values of the parameters of the distribution. The values in θ are fixed but unknown and the problem is to find an estimate for θ given the sample data.
The likelihood function is defined to be the product
The estimator
That is, the estimator that maximizes L also maximizes log L and is the maximum likelihood estimate, or MLE for θ.
The likelihood problem is in general a constrained non-linear optimization problem, where the constraints are determined by the permissible range of θ. In a few situations, the problem has a closed form solution. Otherwise, MLE uses the numerical method as documented in routine BCONF (see Chapter 8, "Optimization" in the Math Library for details) to solve the likelihood problem. If USEMM is .TRUE. (the default), method of moments estimates serve as starting values of the parameters. In some cases, method of moments estimators may not exist, such as when certain moments of the true distribution do not exist; thus it is possible that the starting values are not truly method of moments estimates. If USEMM is set to .FALSE., input values of PARAM are used as starting values.
Upper and lower bounds, when needed for the optimization, have default values for each selection of IPDF (defaults will vary depending on the allowable range of the parameters). It is possible that the optimization will fail. In such cases, the user may try adjusting upper and lower bounds using the optional parameters PARAMLB, PARAMUB, or adjusting up or down the scaling factors in XSCALE, which can sometimes help the optimization converge.
Standard errors and covariances are supplied, in most cases, using the asymptotic properties of ML estimators. Under some general regularity conditions, ML estimates are consistent and asymptotically normally distributed with variance-covariance equal to the inverse Fisher’s Information matrix evaluated at the true value of the parameter, θ0:
MLE approximates the asymptotic variance using the negative inverse Hessian evaluated at the ML estimate:
The Hessian is approximated numerically for all but a few cases where it can be determined in closed form.
In cases when the asymptotic result does not hold, standard errors may be available from the known sampling distribution. For example, the ML estimate of the Pareto distribution location parameter is the minimum of the sample. The variance is estimated using the known sampling distribution of the minium, or first order‑statistic for the Pareto distribution.
For further details regarding the properties of the estimators and the theory of the maximum likelihood method, see Kendall and Stuart (1979). The different probability distributions have wide coverage in the statistical literature. See Johnson & Kotz (1970a, 1970b, or later editions).
Parameter estimation (including maximum likelihoood) for the generalized Pareto distribution is studied in Hosking and Wallis (1987) and Giles and Feng (2009), and estimation for the generalized extreme value distribution is treated in Hosking, Wallis, and Wood (1985).
Comments
1. The location parameter is not estimated for the generalized Pareto distribution (IPDF=15). Instead, the minimum of the sample is subtracted from each observation before the estimation procedure.
2. Only the probability of success parameter is estimated for the binomial and negative binomial distributions, (IPDF=2,3). The number of trials and the number of failures, respectively, must be provided in PARAM(1) on input.
3. MLE issues an error if missing or NaN values are encountered in the input data. Missing or NaN values should be removed before calling MLE.
4. Informational errors
|
Type |
Code |
Description |
|
3 |
1 |
The Hessian is not calculated for the negative binomial distribution. |
|
3 |
2 |
Hessian is not used to calculate the standard errors of the estimates for the continuous uniform distribution. |
|
3 |
3 |
The Hessian is not used to calculate the standard errors of the estimates for the Pareto distribution. |
|
3 |
4 |
For the Pareto distribution, the Hessian cannot be calculated because the parameter estimate is 0. |
Examples
Example 1
The data are N = 100 observations generated from the logistic distribution with location parameter μ = 0.85 and scale parameter σ = 0.5.
use mle_int
implicit none
integer, parameter :: ipdf=18, npar=2
real(kind(1e0)) :: param(npar), stderr(npar), hess(npar,npar)
real(kind(1e0)) :: fval
! Logistic distribution mu = 0.85, sigma=0.5
real(kind(1e0)) :: log1(100)
data log1 /&
2.020394, 2.562315, -0.5453395, 1.258546, 0.7704533, &
0.3662717, 0.6885536, 2.619634, -0.49581, 2.972249, &
0.5356222, 0.4262079, 1.023666, 0.8286033, 1.319018, &
2.123659, 0.3904647, -0.1196832, 1.629261, 1.069602, &
0.9438083, 1.314796, 1.404453, -0.5496156, 0.8326595, &
1.570288, 1.326737, 0.9619384, -0.1795268, 1.330161, &
-0.2916453, 0.7430826, 1.640854, 1.582755, 1.559261, &
0.6177695, 1.739638, 1.308973, 0.568709, 0.2587071, &
0.745583, 1.003815, 1.475413, 1.444586, 0.4515438, &
1.264374, 1.788313, 1.062330, 2.126034, 0.3626510, &
1.365612, 0.5044735, 2.51385, 0.7910572, 0.5932584, &
1.140248, 2.104453, 1.345562, -0.9120445, 0.0006519341, &
1.049729, -0.8246097, 0.8053433, 1.493787, -0.5199705, &
2.285175, 0.9005916, 2.108943, 1.40268, 1.813626, &
1.007817, 1.925250, 1.037391, 0.6767235, -0.3574937, &
0.696697, 1.104745, -0.7691124, 1.554932, 2.090315, &
0.60919, 0.4949385, -2.449544, 0.668952, 0.9480486, &
0.9908558, -1.495384, 2.179275, 0.1858808, -0.3715074, &
0.1447150, 0.857202, 1.805844, 0.405371, 1.425935, &
0.3187476, 1.536181, -0.6352768, 0.5692068, 1.706736/
param = 1.0
stderr = 0.0
hess = 0.0
call mle(log1,ipdf,param,iprint=2,usemm=.true., &
mloglike=fval,se=stderr,hess=hess)
end
Output
Univariate Statistics from UVSTA
Variable Mean Variance Std. Dev. Skewness Kurtosis
1 0.9068 0.8600 0.9274 -0.6251 0.9725
Variable Minimum Maximum Range Coef. Var. Count
1 -2.4495 2.9722 5.4218 1.0227 100.0000
Variable Lower CLM Upper CLM Lower CLV Upper CLV
1 0.7228 1.0908 0.6629 1.1606
Maximum likelihood estimation for the logistic distribution
Starting estimates: 0.90677 0.51128
Initial -log-likelihood: 132.75304
-Log-likelihood 132.61487
MLE for parameter 1 0.95341
MLE for parameter 2 0.50944
Std error for parameter 1 0.08845
Std error for parameter 2 0.04364
Hessian
1 2
1 -127.9 -5.7
2 -5.7 -525.4
Example 2
The data are N = 100 observations generated from the generalized extreme value distribution with location parameter μ = 0, scale parameter σ = 1.0, and shape parameter ξ = ‑0.25.
use mle_int
implicit none
integer, parameter :: ipdf=13, npar=3
real(kind(1e0)) :: param(npar), stderr(npar), &
hess(npar,npar)
real(kind(1e0)) :: fval
! Generalized Extreme Value
! oc = 0, scale =1, shape = =-0.25
real(kind(1e0)) :: gev(100)
data gev/ &
0.7688048, 0.1944504, -0.2992029, -0.3853738, -1.185593, &
0.3056149, -0.4407711, 0.5001115, 0.3635027, -1.058632, &
-0.2927695, -0.3205969, 0.03367599, 0.8850839, 1.860485, &
0.4841038, 0.5421101, 1.883694, 1.707392, 0.2166106, &
1.537204, 1.340291, 0.4589722, 1.616080, -0.8389288, &
0.7057426, 1.532988, 1.161350, 0.9475416, 0.4995294, &
-0.2392898, 0.8167126, 0.992479, -0.8357962, -0.3194499, &
1.233603, 2.321555, -0.3715629, -0.1735171, 0.4624801, &
-0.6249577, 0.7040129, -0.3598889, 0.7121399, -0.5178735, &
-1.069429, 0.7169358, 0.4148059, 1.606248, -0.4640152, &
1.463425, 0.9544342, -1.383239, 0.1393160, 0.622689, &
0.365793, 0.7592438, 0.810005, 0.3483791, 2.375727, &
-0.08124195, -0.4726068, 0.1496043, 0.4961212, 1.532723, &
-0.1106993, 1.028553, 0.856018, -0.6634978, 0.3573150, &
0.06391576, 0.3760349, -0.5998756, 0.4158309, -0.2832369, &
-1.023551, 1.116887, 1.237714, 1.900794, 0.6010037, &
1.599663, -0.3341879, 0.5278575, 0.5497694, 0.6392933, &
0.592865, 1.646261, -1.042950, -1.113611, 1.229645, &
1.655998, 0.6913992, 0.4548073, 0.4982649, -1.073640, &
-0.4765107, -0.8692533, -0.8316462, -0.03609102, 0.655814/
! initialize
param=1.0
stderr=0.0
hess = 0.0
call mle(gev, ipdf, param, iprint=2, usemm=.true., &
mloglike=fval, se=stderr, hess=hess)
end
Output
Univariate Statistics from UVSTA
Variable Mean Variance Std. Dev. Skewness Kurtosis
1 0.3805 0.7484 0.8651 0.05492 -0.6240
Variable Minimum Maximum Range Coef. Var. Count
1 -1.3832 2.3757 3.7590 2.2738 100.0000
Variable Lower CLM Upper CLM Lower CLV Upper CLV
1 0.2088 0.5521 0.5769 1.0100
Maximum likelihood estimation for the generalized extreme value distribution
Starting estimates: -0.00888 0.67451 0.00000
Initial -log-likelihood: 135.43820
-Log-likelihood 126.09403
MLE for parameter 1 0.07500
MLE for parameter 2 0.85115
MLE for parameter 3 -0.27960
Std error for parameter 1 0.09467
Std error for parameter 2 0.07007
Std error for parameter 3 0.06695
Hessian
1 2 3
1 -141.1 -51.3 -112.8
2 -51.3 -337.0 -241.0
3 -112.8 -241.0 -438.8
