Cross Validated is a question and answer site for people interested in statistics, machine learning, data analysis, data mining, and data visualization. It only takes a minute to sign up.
Sign up to join this community
Note: I am posting a question from a former student of mine unable to post on his own for technical reasons.
Given an iid sample $x_1,\ldots,x_n$ from a Weibull distribution with pdf $$ f_k(x) = k x^{k-1} e^{-x^k} \quad x>0 $$ is there a useful missing variable representation $$f_k(x) = \int_\mathcal{Z} g_k(x,z)\,\text{d}z$$and hence an associated EM (expectation-maximisation) algorithm that could be used to find the MLE of $k$, instead of using straightforward numerical optimisation?
I think the answer is yes, if I have understood the question correctly.
Write $z_i = x_i^k$. Then an EM algorithm type of iteration, starting with for example $\hat k = 1$, is
E step: ${\hat z}_i = x_i^{\hat k}$
M step: $\hat k = \frac{n}{\left[\sum({\hat z}_i - 1)\log x_i\right]}$
This is a special case (the case with no censoring and no covariates) of the iteration suggested for Weibull proportional hazards models by Aitkin and Clayton (1980). It can also be found in Section 6.11 of Aitkin et al (1989).
Aitkin, M. and Clayton, D., 1980. The fitting of exponential, Weibull and extreme value distributions to complex censored survival data using GLIM. Applied Statistics, pp.156-163.
Aitkin, M., Anderson, D., Francis, B. and Hinde, J., 1989. Statistical Modelling in GLIM. Oxford University Press. New York.
The Weibull MLE is only numerically solvable:
Let $$ f_{\lambda,\beta}(x) = \begin{cases} \frac{\beta}{\lambda}\left(\frac{x}{\lambda}\right)^{\beta-1}e^{-\left(\frac{x}{\lambda}\right)^{\beta}} & ,\,x\geq0 \\ 0 &,\, x<0 \end{cases} $$ with $\beta,\,\lambda>0$.
1) Likelihoodfunction: $$ \mathcal{L}_{\hat{x}}(\lambda, \beta) =\prod_{i=1}^N f_{\lambda,\beta}(x_i) =\prod_{i=1}^N \frac{\beta}{\lambda}\left(\frac{x_i}{\lambda}\right)^{\beta-1}e^{-\left(\frac{x_i}{\lambda}\right)^{\beta}} = \frac{\beta^N}{\lambda^{N \beta}} e^{-\sum_{i=1}^N\left(\frac{x_i}{\lambda}\right)^{\beta}} \prod_{i=1}^N x_i^{\beta-1} $$
log-Likelihoodfunction: $$ \ell_{\hat{x}}(\lambda, \beta):= \ln \mathcal{L}_{\hat{x}}(\lambda, \beta)=N\ln \beta-N\beta\ln \lambda-\sum_{i=1}^N \left(\frac{x_i}{\lambda}\right)^\beta+(\beta-1)\sum_{i=1}^N \ln x_i $$
2) MLE-Problem: \begin{equation*} \begin{aligned} & & \underset{(\lambda,\beta) \in \mathbb{R}^2}{\text{max}}\,\,\,\,\,\, & \ell_{\hat{x}}(\lambda, \beta) \\ & & \text{s.t.} \,\,\, \lambda>0\\ & & \beta > 0 \end{aligned} \end{equation*} 3) Maximization by $0$-gradients: \begin{align*} \frac{\partial l}{\partial \lambda}&=-N\beta\frac{1}{\lambda}+\beta\sum_{i=1}^N x_i^\beta\frac{1}{\lambda^{\beta+1}}&\stackrel{!}{=} 0\\ \frac{\partial l}{\partial \beta}&=\frac{N}{\beta}-N\ln\lambda-\sum_{i=1}^N \ln\left(\frac{x_i}{\lambda}\right)e^{\beta \ln\left(\frac{x_i}{\lambda}\right)}+\sum_{i=1}^N \ln x_i&\stackrel{!}{=}0 \end{align*} It follows: \begin{align*} -N\beta\frac{1}{\lambda}+\beta\sum_{i=1}^N x_i^\beta\frac{1}{\lambda^{\beta+1}} &= 0\\\\ -\beta\frac{1}{\lambda}N +\beta\frac{1}{\lambda}\sum_{i=1}^N x_i^\beta\frac{1}{\lambda^{\beta}} &= 0\\\\ -1+\frac{1}{N}\sum_{i=1}^N x_i^\beta\frac{1}{\lambda^{\beta}}&=0\\\\ \frac{1}{N}\sum_{i=1}^N x_i^\beta&=\lambda^\beta \end{align*} $$\Rightarrow\lambda^*=\left(\frac{1}{N}\sum_{i=1}^N x_i^{\beta^*}\right)^\frac{1}{\beta^*}$$
Plugging $\lambda^*$ into the second 0-gradient condition:
\begin{align*} \Rightarrow \beta^*=\left[\frac{\sum_{i=1}^N x_i^{\beta^*}\ln x_i}{\sum_{i=1}^N x_i^{\beta^*}}-\overline{\ln x}\right]^{-1} \end{align*}
This equation is only numerically solvable, e.g. Newton-Raphson algorithm. $\hat{\beta}^*$ can then be placed into $\lambda^*$ to complete the ML estimator for the Weibull distribution.
Though this is an old question, it looks like there is an answer in a paper published here: http://home.iitk.ac.in/~kundu/interval-censoring-REVISED-2.pdf
In this work the analysis of interval-censored data, with Weibull distribution as the underlying lifetime distribution has been considered. It is assumed that censoring mechanism is independent and non-informative. As expected, the maximum likelihood estimators cannot be obtained in closed form. In our simulation experiments it is observed that the Newton-Raphson method may not converge many times. An expectation maximization algorithm has been suggested to compute the maximum likelihood estimators, and it converges almost all the times.
In this case the MLE and EM estimators are equivalent, since the MLE estimator is actually just a special case of the EM estimator. (I am assuming a frequentist framework in my answer; this isn't true for EM in a Bayesian context in which we're talking about MAP's). Since there is no missing data (just an unknown parameter), the E step simply returns the log likelihood, regardless of your choice of $k^{(t)}$. The M step then maximizes the log likelihood, yielding the MLE.
EM would be applicable, for example, if you had observed data from a mixture of two Weibull distributions with parameters $k_1$ and $k_2$, but you didn't know which of these two distributions each observation came from.
