进行成对检验和生存趋势检验

Performing Pairwise test and test for survival trend

我正在尝试执行成对测试以确定成对组之间的生存率是否存在差异。

使用的数据:

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)

我尝试使用函数:pairwise_survdiff,但我无法在其上构建代码。

此外,我想进行趋势检验以检验这个有序假设(未经治疗的动物存活率最差,辐射大鼠的存活率略有提高,而辐射大鼠+BPA 的存活率应该最高。 )

这是我对输出所做的,但我不确定哪个是卡方值和 p 值:
这个对吗?

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

请记住在发布问题时指定所有必需的包。您省略了 survminer 包,这是使用 pairwise_survdiff 函数所必需的。

library(survminer)

以下代码适用于您的数据集 myData,因此我不确定您尝试了哪些代码。

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)

因此所有三组的存活率都有显着差异。

应该将组变量转换为一个因子,这不仅是为了在生存曲线上进行标记,而且许多建模函数会假设变量是连续的,否则这里显然不是这种情况。

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

卡方值和p值是多少?他们都是!该功能为您提供六种选择。在这里它们都很重要。有关详细信息,请参阅 package documentation

$tst 部分(未显示)给出了趋势测试的结果。所有趋势比较和检验表明,三组大鼠的存活率存在统计学上的显着差异。未经治疗的大鼠存活率最差(中位数=25 天),其次是辐射大鼠(中位数=30 天)和辐射+BPA(中位数=37 天)。