R Under development (unstable) (2024-03-19 r86153 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

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 an HTML browser interface to help.
Type 'q()' to quit R.

> library(survey)
Loading required package: grid
Loading required package: Matrix
Loading required package: survival

Attaching package: 'survey'

The following object is masked from 'package:graphics':

    dotchart

> library(survival)
> data(nwtco)
> 
> ntwco<-subset(nwtco, !is.na(edrel))
> 
> load("nwtco-subcohort.rda")
> nwtco$subcohort<-subcohort
> 
> d_BorganII <- twophase(id=list(~seqno,~seqno),
+                        strata=list(NULL,~interaction(instit,rel)),
+                        data=nwtco, subset=~I(rel |subcohort))
> 
> ##Coefficient results same as Splus with code from
> ## http://faculty.washington.edu/norm/software.html
> ## SE slightly larger due to using sandwich variance.
> 
> svycoxph(Surv(edrel, rel)~factor(stage)+factor(histol)+I(age/12), design=d_BorganII)
Call:
svycoxph(formula = Surv(edrel, rel) ~ factor(stage) + factor(histol) + 
    I(age/12), design = d_BorganII)

                   coef exp(coef) se(coef) robust se      z        p
factor(stage)2  0.46286   1.58861  0.23762   0.18087  2.559  0.01049
factor(stage)3  0.58309   1.79156  0.23965   0.17848  3.267  0.00109
factor(stage)4  1.05967   2.88541  0.26182   0.20524  5.163 2.43e-07
factor(histol)2 1.59744   4.94035  0.17688   0.13342 11.973  < 2e-16
I(age/12)       0.02994   1.03039  0.02942   0.03337  0.897  0.36972

Likelihood ratio test=  on 5 df, p=
n= 1062, number of events= 571 
> 
> ##
> ## This gives higher standard errors. calibrate() does not recompute the
> ## finite population correction if a calibration variable happens to predict
> ## sampling perfectly. It probably should.
> ##
> d_BorganIIps<-calibrate(twophase(id=list(~seqno,~seqno),
+                            strata=list(NULL,~rel),
+                            data=nwtco, subset=~I(rel |subcohort)), 
+                         phase=2, formula=~interaction(instit,rel),
+                         epsilon=1e-10)
> 
> svycoxph(Surv(edrel, rel)~factor(stage)+factor(histol)+I(age/12), design=d_BorganIIps)
Call:
svycoxph(formula = Surv(edrel, rel) ~ factor(stage) + factor(histol) + 
    I(age/12), design = d_BorganIIps)

                   coef exp(coef) se(coef) robust se      z        p
factor(stage)2  0.46286   1.58861  0.23762   0.18076  2.561  0.01045
factor(stage)3  0.58309   1.79156  0.23965   0.17838  3.269  0.00108
factor(stage)4  1.05967   2.88541  0.26182   0.20513  5.166 2.39e-07
factor(histol)2 1.59744   4.94035  0.17688   0.13044 12.246  < 2e-16
I(age/12)       0.02994   1.03039  0.02942   0.03331  0.899  0.36876

Likelihood ratio test=  on 5 df, p=
n= 1062, number of events= 571 
> 
> proc.time()
   user  system elapsed 
   2.59    0.57    3.15