[R] R-help: Censoring data (actually an optim issue

[John C Nash]

Your function is giving NaN's during the optimization.

The R-forge version of optimx() has functionality specifically intended to deal with this.
NOTE: the CRAN version does not, and the R-forge version still has some glitches!

However, I easily ran the code you supplied by changing optim to optimx in the penultimate
line. Here's the final output.

KKT condition testing
Number of parameters = 2
max abs gradient element = 0.004032794   test tol =  0.1018794
KKT1 result =  TRUE
Hessian eigenvalues:
[1] 7.138974e+02 9.931285e-04
KKT2 result =  TRUE
KKT results: gmax= 0.004032794   evratio= 1.391136e-06   KKT1 & 2:  TRUE TRUE
[1] 7.138974e+02 9.931285e-04
Save results from method  BFGS

> zz
  method                   par  fvalues fns grs hes rs conv KKT1 KKT2
2   BFGS 0.2832468, 52.6096161 100.8794  10   1   0 -2    3 TRUE TRUE
1     NM 0.2827563, 54.3551491 100.8779  53   0   0  1    0 TRUE TRUE
         mtilt xtimes meths
2           NA   0.11  BFGS
1 0.0005130808   0.57    NM
There were 50 or more warnings (use warnings() to see the first 50)
>

Note the Hessian eigenvalues are pretty bad. And the warnings are the NaN's produced.
This is with just the 2 default methods. When I tried "all.methods", one (uobyqa) seemed
to lock up. This is a fairly ill-conditioned problem.


Best, JN


On 04/14/2012 06:00 AM, r-help-request at r-project.org wrote:
> Message: 13
> Date: Fri, 13 Apr 2012 03:54:43 -0700 (PDT)
> From: Christopher Kelvin <chris_kelvin2001 at yahoo.com>
> To: "r-help at r-project.org" <r-help at r-project.org>
> Subject: [R] R-Help: Censoring data
> Message-ID:
> 	<1334314483.27693.YahooMailNeo at web65408.mail.ac4.yahoo.com>
> Content-Type: text/plain; charset=iso-8859-1
> 
> Hello,
> ?I want to estimate weibull parameters with 30% censored data. I have below the code for the censoring
> ?but how it must be put into the likelihood equation to obtain the desire estimate is where i have a problem with,
> ?can some body help?
> ?My likelihood equation is for a random type-I censoring where time for the censored units is different for each unit.
> ?
> n=50;r=35
> p=0.8;b=1.5
> t<-rweibull(50,shape=p,scale=b)
> meantrue<-gamma(1+(1/p))*b
> meantrue
> d<-meantrue/0.30
> 
> cen<- runif(50,min=0,max=d)
> cen
> s<-ifelse(t<=cen,1,0)
> s
> 
> z<-function(p){?
> shape<-p[1]
> scale<-p[2]
> log1<-(r*log(p[1])-r*(p[1])*log(p[2])+(p[1]-1)*sum(log(t))-sum((t/(p[2]))^(p[1])
> )-((n-r)*(sum(cen)/(p[2]))^(p[1])))
> return(-log1)
> }
> 
> start <- c(1,1)
> zz<-optim(start,fn=z,hessian=T)
> zz
> 
> Thanks in anticipation
> 
> Chris Guure
> Researcher
> Institute for Mathematical Research
> UPM