Exo 2 TP Généralités
Olivier Bouaziz
- Les variables d’intérêt sont notées : \(\tilde T_i\)
- Les observations sont : \[\begin{cases} T_i=\min(\tilde T_i,C_i)\\ \Delta_i=I(\tilde T_i\leq C_i). \end{cases}\]
Lois uniformes
On commence par simuler des lois uniformes pour \(\tilde T_i\) et \(C_i\) :
rm(list=ls())
n<-200
trueT<-runif(n,0,5)
cens<-runif(n,0,5)
obs=pmin(trueT,cens)
delta=(trueT<=cens)
mean(delta)## [1] 0.48\(48\%\) d’observations non censurées. Pour voir les données :
head(data.frame(trueT,cens,obs,delta),n=20L)## trueT cens obs delta
## 1 3.4172471 2.1480955 2.1480955 FALSE
## 2 0.5968966 2.6010677 0.5968966 TRUE
## 3 0.6981828 4.2585836 0.6981828 TRUE
## 4 1.0497051 2.6547627 1.0497051 TRUE
## 5 3.0749693 0.4163665 0.4163665 FALSE
## 6 2.3980746 1.2809705 1.2809705 FALSE
## 7 4.1889200 3.9649998 3.9649998 FALSE
## 8 0.9554780 1.4143914 0.9554780 TRUE
## 9 3.4503026 0.4774126 0.4774126 FALSE
## 10 1.3830175 1.2130879 1.2130879 FALSE
## 11 2.0105168 2.0634945 2.0105168 TRUE
## 12 4.0291471 3.2490187 3.2490187 FALSE
## 13 2.5420637 1.9127640 1.9127640 FALSE
## 14 0.8466576 3.2709767 0.8466576 TRUE
## 15 1.6519884 1.7803326 1.6519884 TRUE
## 16 0.8176098 4.9745517 0.8176098 TRUE
## 17 3.3123046 0.4846096 0.4846096 FALSE
## 18 3.5900369 2.2258552 2.2258552 FALSE
## 19 3.8219451 4.1789796 3.8219451 TRUE
## 20 2.5653364 0.6211636 0.6211636 FALSEOn trace les courbes de survie avec les différentes estimations :
plot(c(0,sort(obs)),c(n:0)/n,type="s",col="green",xlab="Temps passé dans l'étude",ylab=
"Survie")
absi<-seq(0,6,by=0.01)
lines(absi,1-punif(absi,0,5),type="l")
lines(c(0,sort(trueT)),c(n:0)/n,type="s",col="red")
segments(max(trueT),0,5,0,col="red")
d<-sum(delta)
lines(c(0,sort(obs[delta==1])),c(d:0)/n,type="s",col="blue")
segments(max(obs[delta==1]),0,5,0,col="blue")
lines(c(0,sort(obs[delta==1])),c(d:0)/d,type="s",col="cyan")
segments(max(obs[delta==1]),0,5,0,col="cyan")
legend("topright",lty=c(1,1,1,1,1),col=c(1,2,3,4,5),c("Vraie fonction de survie",
"Estimation sans biais","Estimation par les observations",
"Estim. par les obs. non-cens. (norm. par n)","Estim. par les obs. non-cens.
(correct. norm.)"),cex=0.5,inset=0.01)
- En noir la vraie fonction de survie d’une loi uniforme de paramètres \(0\) et \(5\).
- En rouge l’estimation de la courbe de survie à partir des vrais temps d’intérêt, ils ne sont pas tous observés en pratique et cet estimateur ne peut pas être utilisé en pratique.
- En vert l’estimation de la courbe de survie par les observations :
\[\hat S(t)=\frac{1}{n}\sum_{i=1}^nI(T_i>t)\]
- En bleu l’estimation de la courbe de survie par les observations non-censurées (renormalisé par \(n\)) :
\[\hat S(t)=\frac{1}{n}\sum_{i=1}^nI(T_i>t,\Delta_i=1)\]
- En bleu clair l’estimation de la courbe de survie par les observations non-censurées (renormalisé par le nombre d’observations non-censurées) :
\[\hat S(t)=\frac{1}{\sum_{i=1}^n \Delta_i}\sum_{i=1}^nI(T_i>t,\Delta_i=1)\]
On voit bien que nos estimateurs n’estiment pas bien \(S(t)\) !! Les données observées ou bien les données des observations non-censurées ont tendance à être plus petites que les temps d’intérêt. On le voit bien regardant les boîtes à moustaches (de gauche à droite) des temps d’intérêt, des temps observés et des temps observés non-censurés :
boxplot(trueT,obs,obs[delta==1])
On estime la moyenne par les différents estimateurs. La vraie moyenne est égale à 2.5 pour une loi uniforme de paramètre 0 et 5.
mean(obs)## [1] 1.648147mean(trueT) ## [1] 2.441301mean(obs[delta==1])## [1] 1.5465sum(obs[delta==1])/n## [1] 0.7423199Lois exponentielles
On simule cette fois-ci des lois exponentielles pour \(\tilde T_i\) et \(C_i\) :
n<-200
trueT<-rexp(n,2)
cens<-rexp(n,3)
obs=pmin(trueT,cens)
delta=(trueT<=cens)
mean(delta)## [1] 0.425plot(c(0,sort(obs)),c(n:0)/n,type="s",col="green",xlab="Temps passé dans l'étude",ylab="Survie")
absi<-seq(0,6,by=0.01)
lines(absi,1-pexp(absi,2),type="l")
lines(c(0,sort(trueT)),c(n:0)/n,type="s",col="red")
segments(max(trueT),0,5,0,col="red")
d<-sum(delta)
lines(c(0,sort(obs[delta==1])),c(d:0)/n,type="s",col="blue")
segments(max(obs[delta==1]),0,5,0,col="blue")
lines(c(0,sort(obs[delta==1])),c(d:0)/d,type="s",col="cyan")
segments(max(obs[delta==1]),0,5,0,col="cyan")
legend("topright",lty=c(1,1,1,1,1),col=c(1,2,3,4,5),c("Vraie fonction de survie",
"Estimation sans biais","Estimation par les observations",
"Estim. par les obs. non-cens. (norm. par n)","Estim. par les obs. non-cens. (correct. norm.)")
,cex=0.5,inset=0.01)
La vraie moyenne pour une loi expoentielle de paramètre \(2\) est égale à \(0.5\).
mean(obs)## [1] 0.1987779mean(trueT)## [1] 0.4720051mean(obs[delta==1])## [1] 0.2205255sum(obs[delta==1])/n## [1] 0.09372336