New rootSolve causes failure

[yfyang86]

Verison: 0.2+
BUG: root solving procedure causes -2LLR to be negative.

The potential issue is that rootSolve package changes in the past.

emplik results:

library(KMsurv)
library(emplik)
data(kidney)
RMSTfun <- function(x) {pmin(x,20) - 15}
el.cen.EM2( x  = kidney$time, d = kidney$delta, fun = RMSTfun, mu = 0, maxit = 50)

reports

$loglik
[1] -131.105

$times
 [1]  0.5  0.5  0.5  0.5  0.5  0.5  1.5  2.5  2.5  3.5  3.5  4.5  4.5  5.5  6.5
[16]  8.5  8.5  9.5 10.5 11.5 15.5 15.5 16.5 18.5 23.5 26.5 28.5

$prob
 [1] 0.01122308 0.01122308 0.01122308 0.01122308 0.01122308 0.01122308
 [7] 0.01243132 0.01279473 0.01279473 0.01381316 0.01381316 0.01533046
[13] 0.01533046 0.01599773 0.01746515 0.02043914 0.02043914 0.02186346
[19] 0.02285866 0.02485554 0.03538362 0.03538362 0.03297763 0.03220609
[25] 0.06200016 0.10114908 0.39333447

$lam
[1] 2.061926

$iters
[1] 50

$`-2LLR`
[1] 2.51874

$Pval
[1] 0.1125004

kmc results

library(kmc)
kmc.solve(kidney$time, kidney$delta, list(RMSTfun))
myfun9 <- function(theta, x, d){
    ff <- function(t) pmin(t, 20) - theta
kmc.solve(x, d, g = list(ff)) 
}
findUL(step = 1, fun = myfun9, MLE = 16, x = kidney$time, d = kidney$delta)

$Low
[1] 12.97085

$Up
[1] 19.56962

$FstepL
[1] 6.103516e-05

$FstepU
[1] 8.69197e-05

$Lvalue
[1] 3.840001

$Uvalue
[1] -198.9877

Warning messages:
1: In kmc.solve(x, d, g = list(ff)) :
The results may be not feasible!

2: In kmc.solve(x, d, g = list(ff)) :
The results may be not feasible!

3: In kmc.solve(x, d, g = list(ff)) :
The results may be not feasible!

4: In kmc.solve(x, d, g = list(ff)) :
The results may be not feasible!

5: In kmc.solve(x, d, g = list(ff)) :
The results may be not feasible!

6: In kmc.solve(x, d, g = list(ff)) :
The results may be not feasible!

7: In kmc.solve(x, d, g = list(ff)) :
The results may be not feasible!

[yfyang86]

After revision, in this case, lambda should meet the following requirements:

1. sum(prob) = 1

therefore, for the rootSolve(.) function. we need to check the root that:

1. in the branch that contains 0
2. satisfies sum(RMSTfun(time) * prob) == 1

An older version checks this, but the newer one fails.

The next version will fix this.

library(KMsurv)
library(emplik)
library(kmc)
library(rootSolve)
data(kidney)
kmc.clean <- function(kmc.time,delta){
  n <- length(kmc.time)
  dataOrder <- order(kmc.time, -delta)
  kmc.time <- kmc.time[dataOrder]
  delta <- delta[dataOrder] 
  FirstUnCenLocation<-which(delta==1)[1];
  if (FirstUnCenLocation==n) {stop('Only one uncensored point.');}
  if (FirstUnCenLocation!=1){
    delta=delta[FirstUnCenLocation:n];
    kmc.time=kmc.time[FirstUnCenLocation:n];
  }
  delta[length(kmc.time)]=1;
  return (list(kmc.time=kmc.time,delta=delta));
}

RMSTfun <- function(x) {pmin(x,20) - 15}

x = kidney$time
d = kidney$delta
inputData = kmc.clean(x, d)
Gmat = RMSTfun(inputData$kmc.time)

##  `kmc_native` is the KMC Solver (1-dim) which returns `prob`
fun.C <- Vectorize(
  function(lam){
    return(sum(Gmat*kmc_native(inputData$delta, lam*Gmat)));
  }
)


All = uniroot.all(fun.C, c(-3, 2), tol = 1e-8)
s = el.cen.EM2( x  = kidney$time, d = kidney$delta, fun = RMSTfun, mu = 0, maxit = 50)

curve(fun.C(x), -3, 5,main = "uniroot.all", ylim=c(-10,10))
points(All, y = rep(0, length(All)), pch = 16, cex = 2, col='steelblue')
points(-s$lam, y = 0, pch = 16, cex = 2, col='red')

[yfyang86]

@yfyang86 26abbb3

Fixed