Comparison to other software"

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

R code and results

library(eventglm)

colon.cifit <- cumincglm(Surv(time, status) ~ rx, time = 2500, data = colon)
summary(colon.cifit)
confint(colon.cifit)

Stata code and results

This uses the st0202_1 package, available from here: https://www.stata-journal.com/article.html?article=st0202_1

```{stata, include = TRUE, eval = FALSE} . import delimited "colon.csv", clear (18 vars, 929 obs)

. . stset time, failure(status==1)

 failure event:  status == 1

obs. time interval: (0, time] exit on or before: failure


    929  total observations
      0  exclusions

    929  observations remaining, representing
    452  failures in single-record/single-failure data

1,551,389 total analysis time at risk and under observation at risk from t = 0 earliest observed entry t = 0 last observed exit t = 3,329

. . // requires st0202_1 install (search stpci) . stpci, at(2500) Pseudo-observations for the cumulative incidence function. Competing risks: (none). Computing pseudo-observations (progress dots indicate percent completed). ----+--- 1 ---+--- 2 ---+--- 3 ---+--- 4 ---+--- 5 .................................................. 50 .................................................. 100 Generated variable: pseudo.

. tabulate rx, gen(rxdum)

     rx |      Freq.     Percent        Cum.

------------+----------------------------------- Lev | 310 33.37 33.37 Lev+5FU | 304 32.72 66.09 Obs | 315 33.91 100.00 ------------+----------------------------------- Total | 929 100.00

. glm pseudo rxdum1 rxdum2, vce(robust)

Iteration 0: log pseudolikelihood = -708.7556

Generalized linear models Number of obs = 929 Optimization : ML Residual df = 926 Scale parameter = .270145 Deviance = 250.154308 (1/df) Deviance = .270145 Pearson = 250.154308 (1/df) Pearson = .270145

Variance function: V(u) = 1 [Gaussian] Link function : g(u) = u [Identity]

                                              AIC             =   1.532305

Log pseudolikelihood = -708.7556003 BIC = -6078.23


         |               Robust
  pseudo |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]

-------------+---------------------------------------------------------------- rxdum1 | -.029075 .0416186 -0.70 0.485 -.1106459 .0524959 rxdum2 | -.1317578 .0417443 -3.16 0.002 -.2135752 -.0499404 _cons | .5434514 .02938 18.50 0.000 .4858676 .6010352 ------------------------------------------------------------------------------

## SAS code and results

Assuming you have loaded the colon data in your workspace.

```{sas, eval = FALSE}
proc lifetest data=colon noprint plots=none timelist=2500 reduceout outsurv=sall;
time time*status(0);
run;

data sall;
set sall;
theta = survival;
keep theta;
run;

data sout;
set colon;
keep id;
run;

%macro pseudosurv;

%do ip=1 %to 929;
data coloni;
set colon;
where id ^= &ip;
run;

proc lifetest data=coloni noprint plots=none timelist=2500 reduceout outsurv=salli;
time time*status(0);
run;

data salli;
set salli;
thetamini = survival;
id = &ip;
keep id thetamini;
run;

data souti;
merge salli sall;
run;

data sout; 
merge sout souti;
by id;
run;
%end;
%mend pseudosurv;

%pseudosurv;


data sout2;
set sout;
pseudoci = 1 - (929 * theta - (929 - 1) * thetamini);
run;

data colon2;
merge colon sout2;
by id;
if rx='Lev' then rxlev = 1;
else rxlev = 0;
if rx='Lev+5FU' then rxlevplus = 1;
else rxlevplus = 0;
run;

proc reg data = colon2;
model pseudoci = rxlev rxlevplus / white;
run;
Parameter Estimates
Variable DF Parameter
Estimate
Standard
Error
t Value Pr > |t| Heteroscedasticity Consistent
Standard
Error
t Value Pr > |t|
Intercept 1 0.54345 0.02928 18.56 <.0001 0.02936 18.51 <.0001
rxlev 1 -0.02907 0.04158 -0.70 0.4846 0.04160 -0.70 0.4847
rxlevplus 1 -0.13176 0.04179 -3.15 0.0017 0.04172 -3.16 0.0016


Try the eventglm package in your browser

Any scripts or data that you put into this service are public.

eventglm documentation built on April 3, 2025, 6:23 p.m.