1### R code from vignette source 'epi.Rnw'
2
3###################################################
4### code chunk number 1: epi.Rnw:45-61
5###################################################
6library(survey)
7load(system.file("doc","nwts.rda",package="survey"))
8nwtsnb<-nwts
9nwtsnb$case<-nwts$case-nwtsb$case
10nwtsnb$control<-nwts$control-nwtsb$control
11
12a<-rbind(nwtsb,nwtsnb)
13a$in.ccs<-rep(c(TRUE,FALSE),each=16)
14
15b<-rbind(a,a)
16b$rel<-rep(c(1,0),each=32)
17b$n<-ifelse(b$rel,b$case,b$control)
18index<-rep(1:64,b$n)
19
20nwt.exp<-b[index,c(1:3,6,7)]
21nwt.exp$id<-1:4088
22
23
24###################################################
25### code chunk number 2: epi.Rnw:65-66
26###################################################
27glm(rel~factor(stage)*factor(histol), family=binomial, data=nwt.exp)
28
29
30###################################################
31### code chunk number 3: epi.Rnw:75-79
32###################################################
33dccs2<-twophase(id=list(~id,~id),subset=~in.ccs,
34                strata=list(NULL,~interaction(instit,rel)),data=nwt.exp)
35
36summary(svyglm(rel~factor(stage)*factor(histol),family=binomial,design=dccs2))
37
38
39###################################################
40### code chunk number 4: epi.Rnw:88-94
41###################################################
42dccs8<-twophase(id=list(~id,~id),subset=~in.ccs,
43                strata=list(NULL,~interaction(instit,stage,rel)),data=nwt.exp)
44gccs8<-calibrate(dccs2,phase=2,formula=~interaction(instit,stage,rel))
45
46summary(svyglm(rel~factor(stage)*factor(histol),family=binomial,design=dccs8))
47summary(svyglm(rel~factor(stage)*factor(histol),family=binomial,design=gccs8))
48
49
50###################################################
51### code chunk number 5: epi.Rnw:122-126
52###################################################
53library(survey)
54library(survival)
55data(nwtco)
56ntwco<-subset(nwtco, !is.na(edrel))
57
58
59###################################################
60### code chunk number 6: epi.Rnw:130-131
61###################################################
62coxph(Surv(edrel, rel)~factor(stage)+factor(histol)+I(age/12),data=nwtco)
63
64
65###################################################
66### code chunk number 7: epi.Rnw:143-155
67###################################################
68(dcch<-twophase(id=list(~seqno,~seqno), strata=list(NULL,~rel),
69                  subset=~I(in.subcohort | rel), data=nwtco))
70svycoxph(Surv(edrel,rel)~factor(stage)+factor(histol)+I(age/12),
71                design=dcch)
72
73subcoh <- nwtco$in.subcohort
74selccoh <- with(nwtco, rel==1|subcoh==1)
75ccoh.data <- nwtco[selccoh,]
76ccoh.data$subcohort <- subcoh[selccoh]
77cch(Surv(edrel, rel) ~ factor(stage) + factor(histol) + I(age/12),
78           data =ccoh.data, subcoh = ~subcohort, id=~seqno,
79           cohort.size=4028, method="LinYing")
80
81
82###################################################
83### code chunk number 8: epi.Rnw:165-176
84###################################################
85nwtco$eventrec<-rep(0,nrow(nwtco))
86nwtco.extra<-subset(nwtco, rel==1)
87nwtco.extra$eventrec<-1
88nwtco.expd<-rbind(subset(nwtco,in.subcohort==1),nwtco.extra)
89nwtco.expd$stop<-with(nwtco.expd,
90                      ifelse(rel & !eventrec, edrel-0.001,edrel))
91nwtco.expd$start<-with(nwtco.expd,
92                       ifelse(rel & eventrec, edrel-0.001, 0))
93nwtco.expd$event<-with(nwtco.expd,
94                       ifelse(rel & eventrec, 1, 0))
95nwtco.expd$pwts<-ifelse(nwtco.expd$event, 1, 1/with(nwtco,mean(in.subcohort | rel)))
96
97
98###################################################
99### code chunk number 9: epi.Rnw:185-189
100###################################################
101(dBarlow<-svydesign(id=~seqno+eventrec, strata=~in.subcohort+rel,
102                    data=nwtco.expd, weight=~pwts))
103svycoxph(Surv(start,stop,event)~factor(stage)+factor(histol)+I(age/12),
104                design=dBarlow)
105
106
107###################################################
108### code chunk number 10: epi.Rnw:194-197
109###################################################
110(dWacholder <- as.svrepdesign(dBarlow,type="bootstrap",replicates=500))
111svycoxph(Surv(start,stop,event)~factor(stage)+factor(histol)+I(age/12),
112                design=dWacholder)
113
114
115###################################################
116### code chunk number 11: epi.Rnw:209-217
117###################################################
118load(system.file("doc","nwtco-subcohort.rda",package="survey"))
119nwtco$subcohort<-subcohort
120
121d_BorganII <- twophase(id=list(~seqno,~seqno),
122                       strata=list(NULL,~interaction(instit,rel)),
123                       data=nwtco, subset=~I(rel |subcohort))
124(b2<-svycoxph(Surv(edrel,rel)~factor(stage)+factor(histol)+I(age/12),
125                design=d_BorganII))
126
127
128###################################################
129### code chunk number 12: epi.Rnw:222-225
130###################################################
131d_BorganIIps <- calibrate(d_BorganII, phase=2, formula=~age+interaction(instit,rel,stage))
132svycoxph(Surv(edrel,rel)~factor(stage)+factor(histol)+I(age/12),
133                  design=d_BorganIIps)
134
135
136