knitr::opts_chunk$set(echo = TRUE)

PHASE II Gastic cancer / OXAliplatin result trt with surgery and metastatic patients

C8H14N2O4Pt

oxyPt geometry

library(asaur)
head(gastricXelox,10)
##    timeWeeks delta
## 1          4     1
## 2          8     1
## 3          8     1
## 4          8     1
## 5          9     1
## 6         11     1
## 7         12     1
## 8         13     1
## 9         16     1
## 10        16     1
table(gastricXelox$delta)#censoring data 1death of Progression 0 censored info xis
## 
##  0  1 
## 16 32
range(gastricXelox$timeWeeks)#length of followup time rounded nearset week
## [1]   4 253
plot(gastricXelox,pch=20)#a lot of deaths occurs soon in the followup cohort
abline(v=max(gastricXelox$timeWeeks[gastricXelox$delta==1]))

#note that censoring clusterr occurs before that timeframe...!!!Suspect Usually expect uniformly distributed censoring if not informative one!

Kaplan Meier estimator

library(survival)
xeloxframe=Surv(gastricXelox$timeWeeks,gastricXelox$delta)
survfit(Surv(gastricXelox$timeWeeks,gastricXelox$delta)~1)#an surv object KM summary
## Call: survfit(formula = Surv(gastricXelox$timeWeeks, gastricXelox$delta) ~ 
##     1)
## 
##       n events median 0.95LCL 0.95UCL
## [1,] 48     32   44.5      28      76
plot(survfit(Surv(gastricXelox$timeWeeks,gastricXelox$delta)~1),xlab="weeks",ylab="prob Survival",main="Xelox Phase II Cohort- Kaplan Meier estim.",sub="conftype=plain",conf.type="plain" )
abline(v=44.5,lty=2,col=2)

plot(survfit(Surv(gastricXelox$timeWeeks,gastricXelox$delta)~1),xlab="weeks",ylab="prob Survival",main="Xelox Phase II Cohort- Kaplan Meier estim.",mark.time = TRUE,conf.type = "log-log",sub="conftype=log-log")
abline(v=44.5,lty=2,col=2)

HAZARD RATE

library(muhaz)
## Warning: le package 'muhaz' a été compilé avec la version R 4.5.1
muhaz(gastricXelox$timeWeeks,gastricXelox$delta,bw.grid=5,b.cor="none")
## $pin
## $pin$times
##  [1]   4   8   8   8   9  11  12  13  16  16  17  17  19  21  24  24  25  28  28
## [20]  30  37  37  42  43  43  46  48  50  51  53  54  57  59  59  60  64  66  76
## [39]  78  91 119 120 139 176 187 203 217 253
## 
## $pin$delta
##  [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 0 0 0 1 0 0 0 1 1 1 1 1
## [39] 1 0 0 0 0 0 0 0 0 0
## 
## $pin$nobs
## [1] 48
## 
## $pin$min.time
## [1] 0
## 
## $pin$max.time
## [1] 78
## 
## $pin$n.min.grid
## [1] 51
## 
## $pin$min.grid
##  [1]  0.00  1.56  3.12  4.68  6.24  7.80  9.36 10.92 12.48 14.04 15.60 17.16
## [13] 18.72 20.28 21.84 23.40 24.96 26.52 28.08 29.64 31.20 32.76 34.32 35.88
## [25] 37.44 39.00 40.56 42.12 43.68 45.24 46.80 48.36 49.92 51.48 53.04 54.60
## [37] 56.16 57.72 59.28 60.84 62.40 63.96 65.52 67.08 68.64 70.20 71.76 73.32
## [49] 74.88 76.44 78.00
## 
## $pin$n.est.grid
## [1] 101
## 
## $pin$bw.pilot
## [1] 4.875
## 
## $pin$bw.smooth
## [1] 24.375
## 
## $pin$method
## [1] 2
## 
## $pin$b.cor
## [1] 0
## 
## $pin$kernel.type
## [1] 1
## 
## 
## $est.grid
##   [1]  0.00  0.78  1.56  2.34  3.12  3.90  4.68  5.46  6.24  7.02  7.80  8.58
##  [13]  9.36 10.14 10.92 11.70 12.48 13.26 14.04 14.82 15.60 16.38 17.16 17.94
##  [25] 18.72 19.50 20.28 21.06 21.84 22.62 23.40 24.18 24.96 25.74 26.52 27.30
##  [37] 28.08 28.86 29.64 30.42 31.20 31.98 32.76 33.54 34.32 35.10 35.88 36.66
##  [49] 37.44 38.22 39.00 39.78 40.56 41.34 42.12 42.90 43.68 44.46 45.24 46.02
##  [61] 46.80 47.58 48.36 49.14 49.92 50.70 51.48 52.26 53.04 53.82 54.60 55.38
##  [73] 56.16 56.94 57.72 58.50 59.28 60.06 60.84 61.62 62.40 63.18 63.96 64.74
##  [85] 65.52 66.30 67.08 67.86 68.64 69.42 70.20 70.98 71.76 72.54 73.32 74.10
##  [97] 74.88 75.66 76.44 77.22 78.00
## 
## $haz.est
##   [1] 0.001125000 0.001828950 0.002380800 0.002780550 0.003492277 0.006329543
##   [7] 0.009402638 0.011819144 0.013768180 0.015575804 0.017413726 0.018913268
##  [13] 0.019861547 0.020148877 0.019272307 0.019209877 0.019321380 0.019542967
##  [19] 0.020368525 0.021969452 0.022086732 0.021897419 0.021666834 0.020918422
##  [25] 0.019950905 0.019601977 0.019188122 0.018107369 0.017591839 0.017718736
##  [31] 0.018785247 0.020155090 0.021113353 0.022137775 0.022579737 0.022061088
##  [37] 0.020147702 0.016839579 0.014576522 0.012669850 0.010683257 0.007966042
##  [43] 0.007583540 0.008266669 0.009088751 0.009336905 0.010365143 0.010862238
##  [49] 0.011798898 0.013251292 0.015018974 0.015682682 0.015242415 0.014554973
##  [55] 0.014174556 0.015594690 0.016124566 0.015764185 0.014513545 0.012372646
##  [61] 0.009341490 0.006836168 0.006164716 0.007139270 0.007412145 0.006983341
##  [67] 0.007165137 0.007721811 0.007894232 0.007682400 0.009342316 0.012430905
##  [73] 0.015909892 0.017996422 0.018690493 0.019650000 0.021716877 0.024637026
##  [79] 0.025987369 0.028675706 0.029998637 0.029143363 0.026109883 0.024077237
##  [85] 0.022856923 0.021551923 0.019076923 0.015431923 0.010616923 0.006651800
##  [91] 0.003680000 0.000099800 0.003830400 0.007106400 0.011577251 0.017541273
##  [97] 0.022111505 0.025287949 0.027070604 0.027459469 0.026454545
## 
## $bw.loc
##  [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [39] 0 0 0 0 0 0 0 0 0 0 0 0 0
## 
## $bw.loc.sm
##   [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
##  [38] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
##  [75] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## 
## $msemin
##  [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [39] 0 0 0 0 0 0 0 0 0 0 0 0 0
## 
## $bias.min
##  [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [39] 0 0 0 0 0 0 0 0 0 0 0 0 0
## 
## $var.min
##  [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [39] 0 0 0 0 0 0 0 0 0 0 0 0 0
## 
## $imse.opt
## [1] 0
## 
## attr(,"class")
## [1] "muhaz"
253/78#kernel range
## [1] 3.24359
plot(muhaz(gastricXelox$timeWeeks,gastricXelox$delta,bw.grid=4,b.cor="none"))

4*78
## [1] 312

It is not clear how the bin grid timeframe is construct What is 78 is weeks or 78 time the bin grid ? TBC soon

from book and R::

data(cancer, package="survival")
attach(ovarian)
summary(ovarian)
##      futime           fustat            age           resid.ds    
##  Min.   :  59.0   Min.   :0.0000   Min.   :38.89   Min.   :1.000  
##  1st Qu.: 368.0   1st Qu.:0.0000   1st Qu.:50.17   1st Qu.:1.000  
##  Median : 476.0   Median :0.0000   Median :56.85   Median :2.000  
##  Mean   : 599.5   Mean   :0.4615   Mean   :56.17   Mean   :1.577  
##  3rd Qu.: 794.8   3rd Qu.:1.0000   3rd Qu.:62.38   3rd Qu.:2.000  
##  Max.   :1227.0   Max.   :1.0000   Max.   :74.50   Max.   :2.000  
##        rx         ecog.ps     
##  Min.   :1.0   Min.   :1.000  
##  1st Qu.:1.0   1st Qu.:1.000  
##  Median :1.5   Median :1.000  
##  Mean   :1.5   Mean   :1.462  
##  3rd Qu.:2.0   3rd Qu.:2.000  
##  Max.   :2.0   Max.   :2.000
range(ovarian$futime)
## [1]   59 1227
plot(ovarian$fustat~ovarian$futime)

#last death occurs at 600 times
fit1 <- muhaz(futime, fustat)
plot(fit1)

#Kernel band with end at around approx last info of death After hazard ahas no reason to be calculated as either censored or Survived anyway
summary(fit1)
## 
## Number of Observations .......... 26
## Censored Observations ........... 14
## Method used ..................... Local
## Boundary Correction Type ........ Left and Right
## Kernel type ..................... Epanechnikov
## Minimum Time .................... 0
## Maximum Time .................... 744
## Number of minimization points ... 51
## Number of estimation points ..... 101
## Pilot Bandwidth ................. 56.58
## Smoothing Bandwidth ............. 282.89
## Minimum IMSE .................... 1e+30
# to compute a globally optimal estimate
fit2 <- muhaz(futime, fustat, bw.method="g")
plot(fit2)

fit2bis <- muhaz(futime, fustat, bw.method="local")
plot(fit2bis)

# to compute an estimate with global bandwidth set to 5
fit3 <- muhaz(futime, fustat, bw.method="g", bw.grid=5)
plot(fit3)

#See that hazard estimate when dt is close to null overstimate end time points due to 1/ censoring 2/decreaswe in population at rate and start point too
fit4 <- muhaz(futime, fustat, bw.method="g", bw.grid=2)
plot(fit4)##BW=2 seems more appropriate see next plot

plot(ovarian$fustat~ovarian$futime)

RAW h(t)

From KM calculation is death given alive until t / persons at risk in Delta Time

K=survfit(Surv(gastricXelox$timeWeeks,gastricXelox$delta)~1)#an surv object KM summary
K
## Call: survfit(formula = Surv(gastricXelox$timeWeeks, gastricXelox$delta) ~ 
##     1)
## 
##       n events median 0.95LCL 0.95UCL
## [1,] 48     32   44.5      28      76
K$n.risk
##  [1] 48 47 44 43 42 41 40 38 36 35 34 32 31 29 28 26 25 23 22 21 20 19 18 17 16
## [26] 14 13 12 11 10  9  8  7  6  5  4  3  2  1
plot(K$time,K$n.risk,type="h",pch=20)
points(K$time,K$n.risk,pch=20)
v=K$n.censor[K$n.censor==0]
points(v,pch="+",col=3)

K$n.event
##  [1] 1 3 1 1 1 1 2 2 1 1 2 1 2 1 2 1 1 1 0 0 0 1 0 0 1 1 1 1 1 1 0 0 0 0 0 0 0 0
## [39] 0
K$n
## [1] 48
gastricXelox$delta/K$n.risk
## Warning in gastricXelox$delta/K$n.risk: la taille d'un objet plus long n'est
## pas multiple de la taille d'un objet plus court
##  [1] 0.02083333 0.02127660 0.02272727 0.02325581 0.02380952 0.02439024
##  [7] 0.02500000 0.02631579 0.02777778 0.02857143 0.02941176 0.03125000
## [13] 0.03225806 0.03448276 0.03571429 0.03846154 0.04000000 0.04347826
## [19] 0.04545455 0.04761905 0.05000000 0.05263158 0.05555556 0.05882353
## [25] 0.00000000 0.07142857 0.00000000 0.00000000 0.00000000 0.10000000
## [31] 0.00000000 0.00000000 0.00000000 0.16666667 0.20000000 0.25000000
## [37] 0.33333333 0.50000000 1.00000000 0.00000000 0.00000000 0.00000000
## [43] 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000
dif=diff(gastricXelox$timeWeeks)
gastricXelox$timeWeeks
##  [1]   4   8   8   8   9  11  12  13  16  16  17  17  19  21  24  24  25  28  28
## [20]  30  37  37  42  43  43  46  48  50  51  53  54  57  59  59  60  64  66  76
## [39]  78  91 119 120 139 176 187 203 217 253
deltat=gastricXelox$delta/K$n.risk
## Warning in gastricXelox$delta/K$n.risk: la taille d'un objet plus long n'est
## pas multiple de la taille d'un objet plus court
dd=deltat[1:47]
gg=(gastricXelox$delta/K$n.risk)[1:47]
## Warning in gastricXelox$delta/K$n.risk: la taille d'un objet plus long n'est
## pas multiple de la taille d'un objet plus court
plot(gg/dif,type="s")#weel we are not far from calculation of Muhaz fx

pe1=pehaz(gastricXelox$timeWeeks,gastricXelox$delta,max.time = 60,width=1)
## 
## max.time= 60
## width= 1
## nbins= 60
pe2=pehaz(gastricXelox$timeWeeks,gastricXelox$delta,max.time = 60,width=10)
## 
## max.time= 60
## width= 10
## nbins= 6
pe1=pehaz(gastricXelox$timeWeeks,gastricXelox$delta,max.time = 60,width=1)
## 
## max.time= 60
## width= 1
## nbins= 60
pe3=pehaz(gastricXelox$timeWeeks,gastricXelox$delta,max.time = 60,width=4)
## 
## max.time= 60
## width= 4
## nbins= 15
plot(pe1)

plot(pe2)

plot(pe3)

library(lattice)
dotplot(dif,xlab="Diff Time intreval Delta T")

Estimating hazard is a very tricky part because is: 1/ it is a (dt) limit when time bin come close to zero therefore some of magnitude change can occur 2/ When constructing bins they is no grantee that the distribution reflect the empirical one (time delay) 3/ Assuming no constant time bins produces different hazard scaling improper to be interpreted 4/ Looking for the best parameters is trial and error and get it according what the data tells (delta==0 position 5/ Hazard is null when not smoothed when is no event (delta==1) and therefoer after last event no hazard is valuable for interpretation