Example 12.2 Maximum Likelihood Weibull Model

This model is obtained from Lawless (1982, pp. 190–194). Suppose you want to find the maximum likelihood estimates for the three-parameter Weibull model. The observed likelihood function is

where is the number of individuals involved in the experiment, is the set of individuals whose lifetimes are observed, and . Then the log-likelihood function is

Note that for , the logarithmic terms become infinite as . That is, is unbounded. Thus our interest is restricted to values greater than or equal to . Further, for the logarithmic terms to be defined, it is required that and .

The following data from Pike (1966) represent the number of days it took the rats painted with the carcinogen DMBA to develop carcinomas:

Note that the last two observations have an asterisk next to them. This indicates that those two rats did not develop carcinomas. These are referred to as censored data.

For these data, a local maximum of occurs at , , and . You can create the input data set by using the following SAS code:

   data pike;
      input days cens @@;
      datalines;
   143  0  164  0  188  0  188  0
   190  0  192  0  206  0  209  0
   213  0  216  0  220  0  227  0
   230  0  234  0  246  0  265  0
   304  0  216  1  244  1
   ;

The first observations in the data set correspond to the rats that develop carcinomas; in this case. Using the preceding data, you can formulate the Weibull model as follows:

proc optmodel;
   set S;  /* set of rats in the experiment */
   set D;  /* set of rats that develop carcinomas */
   number r = card(D);
   number t{S};
   number cens{S};
   var alpha >= 0 init 1;
   var beta  >= 1 init 1;
   var mu >= 0 <= min{i in D} t[i] init 10;
   max logl = r*log(beta) - r*beta*log(alpha)
                + (beta - 1)*sum{i in D} log(t[i] - mu)
                - sum{i in S} ((t[i] - mu)/alpha)^beta;
   read data pike into S=[_n_] t=days cens;
   D = {i in S : cens[i] = 0};
   solve with nlpc / printfreq=1;
   print mu alpha beta;
quit;   

Assume a starting point at , , and . The data for all the rats in the experiment are read using the READ statement. The following statement subsets to the set of rats that develop carcinomas:

      D = {i in S : cens[i] = 0};

The NLPC solver invokes the default optimization technique (the trust region method) to solve this problem.

Output 12.2.1 through Output 12.2.3 show the problem summary, iteration log, solution summary, and solution. As you can see, the solution obtained by the trust region method matches closely with the local maximum given in Lawless (1982, p. 193).