TBSSurvival (version 1.0)

tbs.survreg.mle: MLE of the TBS Model for Failure Data

Description

This function perform the Maximum Likelihhod Estimation of the TBS model. The optimization is done by the function `optim' (or the package Rsolnp when available).

Usage

tbs.survreg.mle(formula,dist="norm",
                 method=c("BFGS","Rsolnp","Nelder-Mead","CG","SANN"),
                 verbose=FALSE,nstart=10,max.time=-1,ncore=1)

Arguments

formula

A formula specification containing a Surv model with right-censored data as in the package survival.

dist

error distribution; dist = "norm", "doubexp", "t", "cauchy" or "logistic".

method

a vector of numerical methods to be used in the optimization. The function try all listed methods and returns the solution with maximal likelihood among them.

verbose

Boolean to indicate the amount of output during the execution of the optimization.

nstart

Number of feasible initial points to guess when performing the optimization.

max.time

Maximum time (in minutes) to run the optimization (

ncore

If parallel computing is available, ncore tells the number of cores to use.

Value

• A list with components:
• par
  The best set of parameters found, the first value of par is the estimate of lambda, the second value is the estimate of xi and the others are the beta parameter.
• std.error
  Standar error for MLE. If it is not possible to evaluate return NA
• log.lik
  The log-likelihood at parameters par.
• error.dist
  The error distribution chosen.
• AIC
  Akaike Information Criterion.
• AICc
  AICc is AIC with a second order correction for small sample sizes.
• BIC
  Bayesian Information Criterion.
• KS
  Komogorov-Simirnov Statistic. It is not evaluated with matrix x of co-variables.
• method
  Numerical method used to achive the MLE.
• convergence
  If convergence is FALSE then it was not possible to find the MLE.
• time
  observed survival times.
• error
  error of the estimated model.
• call
  function evaluated.
• formula
  formula entered by user.
• run.time
  Time spent with the function.

Details

This function calls numerical optimization methods to maximize the likelihood of the TBS model, according to the given error distribution, method of optimization, and formula. The formula is supposed to have a Surv object and possibility co-variates, just as the standard specification of R formulas. The optimizers are going to do their best to find high likelihood estimates, but as in most estimation methods that need a numerical optimization procedure, the obtained estimate cannot be guaranteed to be a global optimal solution, but instead is dependent on the initial points, and thus on the seed of the random number generation.

References

Meeker, W. and Escobar, L. (1998) Statistical Methods for Reliability Data. Willey, ISBN 0-471-14328-6.

See Also

tbs.survreg.be,dtbs,ptbs,qtbs,rtbs.

Examples

# Alloy - T7987: data extracted from Meeker and Escobar (1998), pp. 131.
data(alloyT7987)
alloyT7987$time <- as.double(alloyT7987$time)
alloyT7987$delta <- as.double(alloyT7987$delta)

# MLE estimation with logistic error
formula <- Surv(alloyT7987$time,alloyT7987$delta == 1) ~ 1
tbs.mle <- tbs.survreg.mle(formula,dist="logistic",method="Nelder-Mead",nstart=3)

# Kaplan-Meier estimation
km <- survfit(formula)

plot(km,ylab="",xlab="",xlim=c(min(alloyT7987$time),max(alloyT7987$time)),
     conf.int=FALSE,axes=FALSE,lty=1,lwd=1)
t <- seq(min(alloyT7987$time),max(alloyT7987$time),
         (max(alloyT7987$time)-min(alloyT7987$time)-0.01)/1000)
title(ylab="R(t)",xlab="t: number of cycles (in thousands)",
      main="Reliability function (MLE)",
      sub="Alloy - T7987 (cf. Meeker and Escobar (1998), pp. 131)",
      cex.lab=1.2)
axis(1,at=c(93,100,150,200,250,300),lwd=1,lwd.ticks=1,pos=0)
axis(2,lwd=1,lwd.ticks=1,pos=min(alloyT7987$time))
legend(200,0.95,c("Kaplan-Meier",
                  expression(textstyle(paste("TBS / ",sep="")) ~ epsilon
                             ~ textstyle(paste("~",sep="")) ~ Logistic)),
       col=c(1,2),lty=c(1,1),cex=1.1,lwd=c(1,2),bg="white")
lines(t,1-ptbs(t,lambda=tbs.mle$par[1],xi=tbs.mle$par[2],
               beta=tbs.mle$par[3],dist=tbs.mle$error.dist),type="l",
               lwd=2,col=2,lty=1)