NAG Library Function Document
nag_censored_normal (g07bbc)
1
Purpose
nag_censored_normal (g07bbc) computes maximum likelihood estimates and their standard errors for arguments of the Normal distribution from grouped and/or censored data.
2
Specification
#include <nag.h> |
#include <nagg07.h> |
void |
nag_censored_normal (Nag_CEMethod method,
Integer n,
const double x[],
const double xc[],
const Integer ic[],
double *xmu,
double *xsig,
double tol,
Integer maxit,
double *sexmu,
double *sexsig,
double *corr,
double *dev,
Integer nobs[],
Integer *nit,
NagError *fail) |
|
3
Description
A sample of size
is taken from a Normal distribution with mean
and variance
and consists of grouped and/or censored data. Each of the
observations is known by a pair of values
such that:
The data is represented as particular cases of this form:
- exactly specified observations occur when ,
- right-censored observations, known only by a lower bound, occur when ,
- left-censored observations, known only by a upper bound, occur when ,
- and interval-censored observations when .
Let the set
identify the exactly specified observations, sets
and
identify the observations censored on the right and left respectively, and set
identify the observations confined between two finite limits. Also let there be
exactly specified observations, i.e., the number in
. The probability density function for the standard Normal distribution is
and the cumulative distribution function is
The log-likelihood of the sample can be written as:
where
and
.
Let
and
then the first derivatives of the log-likelihood can be written as:
and
The maximum likelihood estimates,
and
, are the solution to the equations:
and
and if the second derivatives
,
and
are denoted by
,
and
respectively, then estimates of the standard errors of
and
are given by:
and an estimate of the correlation of
and
is given by:
To obtain the maximum likelihood estimates the equations
(1) and
(2) can be solved using either the Newton–Raphson method or the Expectation-maximization
algorithm of
Dempster et al. (1977).
Newton–Raphson Method
This consists of using approximate estimates
and
to obtain improved estimates
and
by solving
for the corrections
and
.
EM Algorithm
The expectation step consists of constructing the variable
as follows:
the maximization step consists of substituting
(3),
(4),
(5) and
(6) into
(1) and
(2) giving:
and
where
and where
,
and
are
,
and
evaluated at
and
. Equations
(3) to
(8) are the basis of the
iterative procedure for finding
and
. The procedure consists of alternately estimating
and
using
(7) and
(8) and estimating
using
(3) to
(6).
In choosing between the two methods a general rule is that the Newton–Raphson method converges more quickly but requires good initial estimates whereas the algorithm converges slowly but is robust to the initial values. In the case of the censored Normal distribution, if only a small proportion of the observations are censored then estimates based on the exact observations should give good enough initial estimates for the Newton–Raphson method to be used. If there are a high proportion of censored observations then the algorithm should be used and if high accuracy is required the subsequent use of the Newton–Raphson method to refine the estimates obtained from the algorithm should be considered.
4
References
Dempster A P, Laird N M and Rubin D B (1977) Maximum likelihood from incomplete data via the algorithm (with discussion) J. Roy. Statist. Soc. Ser. B 39 1–38
Swan A V (1969) Algorithm AS 16. Maximum likelihood estimation from grouped and censored normal data Appl. Statist. 18 110–114
Wolynetz M S (1979) Maximum likelihood estimation from confined and censored normal data Appl. Statist. 28 185–195
5
Arguments
- 1:
– Nag_CEMethodInput
-
On entry: indicates whether the Newton–Raphson or
algorithm should be used.
If , the Newton–Raphson algorithm is used.
If , the algorithm is used.
Constraint:
or .
- 2:
– IntegerInput
-
On entry: , the number of observations.
Constraint:
.
- 3:
– const doubleInput
-
On entry: the observations
,
or
, for
.
If the observation is exactly specified – the exact value, .
If the observation is right-censored – the lower value, .
If the observation is left-censored – the upper value, .
If the observation is interval-censored – the lower or upper value,
or
, (see
xc).
- 4:
– const doubleInput
-
On entry: if the
th observation, for is an interval-censored observation then should contain the complementary value to , that is, if , then contains upper value, , and if , then contains lower value, . Otherwise if the th observation is exact or right- or left-censored need not be set.
Note: if then the observation is ignored.
- 5:
– const IntegerInput
-
On entry:
contains the censoring codes for the
th observation, for
.
If , the observation is exactly specified.
If , the observation is right-censored.
If , the observation is left-censored.
If , the observation is interval-censored.
Constraint:
, , or , for .
- 6:
– double *Input/Output
-
On entry: if
the initial estimate of the mean,
; otherwise
xmu need not be set.
On exit: the maximum likelihood estimate, , of .
- 7:
– double *Input/Output
-
On entry: specifies whether an initial estimate of
and
are to be supplied.
- xsig is the initial estimate of and xmu must contain an initial estimate of .
- Initial estimates of xmu and xsig are calculated internally from:
(a) |
the exact observations, if the number of exactly specified observations is ; or |
(b) |
the interval-censored observations; if the number of interval-censored observations is ; or |
(c) |
they are set to and respectively. |
On exit: the maximum likelihood estimate, , of .
- 8:
– doubleInput
-
On entry: the relative precision required for the final estimates of
and
. Convergence is assumed when the absolute relative changes in the estimates of both
and
are less than
tol.
If , a relative precision of is used.
Constraint:
or .
- 9:
– IntegerInput
-
On entry: the maximum number of iterations.
If , a value of is used.
- 10:
– double *Output
-
On exit: the estimate of the standard error of .
- 11:
– double *Output
-
On exit: the estimate of the standard error of .
- 12:
– double *Output
-
On exit: the estimate of the correlation between and .
- 13:
– double *Output
-
On exit: the maximized log-likelihood, .
- 14:
– IntegerOutput
-
On exit: the number of the different types of each observation;
contains number of right-censored observations.
contains number of left-censored observations.
contains number of interval-censored observations.
contains number of exactly specified observations.
- 15:
– Integer *Output
-
On exit: the number of iterations performed.
- 16:
– NagError *Input/Output
-
The NAG error argument (see
Section 3.7 in How to Use the NAG Library and its Documentation).
6
Error Indicators and Warnings
- NE_ALLOC_FAIL
-
Dynamic memory allocation failed.
See
Section 2.3.1.2 in How to Use the NAG Library and its Documentation for further information.
- NE_BAD_PARAM
-
On entry, argument had an illegal value.
- NE_CONVERGENCE
-
Method has not converged in iterations.
- NE_DIVERGENCE
-
Process has diverged.
- NE_EM_PROCESS
-
The EM process has failed.
- NE_INT
-
On entry, .
Constraint: .
- NE_INT_ARRAY
-
On entry, is invalid, it contains .
- NE_INTERNAL_ERROR
-
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact
NAG for assistance.
See
Section 2.7.6 in How to Use the NAG Library and its Documentation for further information.
- NE_NO_LICENCE
-
Your licence key may have expired or may not have been installed correctly.
See
Section 2.7.5 in How to Use the NAG Library and its Documentation for further information.
- NE_OBSERVATIONS
-
On entry, effective number of observations .
- NE_REAL
-
On entry,
tol is invalid:
.
- NE_STANDARD_ERRORS
-
Standard errors cannot be computed.
7
Accuracy
The accuracy is controlled by the argument
tol.
If high precision is requested with the
algorithm then there is a possibility that, due to the slow convergence, before the correct solution has been reached the increments of
and
may be smaller than
tol and the process will prematurely assume convergence.
8
Parallelism and Performance
nag_censored_normal (g07bbc) is not threaded in any implementation.
The process is deemed divergent if three successive increments of or increase.
10
Example
A sample of observations and their censoring codes are read in and the Newton–Raphson method used to compute the estimates.
10.1
Program Text
Program Text (g07bbce.c)
10.2
Program Data
Program Data (g07bbce.d)
10.3
Program Results
Program Results (g07bbce.r)