1
2R version 3.1.0 (2014-04-10) -- "Spring Dance"
3Copyright (C) 2014 The R Foundation for Statistical Computing
4Platform: x86_64-apple-darwin13.1.0 (64-bit)
5
6R is free software and comes with ABSOLUTELY NO WARRANTY.
7You are welcome to redistribute it under certain conditions.
8Type 'license()' or 'licence()' for distribution details.
9
10R is a collaborative project with many contributors.
11Type 'contributors()' for more information and
12'citation()' on how to cite R or R packages in publications.
13
14Type 'demo()' for some demos, 'help()' for on-line help, or
15'help.start()' for an HTML browser interface to help.
16Type 'q()' to quit R.
17
18> library(survey)
19
20Attaching package: 'survey'
21
22The following object is masked from 'package:graphics':
23
24    dotchart
25
26> library(survival)
27Loading required package: splines
28>
29> pbc2<-rbind(pbc,pbc)
30> pbc2$id<-rep(1:418,2)
31>
32> dpbc1<-svydesign(id=~1, data=pbc)
33Warning message:
34In svydesign.default(id = ~1, data = pbc) :
35  No weights or probabilities supplied, assuming equal probability
36> dpbc2<-svydesign(id=~id, data=pbc2)
37Warning message:
38In svydesign.default(id = ~id, data = pbc2) :
39  No weights or probabilities supplied, assuming equal probability
40>
41> s1<-svykm(Surv(time,status>0)~1, subset(dpbc1, bili>6), se=TRUE)
42> s2<-svykm(Surv(time,status>0)~1, subset(dpbc2, bili>6), se=TRUE)
43>
44> (c1<-confint(s1,(1:5)*365))
45         0.025     0.975
46365  0.6446215 0.8594153
47730  0.5410938 0.7766848
481095 0.2683127 0.5103356
491460 0.1444731 0.3722001
501825 0.1009672 0.3204713
51> (c2<-confint(s2,(1:5)*365))
52         0.025     0.975
53365  0.6446215 0.8594153
54730  0.5410938 0.7766848
551095 0.2683127 0.5103356
561460 0.1444731 0.3722001
571825 0.1009672 0.3204713
58> all.equal(c1, c2)
59[1] TRUE
60>
61> m1<-svycoxph(Surv(time,status>0)~log(bili), design=dpbc1)
62> m2<-svycoxph(Surv(time,status>0)~log(bili), design=dpbc2)
63>
64> d<-data.frame(bili=c(5,10))
65> p1<-predict(m1, se=TRUE, newdata=d,type="curve")
66> p2<-predict(m2, se=TRUE, newdata=d,type="curve")
67>
68> (pc1<-confint(p1[[1]],(1:5)*365))
69         0.025     0.975
70365  0.8410027 0.9266263
71730  0.7371114 0.8548312
721095 0.5517779 0.7018583
731460 0.4335073 0.5992819
741825 0.3260899 0.5046241
75> (pc2<-confint(p2[[1]],(1:5)*365))
76         0.025     0.975
77365  0.8409490 0.9267054
78730  0.7370152 0.8549432
791095 0.5515848 0.7019513
801460 0.4332252 0.5992968
811825 0.3257172 0.5045795
82> all.equal(pc1, pc2)
83[1] "Mean relative difference: 0.0002070722"
84>
85> (q1<-quantile(p1[[2]]))
860.75  0.5 0.25
87 489  930 1492
88> (q2<-quantile(p2[[2]]))
890.75  0.5 0.25
90 489  930 1492
91> all.equal(q1,q2)
92[1] TRUE
93>
94> proc.time()
95   user  system elapsed
96  3.410   0.099   3.519
97