Performing Pairwise test and test for survival trend

Question

I am trying to perform pairwise tests to determine if there is any difference in survival between pairs of groups.

the data used:

time_Untreated<- c(20, 21, 23, 24, 24, 26, 26, 27, 28, 30)
censor_Untreated<- c(rep(1,10), rep(0,0))
censor_Untreated

time_Radiated<- c(26,28, 29, 29, 30, 30, 31, 31, 32, 35)
censor_Radiated<- c(rep(1,9), rep(0,1))
censor_Radiated

time_Radiated_BPA <- c(31, 32, 34, 35, 36, 38, 38, 39, 42, 42)
censor_Radiated_BPA <- c(rep(1,8), rep(0,2))
censor_Radiated_BPA

myData <- data.frame(time=c(time_Untreated, time_Radiated, time_Radiated_BPA),
                     status=c(censor_Untreated, censor_Radiated, censor_Radiated_BPA),
                     group= rep(1:3, each=10))

library(KMsurv)
library(survival)

I have tried to use the function: pairwise_survdiff but I could not build a code on it.

Also, I want to perform the test for trend which would test this ordered hypothesis (untreated animals will have the worst survival, radiated rats will have slightly improved survival, and the radiated rats+BPA should have the best survival.)

Here is what I have done with the output, but I am unsure which is the value of Chi square and for p-value:
Is this correct?

KM.fit<-survfit(Surv(time,status)~group, conf.type="none", data=myData)
KM.fit

Call: survfit(formula = Surv(time, status) ~ group, data = myData, conf.type = "none") 

         n events median
group=1 10     10     25
group=2 10      9     30
group=3 10      8     37

myData.fit<-ten(Surv(time,status)~group, data=myData)
comp(myData.fit, p=0, q=0,scores =c(1,2,3))
            chiSq df     pChisq    
1          33.380  2 5.6436e-08 ***
n          30.255  2 2.6925e-07 ***
sqrtN      32.037  2 1.1046e-07 ***
S1         29.657  2 3.6307e-07 ***
S2         29.496  2 3.9349e-07 ***
FH_p=0_q=0 33.380  2 5.6436e-08 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
$tft
                   Q       Var      Z      pNorm    
1            16.0869    8.6116 5.4819 4.2081e-08 ***
n           364.0000 4741.0509 5.2864 1.2471e-07 ***
sqrtN        76.0224  196.2111 5.4272 5.7230e-08 ***
S1           11.1539    4.5558 5.2257 1.7351e-07 ***
S2           10.6871    4.2060 5.2110 1.8779e-07 ***
FH_p=0_q=0   16.0869    8.6116 5.4819 4.2081e-08 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

$scores
[1] 1 2 3

Answer 1

Please remember to specify all required packages when posting a question. You've omitted the survminer package, which is required for using the pairwise_survdiff function.

library(survminer)

The following code works on your dataset, myData, so I'm not sure what code you tried.

pairwise_survdiff(Surv(time,status)~group, data=myData)

---

    Pairwise comparisons using Log-Rank test 

data:  myData and group 

  1       2     
2 0.0011  -     
3 9.7e-06 0.0014

P value adjustment method: BH  # Bonferroni-Holm method of adjustment (default)

---

So all three groups have a significantly different survival.

The group variable should be converted into a factor, not just for labeling purposes on survival curves, but many modelling functions will assume the variable is continuous otherwise, which is clearly not the case here.

myData$group <- factor(myData$group, labels=c("Untreated","Radiated","Rad+BPA"))

KM.fit <- survfit(Surv(time,status)~group, data=myData)

ggsurvplot(KM.fit, data=myData,

       # Add median survival times (horizontal and vertical)
       surv.median.line = "hv", 

       # Legend placement, title, and labels
       legend = c(0.25, 0.75),
       legend.title = "Treatment Group",
       legend.labs = c("Untreated","Radiated","Rad+BPA"),
       # Add p-value and 95% confidence intervals
       pval = TRUE,
       conf.int = TRUE,

       # Add risk table
       risk.table = TRUE,
       tables.theme = theme_cleantable(),

       # Color palette
       palette = c("green4", "orange", "red"),
       ggtheme = theme_bw()
)

---

myData.fit <- ten(Surv(time,status)~group, data=myData)
comp(myData.fit, p=0, q=0, scores=1:3) # Scores 1, 2, 3 will be the default 

            chiSq df     pChisq    
1          33.380  2 5.6436e-08 ***  # log-rank test
n          30.255  2 2.6925e-07 ***  # Gehan-Breslow generalized Wilcoxon
sqrtN      32.037  2 1.1046e-07 ***  # Tarone-Ware
S1         29.657  2 3.6307e-07 ***  # Peto-Peto’s modified survival estimate
S2         29.496  2 3.9349e-07 ***  # modified Peto-Peto (by Andersen)
FH_p=0_q=0 33.380  2 5.6436e-08 ***  # Fleming-Harrington

Which is the value of the Chi-square and p-value? All of them are! The function gives you six to choose from. Here they are all significant. See the package documentation for more details.

The $tst section (not shown) gives the results of the test for trends. All comparisons and tests for trends indicate that there is a statistically significant difference in the survival of the rats in the three groups. Untreated rats have the worst survival (median=25 days), followed by radiated rats (median=30 days) and radiated+BPA (median=37 days).

---