Prova 1 de análise de sobrevivência
Arthur Henrique Elias
23/03/2023
library (asaur)
library(survival)
dados=asaur::pharmacoSmoking
attach(dados)
cens_smoking<-dados[,3]
tempo_smoking=ttr
Surv(tempo_smoking,cens_smoking)## [1] 182+ 14 5 16 0 182+ 14 77 2 0 12 182+ 21 3 170
## [16] 25 4 182+ 140 63 15 140 110 182+ 0 182+ 15 182+ 4 56
## [31] 2 80 182+ 56 0 14 14 28 182+ 6 182+ 14 15 182+ 75
## [46] 30 4 56 182+ 182+ 5 8 140 20 63 30 8 50 14 0
## [61] 84 0 105 182+ 182+ 182+ 7 182+ 0 8 1 182+ 12 182+ 49
## [76] 182+ 182+ 2 182+ 56 182+ 0 28 155 2 0 0 1 140 1
## [91] 28 1 182+ 77 56 182+ 182+ 182+ 182+ 182+ 21 60 0 182+ 65
## [106] 182+ 182+ 182+ 182+ 2 40 100 1 45 14 30 42 2 182+ 60
## [121] 10 0 170 15 182+- A
ekmE<- survfit(Surv(tempo_smoking,cens_smoking)~1)
summary(ekmE)## Call: survfit(formula = Surv(tempo_smoking, cens_smoking) ~ 1)
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 0 125 12 0.904 0.0263 0.854 0.957
## 1 113 5 0.864 0.0307 0.806 0.926
## 2 108 6 0.816 0.0347 0.751 0.887
## 3 102 1 0.808 0.0352 0.742 0.880
## 4 101 3 0.784 0.0368 0.715 0.860
## 5 98 2 0.768 0.0378 0.697 0.846
## 6 96 1 0.760 0.0382 0.689 0.839
## 7 95 1 0.752 0.0386 0.680 0.832
## 8 94 3 0.728 0.0398 0.654 0.810
## 10 91 1 0.720 0.0402 0.645 0.803
## 12 90 2 0.704 0.0408 0.628 0.789
## 14 88 7 0.648 0.0427 0.569 0.737
## 15 81 4 0.616 0.0435 0.536 0.707
## 16 77 1 0.608 0.0437 0.528 0.700
## 20 76 1 0.600 0.0438 0.520 0.692
## 21 75 2 0.584 0.0441 0.504 0.677
## 25 73 1 0.576 0.0442 0.496 0.669
## 28 72 3 0.552 0.0445 0.471 0.646
## 30 69 3 0.528 0.0447 0.447 0.623
## 40 66 1 0.520 0.0447 0.439 0.615
## 42 65 1 0.512 0.0447 0.431 0.608
## 45 64 1 0.504 0.0447 0.424 0.600
## 49 63 1 0.496 0.0447 0.416 0.592
## 50 62 1 0.488 0.0447 0.408 0.584
## 56 61 5 0.448 0.0445 0.369 0.544
## 60 56 2 0.432 0.0443 0.353 0.528
## 63 54 2 0.416 0.0441 0.338 0.512
## 65 52 1 0.408 0.0440 0.330 0.504
## 75 51 1 0.400 0.0438 0.323 0.496
## 77 50 2 0.384 0.0435 0.308 0.479
## 80 48 1 0.376 0.0433 0.300 0.471
## 84 47 1 0.368 0.0431 0.292 0.463
## 100 46 1 0.360 0.0429 0.285 0.455
## 105 45 1 0.352 0.0427 0.277 0.447
## 110 44 1 0.344 0.0425 0.270 0.438
## 140 43 4 0.312 0.0414 0.240 0.405
## 155 39 1 0.304 0.0411 0.233 0.396
## 170 38 2 0.288 0.0405 0.219 0.379plot(ekmE, xlab="Tempo (semanas)",ylab="S(t) estimada", main = "Curva de Sobrevivência Estimada para dados de Smoking")
Nesse exercicio, sobreviver significa não voltar a fumar.
- B
ekm2<- survfit(Surv(tempo_smoking,cens_smoking)~1)
t<- tempo_smoking[cens_smoking==1]
tj<-c(0,as.numeric(levels(as.factor(t))))
surv<-c(1,as.numeric(levels(as.factor(ekm2$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] 73.984tempo medio de 147.9398 e a mediana é 49.
- C Para questão 3 letra C, usaremos as variaveis grp e levelsmoking para comparamos.
ekm4E<- survfit(Surv(tempo_smoking,cens_smoking)~levelSmoking)
plot(ekm4E,xlab="Tempo (dias)",ylab="S(t) estimada",col=c(1,2))
legend("topright",col=c(2,1),c("leve","forte"),lwd=1, bty="n")
As curvas se tocam em mais de um ponto, logo o teste de Peto deve ser o mais indicado. Vamos realizar os dois apenas por fins interpretativos.
survdiff(Surv(tempo_smoking,cens_smoking)~levelSmoking,rho=0)## Call:
## survdiff(formula = Surv(tempo_smoking, cens_smoking) ~ levelSmoking,
## rho = 0)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## levelSmoking=heavy 89 62 62.7 0.00842 0.0297
## levelSmoking=light 36 27 26.3 0.02010 0.0297
##
## Chisq= 0 on 1 degrees of freedom, p= 0.9survdiff(Surv(tempo_smoking,cens_smoking)~levelSmoking,rho=1)## Call:
## survdiff(formula = Surv(tempo_smoking, cens_smoking) ~ levelSmoking,
## rho = 1)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## levelSmoking=heavy 89 42.1 41.4 0.0112 0.054
## levelSmoking=light 36 16.8 17.5 0.0265 0.054
##
## Chisq= 0.1 on 1 degrees of freedom, p= 0.8A um nivel de significancia de 5%, não rejeito a hipotese de que existe diferença nos grupos.
ekm6E<- survfit(Surv(tempo_smoking,cens_smoking)~grp)
plot(ekm6E,xlab="Tempo (dias)",ylab="S(t) estimada",col=c(1,2))
legend("topright",col=c(2,1),c("combination","patchOnly"),lwd=1, bty="n")
As curvas não tocam em mais de um ponto, logo podemos usar os dois testes.
survdiff(Surv(tempo_smoking,cens_smoking)~grp,rho=0)## Call:
## survdiff(formula = Surv(tempo_smoking, cens_smoking) ~ grp, rho = 0)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## grp=combination 61 37 49.9 3.36 8.03
## grp=patchOnly 64 52 39.1 4.29 8.03
##
## Chisq= 8 on 1 degrees of freedom, p= 0.005survdiff(Surv(tempo_smoking,cens_smoking)~grp,rho=1)## Call:
## survdiff(formula = Surv(tempo_smoking, cens_smoking) ~ grp, rho = 1)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## grp=combination 61 23.1 32.1 2.53 8.01
## grp=patchOnly 64 35.8 26.8 3.04 8.01
##
## Chisq= 8 on 1 degrees of freedom, p= 0.005A um nivel de significancia de 5%, rejeito a hipotese de que existe diferença nos grupos.