R : Copyright 2005, The R Foundation for Statistical Computing Version 2.1.1 (2005-06-20), ISBN 3-900051-07-0 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for a HTML browser interface to help. Type 'q()' to quit R. > ### *
> ### > attach(NULL, name = "CheckExEnv") > assign(".CheckExEnv", as.environment(2), pos = length(search())) # base > ## add some hooks to label plot pages for base and grid graphics > setHook("plot.new", ".newplot.hook") > setHook("persp", ".newplot.hook") > setHook("grid.newpage", ".gridplot.hook") > > assign("cleanEx", + function(env = .GlobalEnv) { + rm(list = ls(envir = env, all.names = TRUE), envir = env) + RNGkind("default", "default") + set.seed(1) + options(warn = 1) + delayedAssign("T", stop("T used instead of TRUE"), + assign.env = .CheckExEnv) + delayedAssign("F", stop("F used instead of FALSE"), + assign.env = .CheckExEnv) + sch <- search() + newitems <- sch[! sch %in% .oldSearch] + for(item in rev(newitems)) + eval(substitute(detach(item), list(item=item))) + missitems <- .oldSearch[! .oldSearch %in% sch] + if(length(missitems)) + warning("items ", paste(missitems, collapse=", "), + " have been removed from the search path") + }, + env = .CheckExEnv) > assign("..nameEx", "__{must remake R-ex/*.R}__", env = .CheckExEnv) # for now > assign("ptime", proc.time(), env = .CheckExEnv) > grDevices::postscript("dblcens-Examples.ps") > assign("par.postscript", graphics::par(no.readonly = TRUE), env = .CheckExEnv) > options(contrasts = c(unordered = "contr.treatment", ordered = "contr.poly")) > options(warn = 1) > library('dblcens') > > assign(".oldSearch", search(), env = .CheckExEnv) > assign(".oldNS", loadedNamespaces(), env = .CheckExEnv) > cleanEx(); ..nameEx <- "d011" > > ### * d011 > > flush(stderr()); flush(stdout()) > > ### Name: d011 > ### Title: Compute NPMLE of CDF from doubly censored data > ### Aliases: d011 > ### Keywords: survival nonparametric > > ### ** Examples > > d011(z=c(1,2,3,4,5), d=c(1,0,2,2,1)) $time [1] 1.0 2.0 2.5 5.0 $status [1] 1 0 -1 1 $surv [1] 0.5000351 0.5000351 0.3333177 0.0000000 $jump [1] 0.4999649 0.0000000 0.1667174 0.3333177 $exttime [1] 1.0 2.0 2.5 3.0 4.0 5.0 $extstatus [1] 1 0 -1 2 2 1 $extweight [1] 1 1 0 1 1 1 $extjump [1] 0.4999649 0.0000000 0.1667174 0.0000000 0.0000000 0.3333177 $extsurv.Sx [1] 0.5000351 0.5000351 0.3333177 0.3333177 0.3333177 0.0000000 $surv0.Sy [1] 1.0000000 0.6000281 0.6000281 0.6000281 0.6000281 0.6000281 $jump0 [1] 0.0000000 0.3999719 0.0000000 0.0000000 0.0000000 0.0000000 $surv2.Sz [1] 0.599986 0.599986 0.599986 0.599986 0.299993 0.000000 $jump2 [1] 0.000000 0.000000 0.000000 0.299993 0.299993 0.000000 $conv [1] 3.300000e+01 8.788214e-06 > # > # you should get something like below (and more) > # > # $time: > # [1] 1.0 2.0 2.5 5.0 (notice the times, (3,4), corresponding > # to d=2 are removed, and time 2.5 added > # $status: since there is a (0,2) pattern at > # [1] 1 0 -1 1 times 2, 3. The status indicator of -1 > # show that it is an added time ) > # $surv > # [1] 0.5000351 0.5000351 0.3333177 0.0000000 > # > # $jump > # [1] 0.4999649 0.0000000 0.1667174 0.3333177 > # > # $exttime > # [1] 1.0 2.0 2.5 3.0 4.0 5.0 > # > # $extstatus > # [1] 1 0 -1 2 2 1 > # > # ...... > # > # $conv > # [1] 3.300000e+01 8.788214e-06 ### did 33 iterations > # > # BTW, the true NPMLE of surv is (1/2, 1/2, 1/3, 0) at times (1,2,2.5,5). > ###### Example 2. > d011(c(1,2,3,4,5), c(1,2,1,0,1),influence.fun=TRUE) $time [1] 1 3 4 5 $status [1] 1 1 0 1 $surv [1] 0.6 0.4 0.4 0.0 $jump [1] 0.4 0.2 0.0 0.4 $exttime [1] 1 2 3 4 5 $extstatus [1] 1 2 1 0 1 $extweight [1] 1 1 1 1 1 $extjump [1] 0.4 0.0 0.2 0.0 0.4 $extsurv.Sx [1] 0.6 0.6 0.4 0.4 0.0 $surv0.Sy [1] 1.0 1.0 1.0 0.5 0.5 $jump0 [1] 0.0 0.0 0.0 0.5 0.0 $surv2.Sz [1] 0.5 0.5 0.0 0.0 0.0 $jump2 [1] 0.0 0.5 0.0 0.0 0.0 $conv [1] 3 0 $Nodes [1] 2 4 $NodeStatus [1] 2 0 $IC1tu [,1] [,2] [1,] -1 0 [2,] -1 -2 $IC1tu2 [,1] [,2] [,3] [1,] -1 0 0 [2,] -1 -1 0 $IC2tu [,1] [,2] [1,] 0.0000000 0 [2,] -0.3333333 0 $IC3tu [,1] [,2] [1,] -1 -0.6666667 [2,] -1 -1.0000000 $VarFt [1] 0.24 0.24 > # we get > # ...... > #$conv: > #[1] 3 0 > # > #$Nodes: > #[1] 2 4 > # > #$IC1tu: > # [,1] [,2] > #[1,] -1 0 > #[2,] -1 -2 > # > #$IC2tu: > # [,1] [,2] > #[1,] 0.0000000 0 > #[2,] -0.3333333 0 > # > #$IC3tu: > # [,1] [,2] > #[1,] -1 -0.6666667 > #[2,] -1 -1.0000000 > # > #$VarFt: > #[1] 0.24 0.24 ## est var of F(t) at t=nodes > ####################################################### > > > > cleanEx(); ..nameEx <- "d011ch" > > ### * d011ch > > flush(stderr()); flush(stdout()) > > ### Name: d011ch > ### Title: Compute NPMLE of CDF from doubly censored data, with a > ### constraint > ### Aliases: d011ch > ### Keywords: survival nonparametric > > ### ** Examples > > d011ch(z=c(1,2,3,4,5), d=c(1,0,2,2,1), K=3.5, konst=0.6) $time [1] 1.0 2.0 2.5 5.0 $status [1] 1 0 -1 1 $surv [1] 0.5000351 0.5000351 0.3333177 0.0000000 $jump [1] 0.4999649 0.0000000 0.1667174 0.3333177 $exttime [1] 1.0 2.0 2.5 3.0 4.0 5.0 $extstatus [1] 1 0 -1 2 2 1 $extjump [1] 0.4999649 0.0000000 0.1667174 0.0000000 0.0000000 0.3333177 $extsurv.Sx [1] 0.5000351 0.5000351 0.3333177 0.3333177 0.3333177 0.0000000 $konstdist [1] 0.4999365 0.4999365 0.6000000 0.6000000 0.6000000 1.0000000 $konstjump [1] 0.4999365 0.0000000 0.1000635 0.0000000 0.0000000 0.4000000 $konsttime [1] 3.5 $theta [1] 0.6 $"-2loglikR" [1] 0.05679897 $maxiter [1] 33 > # > # Here we are testing Ho: F(3.5) = 0.6 with a two-sided alternative > # you should get something like > # > # $time: > # [1] 1.0 2.0 2.5 5.0 (notice the times, (3,4), corresponding > # to d=2 are removed, and time 2.5 added > # $status: since there is a (0,2) pattern at > # [1] 1 0 -1 1 times 2, 3. The status indicator of -1 > # show that it is an added time ) > # $surv > # [1] 0.5000351 0.5000351 0.3333177 0.0000000 > # > # $jump > # [1] 0.4999649 0.0000000 0.1667174 0.3333177 > # > # $exttime > # [1] 1.0 2.0 2.5 3.0 4.0 5.0 (exttime include all the times, > # censor or not, plus the added time) > # $extstatus > # [1] 1 0 -1 2 2 1 > # > # $extjump > # [1] 0.4999649 0.0000000 0.1667174 0.0000000 0.0000000 0.3333177 > # > # $extsurv.Sx > # [1] 0.5000351 0.5000351 0.3333177 0.3333177 0.3333177 0.0000000 > # > # $konstdist > # [1] 0.4999365 0.4999365 0.6000000 0.6000000 0.6000000 1.0000000 > # > # $konstjump > # [1] 0.4999365 0.0000000 0.1000635 0.0000000 0.0000000 0.4000000 > # > # $konsttime > # [1] 3.5 > # > # $theta > # [1] 0.6 > # > # $"-2loglikR" (the Wilks statistics to test Ho: > # [1] 0.05679897 F(K)=konst) > # > # $maxiter > # [1] 33 > # > # The Wilks statistic is only 0.05679897, there is no evidence agaist Ho > > > > ### *