Find the MLEs of a data set with censoring by using normfit, and then find the negative loglikelihood of the MLEs by using normlike.

Load the sample data.

load lightbulb

The first column of the data contains the lifetime (in hours) of two types of bulbs. The second column contains the binary variable indicating whether the bulb is fluorescent or incandescent. 1 indicates that the bulb is fluorescent, and 0 indicates that the bulb is incandescent. The third column contains the censorship information, where 0 indicates the bulb is observed until failure, and 1 indicates the item (bulb) is censored.

Find the indices for fluorescent bulbs.

idx = find(lightbulb(:,2) == 0);

Find the MLEs of the normal distribution parameters. The second input argument of normfit specifies the confidence level. Pass in [] to use its default value 0.05. The third input argument specifies the censorship information.

censoring = lightbulb(idx,3) == 1;
[muHat,sigmaHat] = normfit(lightbulb(idx,1),[],censoring)

muHat = 9.4966e+03

sigmaHat = 3.0640e+03

Find the negative loglikelihood of the MLEs.

nlogL = normlike([muHat,sigmaHat],lightbulb(idx,1),censoring)

nlogL = 376.2305

The function normfit finds the sample mean and the square root of the unbiased estimator of the variance with no censoring. The sample mean is equal to the MLE of the mean parameter, but the square root of the unbiased estimator of the variance is not equal to the MLE of the standard deviation parameter.

Find the normal distribution parameters by using normfit, convert them into MLEs, and then compare the negative log likelihoods of the estimates by using normlike.

Generate 100 normal random numbers from the standard normal distribution.

rng('default') % For reproducibility
n = 100;
x = normrnd(0,1,[n,1]);

Find the sample mean and the square root of the unbiased estimator of the variance.

[muHat,sigmaHat] = normfit(x)

muHat = 0.1231

sigmaHat = 1.1624

Convert the square root of the unbiased estimator of the variance into the MLE of the standard deviation parameter.

sigmaHat_MLE = sqrt((n-1)/n)*sigmaHat

sigmaHat_MLE = 1.1566

The difference between sigmaHat and sigmaHat_MLE is negligible for large n.

Alternatively, you can find the MLEs by using the function mle.

phat = mle(x)

phat = 1×2

    0.1231    1.1566

phat(1) and phat(2) are the MLEs of the mean and the standard deviation parameter, respectively.

Confirm that the log likelihood of the MLEs (muHat and sigmaHat_MLE) is greater than the log likelihood of the unbiased estimators (muHat and sigmaHat) by using the normlike function.

logL = -normlike([muHat,sigmaHat],x)

logL = -156.4424

logL_MLE = -normlike([muHat,sigmaHat_MLE],x)

logL_MLE = -156.4399

Find the maximum likelihood estimates (MLEs) of the normal distribution parameters, and then find the confidence interval of the corresponding inverse cdf value.

Generate 1000 normal random numbers from the normal distribution with mean 5 and standard deviation 2.

rng('default') % For reproducibility
n = 1000; % Number of samples
x = normrnd(5,2,[n,1]);

Find the MLEs for the distribution parameters (mean and standard deviation) by using mle.

phat = mle(x)

phat = 1×2

    4.9347    1.9969

muHat = phat(1);
sigmaHat = phat(2);

Estimate the covariance of the distribution parameters by using normlike. The function normlike returns an approximation to the asymptotic covariance matrix if you pass the MLEs and the samples used to estimate the MLEs.

[~,pCov] = normlike([muHat,sigmaHat],x)

pCov = 2×2

    0.0040   -0.0000
   -0.0000    0.0020

Find the inverse cdf value at 0.5 and its 99% confidence interval.

[x,xLo,xUp] = norminv(0.5,muHat,sigmaHat,pCov,0.01)

x = 4.9347

xLo = 4.7721

xUp = 5.0974

x is the inverse cdf value using the normal distribution with the parameters muHat and sigmaHat. The interval [xLo,xUp] is the 99% confidence interval of the inverse cdf value evaluated at 0.5, considering the uncertainty of muHat and sigmaHat using pCov. The 99% confidence interval means the probability that [xLo,xUp] contains the true inverse cdf value is 0.99.