A análise de Sobrevivência é, portanto, um conjunto de técnicas estatísticas desenvolvidas para estudar dados censurados. Seu principal campo de aplicação é em situações médicas envolvendo dados censurados. No entanto, ocorrem situações similares em que a aplicação destas técnicas é adequada, especialmente em engenharia, ciências sociais, demografia.
A grosso modo, a análise descritiva consiste em apresentar medidas de tendência central e variabilidade dos dados utilizados no estudo. No entanto, a presença de censura nos dados inviabiliza este procedimento. Assim, em dados de tempo de vida, o principal método de análise exploratória é a representação da função de sobrevivência.
Para dados não censurados, a função de sobrevivência empírica é dada pela divisão entre o número de observações que não falharam até o tempo t e o total de observações no estudo. Para dados censurados, podemos utilizar a técnica não-paramétrica conhecida como Estimador de Kaplan-Meier.
Os dados da tabela abaixo referem-se aos tempos de sobrevivência (em dias) de pacientes com câncer submetidos à radioterapia (o símbolo + indica censura).
dados_t=c(7,34,42,63,64,-74,83,84,91,108,112,129,133,133,139,140,
140,146,149,154,160,160,165,173,176,-185,218,225,241,
248,273,277,-279,297,-319,405,417,420,440,-523,
583,594,1101,-1116,1146,-1226,-1349,-1412,1417)
cens = ifelse(dados_t<= 0,0,1)
tempos = abs(dados_t)
library(survival)
ekm=survfit(Surv(tempos,cens)~1)
summary(ekm)
## Call: survfit(formula = Surv(tempos, cens) ~ 1)
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 7 49 1 0.980 0.0202 0.9408 1.000
## 34 48 1 0.959 0.0283 0.9054 1.000
## 42 47 1 0.939 0.0342 0.8740 1.000
## 63 46 1 0.918 0.0391 0.8448 0.998
## 64 45 1 0.898 0.0432 0.8171 0.987
## 83 43 1 0.877 0.0470 0.7896 0.974
## 84 42 1 0.856 0.0503 0.7630 0.961
## 91 41 1 0.835 0.0532 0.7372 0.946
## 108 40 1 0.814 0.0559 0.7120 0.932
## 112 39 1 0.794 0.0582 0.6873 0.916
## 129 38 1 0.773 0.0603 0.6631 0.900
## 133 37 2 0.731 0.0639 0.6159 0.867
## 139 35 1 0.710 0.0654 0.5928 0.850
## 140 34 2 0.668 0.0679 0.5476 0.815
## 146 32 1 0.647 0.0689 0.5255 0.797
## 149 31 1 0.626 0.0698 0.5037 0.779
## 154 30 1 0.606 0.0705 0.4821 0.761
## 160 29 2 0.564 0.0715 0.4397 0.723
## 165 27 1 0.543 0.0719 0.4189 0.704
## 173 26 1 0.522 0.0721 0.3983 0.684
## 176 25 1 0.501 0.0722 0.3780 0.665
## 218 23 1 0.479 0.0722 0.3568 0.644
## 225 22 1 0.458 0.0722 0.3359 0.623
## 241 21 1 0.436 0.0719 0.3153 0.602
## 248 20 1 0.414 0.0716 0.2950 0.581
## 273 19 1 0.392 0.0710 0.2750 0.559
## 277 18 1 0.370 0.0704 0.2553 0.538
## 297 16 1 0.347 0.0697 0.2344 0.515
## 405 14 1 0.322 0.0690 0.2121 0.490
## 417 13 1 0.298 0.0680 0.1903 0.466
## 420 12 1 0.273 0.0667 0.1690 0.441
## 440 11 1 0.248 0.0651 0.1483 0.415
## 583 9 1 0.221 0.0634 0.1255 0.387
## 594 8 1 0.193 0.0612 0.1036 0.359
## 1101 7 1 0.165 0.0583 0.0828 0.330
## 1146 5 1 0.132 0.0552 0.0584 0.300
## 1417 1 1 0.000 NaN NA NA
library(ggplot2)
ggplot(,aes(x=ekm$time,y=ekm$surv))+
geom_step()+
labs(x = "Tempo", y = "Sobrevivência")+
theme(panel.grid = element_blank())
ena=survfit(coxph(Surv(tempos,cens)~1,method="breslow"))
summary(ena)
## Call: survfit(formula = coxph(Surv(tempos, cens) ~ 1, method = "breslow"))
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 7 49 1 0.9798 0.0200 0.94138 1.000
## 34 48 1 0.9596 0.0280 0.90629 1.000
## 42 47 1 0.9394 0.0339 0.87523 1.000
## 63 46 1 0.9192 0.0387 0.84633 0.998
## 64 45 1 0.8990 0.0428 0.81885 0.987
## 83 43 1 0.8783 0.0466 0.79165 0.975
## 84 42 1 0.8577 0.0498 0.76533 0.961
## 91 41 1 0.8370 0.0528 0.73974 0.947
## 108 40 1 0.8163 0.0553 0.71475 0.932
## 112 39 1 0.7957 0.0577 0.69028 0.917
## 129 38 1 0.7750 0.0598 0.66628 0.901
## 133 37 2 0.7342 0.0632 0.62024 0.869
## 139 35 1 0.7135 0.0647 0.59734 0.852
## 140 34 2 0.6728 0.0671 0.55328 0.818
## 146 32 1 0.6521 0.0682 0.53126 0.800
## 149 31 1 0.6314 0.0691 0.50951 0.782
## 154 30 1 0.6107 0.0698 0.48804 0.764
## 160 29 2 0.5700 0.0709 0.44671 0.727
## 165 27 1 0.5493 0.0713 0.42593 0.708
## 173 26 1 0.5285 0.0715 0.40540 0.689
## 176 25 1 0.5078 0.0717 0.38511 0.670
## 218 23 1 0.4862 0.0718 0.36403 0.649
## 225 22 1 0.4646 0.0718 0.34322 0.629
## 241 21 1 0.4430 0.0716 0.32269 0.608
## 248 20 1 0.4214 0.0713 0.30244 0.587
## 273 19 1 0.3998 0.0708 0.28247 0.566
## 277 18 1 0.3782 0.0702 0.26279 0.544
## 297 16 1 0.3553 0.0696 0.24197 0.522
## 405 14 1 0.3308 0.0690 0.21979 0.498
## 417 13 1 0.3063 0.0681 0.19811 0.474
## 420 12 1 0.2818 0.0669 0.17695 0.449
## 440 11 1 0.2573 0.0654 0.15634 0.423
## 583 9 1 0.2302 0.0639 0.13367 0.397
## 594 8 1 0.2032 0.0618 0.11192 0.369
## 1101 7 1 0.1761 0.0592 0.09114 0.340
## 1146 5 1 0.1442 0.0564 0.06700 0.310
## 1417 1 1 0.0531 0.0570 0.00647 0.435
a =(218-176)
b =(0.497-0.501)
c =(0.5-0.501)
md = ((c*a)/b)+176
md
## [1] 186.5
t<- tempos[cens==1]
tj<-c(0,as.numeric(levels(as.factor(t))))
surv<-c(1,as.numeric(levels(as.factor(ekm$surv))))
surv<-sort(surv, decreasing=T)
k<-length(tj)-1
prod<-matrix(0,k,1)
for(j in 1:k){
prod[j]<-(tj[j+1]-tj[j])*surv[j]
}
tm<-sum(prod)
tm
## [1] 427.1155
time=tempos
cens=cens
TTTplot <- function(time,censur=NULL)
{
n <- length(time)
oo <- order(time)
if(is.null(censur)) censur = rep(1,n)
sorttime <- time[oo]
sortcens <- cens[oo]
ttt <- numeric(n)
ttt[1] <- n*sorttime[1]
for(i in 2:n)
{
ttt[i] <- ttt[i-1]+(n-i+1)*(sorttime[i]-sorttime[i-1])
}
ttt <- data.frame(time=sorttime,censur = sortcens,
rank=1:n,"cum total time"=ttt)
colnames(ttt) <- c("time","censur","rank","Cum.Total.Time")#,"i/n","TTT")
TTT <- ttt[ttt$censur==1,]
nfail <- dim(TTT)[1]
TTT <- cbind(TTT[,-3],TTT$Cum.Total.Time/TTT$Cum.Total.Time[nfail],
(1:nfail)/nfail)
colnames(TTT) <- c("time","censur","Cum.Total.Time","TTT","i/n")
W <- sum(TTT$"Cum.Total.Time"[1:(nfail-1)])/TTT$"Cum.Total.Time"[nfail]
Z <- (W-(nfail-1)/2)/sqrt((nfail-1)/12)
print(TTT)
cat("Barlow-Proschan's test\n")
cat(paste("W=",round(W,4)," k=",nfail,"\n",sep=""))
cat(paste("Z=",round(Z,4)," p.value=",round(pnorm(-abs(Z)),4),"\n",sep=""))
plot(TTT$"i/n",TTT$TTT,xlab="",ylab="TTT-plot",xlim=c(0,1),
ylim=c(0,1),pch=20, main="",sub=paste("nº. falhas = ",nfail," nº. censuras = ", n-nfail,sep=""),type="l",lwd=1)
lines(c(0,1),c(0,1),lty=2,lwd=1)
# return(list(TTT=TTT,W=W,Z=Z))
}
TTTplot(time,cens)
## time censur Cum.Total.Time TTT i/n
## 1 7 1 343 0.01952191 0.025
## 2 34 1 1639 0.09328401 0.050
## 3 42 1 2015 0.11468412 0.075
## 4 63 1 2981 0.16966420 0.100
## 5 64 1 3026 0.17222538 0.125
## 7 83 1 3853 0.21929425 0.150
## 8 84 1 3895 0.22168469 0.175
## 9 91 1 4182 0.23801935 0.200
## 10 108 1 4862 0.27672168 0.225
## 11 112 1 5018 0.28560046 0.250
## 12 129 1 5664 0.32236767 0.275
## 13 133 1 5812 0.33079112 0.300
## 14 133 1 5812 0.33079112 0.325
## 15 139 1 6022 0.34274331 0.350
## 16 140 1 6056 0.34467843 0.375
## 17 140 1 6056 0.34467843 0.400
## 18 146 1 6248 0.35560615 0.425
## 19 149 1 6341 0.36089926 0.450
## 20 154 1 6491 0.36943654 0.475
## 21 160 1 6665 0.37933978 0.500
## 22 160 1 6665 0.37933978 0.525
## 23 165 1 6800 0.38702334 0.550
## 24 173 1 7008 0.39886170 0.575
## 25 176 1 7083 0.40313034 0.600
## 27 218 1 8058 0.45862265 0.625
## 28 225 1 8212 0.46738759 0.650
## 29 241 1 8548 0.48651110 0.675
## 30 248 1 8688 0.49447923 0.700
## 31 273 1 9163 0.52151394 0.725
## 32 277 1 9235 0.52561184 0.750
## 34 297 1 9557 0.54393853 0.775
## 36 405 1 11091 0.63124644 0.800
## 37 417 1 11247 0.64012521 0.825
## 38 420 1 11283 0.64217416 0.850
## 39 440 1 11503 0.65469550 0.875
## 41 583 1 12873 0.73266932 0.900
## 42 594 1 12961 0.73767786 0.925
## 43 1101 1 16510 0.93966989 0.950
## 45 1146 1 16750 0.95332954 0.975
## 49 1417 1 17570 1.00000000 1.000
## Barlow-Proschan's test
## W=16.29 k=40
## Z=-1.7806 p.value=0.0375
Os dados representados na tabela abaixo representam o tempo (em dias) até a morte de pacientes com câncer de ovário tratados na Mayo Clinic Obtenha: o estimador de Kaplan-Meier (estratificado), de Nelson-Aalen, o tempo mediano e o tempo médio para cada grupo, compare as curvas utilizando os testes de Log-rank e Peto.
T1<-c(28,89,175,195,309,-377,-393,-421,-447,462,-709,-744,-770,-1106,-1206)
C1<-ifelse(T1<=0,0,1)
T1<-abs(T1)
T2<-c(34,88,137,199,280,291,-299,-300,309,351,358,369,369,370,375,382,392,-429,451,-1119)
C2<-ifelse(T2<=0,0,1)
T2<-abs(T2)
Data<-cbind(T2,C2) # tumor P
Data2<-cbind(T1,C1) # tumor G
ekm1<-survfit(Surv(Data[,1],Data[,2])~1)
ekm2<-survfit(Surv(Data2[,1],Data2[,2])~1)
ggplot()+
geom_step(data = data.frame(time = ekm1$time, surv = ekm1$surv, group = "Tumor pequeno"),
aes(x = time, y = surv, color = group))+
geom_step(data = data.frame(time = ekm2$time, surv = ekm2$surv, group = "Tumor grande"), aes(x = time, y = surv, color = group),size = 0.75) +
theme_minimal()+
labs(x = "Tempo", y = "Sobrevivência", color = "Grupos")
ENA1<-survfit(coxph(Surv(Data[,1],Data[,2])~1,method="breslow"))
ENA2<-survfit(coxph(Surv(Data2[,1],Data2[,2])~1,method="breslow"))
summary(ENA1)
## Call: survfit(formula = coxph(Surv(Data[, 1], Data[, 2]) ~ 1, method = "breslow"))
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 34 20 1 0.951 0.0476 0.8624 1.000
## 88 19 1 0.902 0.0655 0.7828 1.000
## 137 18 1 0.854 0.0780 0.7137 1.000
## 199 17 1 0.805 0.0875 0.6505 0.996
## 280 16 1 0.756 0.0948 0.5914 0.967
## 291 15 1 0.707 0.1005 0.5355 0.934
## 309 12 1 0.651 0.1072 0.4713 0.899
## 351 11 1 0.594 0.1118 0.4110 0.859
## 358 10 1 0.538 0.1145 0.3542 0.816
## 369 9 2 0.431 0.1140 0.2563 0.723
## 370 7 1 0.373 0.1123 0.2070 0.673
## 375 6 1 0.316 0.1086 0.1610 0.620
## 382 5 1 0.259 0.1029 0.1186 0.564
## 392 4 1 0.201 0.0947 0.0802 0.506
## 451 2 1 0.122 0.0838 0.0318 0.469
summary(ENA2)
## Call: survfit(formula = coxph(Surv(Data2[, 1], Data2[, 2]) ~ 1, method = "breslow"))
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 28 15 1 0.936 0.0624 0.821 1.000
## 89 14 1 0.871 0.0851 0.719 1.000
## 175 13 1 0.807 0.1003 0.632 1.000
## 195 12 1 0.742 0.1111 0.553 0.995
## 309 11 1 0.678 0.1187 0.481 0.955
## 462 6 1 0.574 0.1387 0.357 0.921
ggplot() +
geom_step(data = data.frame(time = ENA1$time, surv = ENA1$surv, group = "Tumor pequeno"), aes(x = time, y = surv, color = group),size = 0.75) +
geom_step(data = data.frame(time = ENA2$time, surv = ENA2$surv, group = "Tumor grande"), aes(x = time, y = surv, color = group),size = 0.75) +
labs(x = "Tempo", y = "Sobrevivência", color = "Grupos")+
theme_light()
a = (369-358)
b = (0.4083 - 0.5250)
c = (0.5-0.5250)
md = ((c*a)/b)+358
md
## [1] 360.3565
tempos2 = c(28,89,175,195,309,377,393,421,447,462,709,744,770,1106,1206,34,88,137,199,280,291,299,300,309,351,358,369,369,370,375,382,392,429,451,1119)
cens2= c(1,1, 1, 1, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0,
1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0)
grupos = c(rep(1,15),rep(2,20))
ekm3=survfit(Surv(tempos2,cens2)~grupos)
t<- tempos2[cens2==1]
tj<-c(0,as.numeric(levels(as.factor(t))))
surv<-c(1,as.numeric(levels(as.factor(ekm3$surv))))
surv<-sort(surv, decreasing=T)
k<-length(tj)-1
prod<-matrix(0,k,1)
for(j in 1:k){
prod[j]<-(tj[j+1]-tj[j])*surv[j]
}
tm<-sum(prod)
tm
## [1] 317.7694
survdiff(Surv(tempos2,cens2)~grupos,rho=0)
## Call:
## survdiff(formula = Surv(tempos2, cens2) ~ grupos, rho = 0)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## grupos=1 15 6 11.3 2.51 5.57
## grupos=2 20 16 10.7 2.67 5.57
##
## Chisq= 5.6 on 1 degrees of freedom, p= 0.02
survdiff(Surv(tempos2,cens2)~grupos,rho=1)
## Call:
## survdiff(formula = Surv(tempos2, cens2) ~ grupos, rho = 1)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## grupos=1 15 4.66 7.40 1.013 2.74
## grupos=2 20 10.60 7.87 0.953 2.74
##
## Chisq= 2.7 on 1 degrees of freedom, p= 0.1
Dados de tempo até a ocorrência de cio em ovelhas submetidas a três tratamentos.
TratA<-c(57,-60,-60,-60,48,40,49,45,25,-60,51)
TratB<-c(32,-60,-60,22,-60,54,-60,-60,-60,51)
TratC<-c(28,32,-60,41,-60,54,59,40,32,-60,-60)
Censura_TratA<-ifelse(TratA<=0,0,1)
Censura_TratB<-ifelse(TratB<=0,0,1)
Censura_TratC<-ifelse(TratC<=0,0,1)
TratA<-abs(TratA)
TratB<-abs(TratB)
TratC<-abs(TratC)
Data_TratA<-cbind(TratA,Censura_TratA)
Data_TratB<-cbind(TratB,Censura_TratB)
Data_TratC<-cbind(TratC,Censura_TratC)
EKM1<-survfit(Surv(Data_TratA[,1],Data_TratA[,2])~1)
EKM2<-survfit(Surv(Data_TratB[,1],Data_TratB[,2])~1)
EKM3<-survfit(Surv(Data_TratC[,1],Data_TratC[,2])~1)
summary(EKM1)
## Call: survfit(formula = Surv(Data_TratA[, 1], Data_TratA[, 2]) ~ 1)
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 25 11 1 0.909 0.0867 0.754 1.000
## 40 10 1 0.818 0.1163 0.619 1.000
## 45 9 1 0.727 0.1343 0.506 1.000
## 48 8 1 0.636 0.1450 0.407 0.995
## 49 7 1 0.545 0.1501 0.318 0.936
## 51 6 1 0.455 0.1501 0.238 0.868
## 57 5 1 0.364 0.1450 0.166 0.795
summary(EKM2)
## Call: survfit(formula = Surv(Data_TratB[, 1], Data_TratB[, 2]) ~ 1)
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 22 10 1 0.9 0.0949 0.732 1.000
## 32 9 1 0.8 0.1265 0.587 1.000
## 51 8 1 0.7 0.1449 0.467 1.000
## 54 7 1 0.6 0.1549 0.362 0.995
summary(EKM3)
## Call: survfit(formula = Surv(Data_TratC[, 1], Data_TratC[, 2]) ~ 1)
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 28 11 1 0.909 0.0867 0.754 1.000
## 32 10 2 0.727 0.1343 0.506 1.000
## 40 8 1 0.636 0.1450 0.407 0.995
## 41 7 1 0.545 0.1501 0.318 0.936
## 54 6 1 0.455 0.1501 0.238 0.868
## 59 5 1 0.364 0.1450 0.166 0.795
ggplot() +
geom_step(data = data.frame(time = EKM1$time, surv = EKM1$surv, group = "Tratamento A"), aes(x = time, y = surv, color = group),size = 0.75) +
geom_step(data = data.frame(time = EKM2$time, surv = EKM2$surv, group = "Tratamento B"), aes(x = time, y = surv, color = group),size = 0.75) +
geom_step(data = data.frame(time = EKM3$time, surv = EKM3$surv, group = "Tratamento C"), aes(x = time, y = surv, color = group),size = 0.75) +
labs(x = "Tempo", y = "Sobrevivência", color = "Grupos")+
theme_light()
EnaA<-survfit(coxph(Surv(Data_TratA[,1],Data_TratA[,2])~1,method="breslow"))
EnaB<-survfit(coxph(Surv(Data_TratB[,1],Data_TratB[,2])~1,method="breslow"))
EnaC<-survfit(coxph(Surv(Data_TratC[,1],Data_TratC[,2])~1,method="breslow"))
summary(EnaA)
## Call: survfit(formula = coxph(Surv(Data_TratA[, 1], Data_TratA[, 2]) ~
## 1, method = "breslow"))
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 25 11 1 0.913 0.083 0.764 1.000
## 40 10 1 0.826 0.112 0.634 1.000
## 45 9 1 0.739 0.129 0.525 1.000
## 48 8 1 0.652 0.140 0.428 0.994
## 49 7 1 0.566 0.146 0.341 0.938
## 51 6 1 0.479 0.147 0.262 0.874
## 57 5 1 0.392 0.144 0.191 0.804
summary(EnaA)
## Call: survfit(formula = coxph(Surv(Data_TratA[, 1], Data_TratA[, 2]) ~
## 1, method = "breslow"))
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 25 11 1 0.913 0.083 0.764 1.000
## 40 10 1 0.826 0.112 0.634 1.000
## 45 9 1 0.739 0.129 0.525 1.000
## 48 8 1 0.652 0.140 0.428 0.994
## 49 7 1 0.566 0.146 0.341 0.938
## 51 6 1 0.479 0.147 0.262 0.874
## 57 5 1 0.392 0.144 0.191 0.804
summary(EnaB)
## Call: survfit(formula = coxph(Surv(Data_TratB[, 1], Data_TratB[, 2]) ~
## 1, method = "breslow"))
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 22 10 1 0.905 0.0905 0.744 1.000
## 32 9 1 0.810 0.1210 0.604 1.000
## 51 8 1 0.715 0.1392 0.488 1.000
## 54 7 1 0.619 0.1497 0.386 0.995
summary(EnaC)
## Call: survfit(formula = coxph(Surv(Data_TratC[, 1], Data_TratC[, 2]) ~
## 1, method = "breslow"))
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 28 11 1 0.913 0.083 0.764 1.000
## 32 10 2 0.748 0.126 0.538 1.000
## 40 8 1 0.660 0.138 0.438 0.995
## 41 7 1 0.572 0.145 0.348 0.940
## 54 6 1 0.484 0.147 0.267 0.877
## 59 5 1 0.396 0.144 0.194 0.808
ggplot() +
geom_step(data = data.frame(time = EnaA$time, surv = EnaA$surv, group = "Tratamento A"), aes(x = time, y = surv, color = group),size = 0.75) +
geom_step(data = data.frame(time = EnaB$time, surv = EnaB$surv, group = "Tratamento B"), aes(x = time, y = surv, color = group),size = 0.75) +
geom_step(data = data.frame(time = EnaC$time, surv = EnaC$surv, group = "Tratamento C"), aes(x = time, y = surv, color = group),size = 0.75) +
labs(x = "Tempo", y = "Sobrevivência", color = "Grupos")+
theme_light()
tempos3 = c(57, 60, 60, 60, 48, 40, 49, 45, 25, 60, 51,
32, 60, 60, 22, 60, 54, 60, 60, 60, 51,
28, 32, 60, 41, 60, 54, 59, 40, 32, 60, 60)
cens3 = c(1, 0, 0, 0, 1, 1, 1, 1, 1, 0, 1,
1, 0, 0, 1, 0, 1, 0, 0, 0, 1,
1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 0)
grupos3 = c(1,1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2,3,3,3,3,3,3,3,3,3,3,3)
t<- tempos3[cens3==1]
tj<-c(0,as.numeric(levels(as.factor(t))))
surv<-c(1,as.numeric(levels(as.factor(ekm3$surv))))
surv<-sort(surv, decreasing=T)
k<-length(tj)-1
prod<-matrix(0,k,1)
for(j in 1:k){
prod[j]<-(tj[j+1]-tj[j])*surv[j]
}
tm<-sum(prod)
tm
## [1] 51.70833
survdiff(Surv(tempos3,cens3)~grupos3,rho=0)
## Call:
## survdiff(formula = Surv(tempos3, cens3) ~ grupos3, rho = 0)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## grupos3=1 11 7 6.03 0.155 0.241
## grupos3=2 10 4 6.23 0.797 1.263
## grupos3=3 11 7 5.74 0.277 0.421
##
## Chisq= 1.3 on 2 degrees of freedom, p= 0.5
survdiff(Surv(tempos3,cens3)~grupos3,rho=1)
## Call:
## survdiff(formula = Surv(tempos3, cens3) ~ grupos3, rho = 1)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## grupos3=1 11 4.97 4.56 0.0362 0.0723
## grupos3=2 10 3.09 4.53 0.4560 0.9241
## grupos3=3 11 5.34 4.31 0.2466 0.4792
##
## Chisq= 1 on 2 degrees of freedom, p= 0.6