Tutorial on survival modeling with applications to omics data
关于生存建模与组学数据应用的教程

Abstract 抽象

1 Introduction 1 引言

Personalized medicine improves patient diagnosis and treatment by making use of patient-specific genomic and molecular markers that are indicative of disease development. Time to an event of interest (e.g. time to death or disease progression) is a widely used end point and patient outcome for many diseases, and therefore it has become popular to identify genomic, molecular and clinical markers for survival or progression prediction of patients suffering from complex diseases such as cancer. In this tutorial, we will only consider so-called right-censored survival data, where a patient has been followed for a certain time period and the event of interest is either observed in this time period or might occur at a later (as yet unobserved) time point. Right-censored survival data include both a time and the status of each patient at that time as joint outcomes, where time is a continuous variable and status is a binary variable indicating whether the event of interest has been observed up to the given time or not; in the latter case, we refer to this observation as censored at the observed time. See Table 1 for an example illustration, and more details of the survival data in Section 2. When using the status label as an outcome in an ordinary logistic regression, the regression coefficients will become increasingly uncertain and less reliable with increasing follow-up time (Green and Symons 1983). When using the observed time (or its transformation, e.g. logarithm of time) as an outcome in an ordinary linear regression, the presence of censored observations (i.e. patients still alive by the end of follow-up period) causes considerable difficulties for assessing the accuracy of predictions (Henderson 1995). For example, Table 1 shows an example where the proportion of patients surviving past 10 years is $1/5=20%$ based on the observed data, but the (unobserved) factual proportion surviving past 10 years is $3/5=60%$⁠.
个性化医疗通过利用指示疾病发展的患者特异性基因组和分子标志物来改善患者的诊断和治疗。发生感兴趣事件的时间（例如死亡或疾病进展的时间）是许多疾病广泛使用的终点和患者结果，因此识别基因组、分子和临床标志物以预测患有癌症等复杂疾病的患者的生存或进展已变得很流行。在本教程中，我们将只考虑所谓的右删失生存数据，其中患者已经被跟踪了一段时间，并且感兴趣的事件要么在该时间段内观察到，要么可能发生在较晚的（尚未观察到的）时间点。右删失生存数据包括时间和每个患者在该时间的状态作为联合结果，其中时间是一个连续变量，状态是一个二进制变量，表明在给定时间之前是否观察到感兴趣的事件;在后一种情况下，我们将此观测称为在观测时间删失。有关示例说明，请参见表 1，并在第 2 节中有关生存数据的更多详细信息。当在普通 logistic 回归中使用状态标签作为结果时，随着随访时间的增加，回归系数将变得越来越不确定和不可靠 （Green 和 Symons 1983）。当使用观测时间（或其变换，例如 时间的对数）作为普通线性回归的结果时，删失观察的存在（即 患者在随访期结束时仍然活着）对评估预测的准确性造成相当大的困难（Henderson 1995）。 例如，表 1 显示了一个例子，根据观察到的数据，存活超过 10 年的患者比例为 $1/5=20%，$但（未观察到的）存活超过 10 年的事实比例为 $3/5=60%。$

The availability of multiple types of genomic and molecular data poses great opportunities but also further challenges for building effective statistical models to identify biomarkers that are prognostic of patient survival. For example, omics profiles, such as those from mRNA expression, DNA copy number, single-point and other genetic mutations, may be available from the same patient, and these high-dimensional data come with intra- and inter-dataset correlations, heterogeneous measurement scales, missing values, technical variability, and other background noise (Hasin et al. 2017). Rahnenführer et al. (2023) provided a general guideline for high-dimensional omics and electronic health records data analysis, and discussed statistical challenges and opportunities for survival modeling. Bøvelstad et al. (2007) and Bøvelstad et al. (2009) compared various machine learning methods and showed improved survival prediction performance by coefficient shrinkage methods that combine several data sources, in their case clinical and gene expression data. Another recent study by Bommert et al. (2022) performed a benchmark of filter methods for feature selection in high-dimensional gene expression data for survival prediction, and recommended using the simple variance filter before fitting a $ℓ2$-regularized Cox regression model (Simon et al. 2011) for accurate survival prediction and stable feature selection. Vinga (2021) reviewed structured penalized regressions for analyzing high-dimensional omics data for survival prediction or evaluation. Multiple studies have shown that it is possible to further improve the prediction accuracy and feature selection by considering more complex structures, such as biological pathways, or by identifying significant features among functional relationships between the omics features (Chekouo et al. 2015, Kundu et al. 2018, Wang and Liu 2020, Madjar et al. 2021).
多种类型的基因组和分子数据的可用性带来了巨大的机遇，但也为构建有效的统计模型以识别预测患者生存的生物标志物带来了进一步的挑战。例如，组学图谱，例如来自 mRNA 表达、DNA 拷贝数、单点和其他基因突变的组学图谱，可能来自同一患者，并且这些高维数据带有数据集内和数据集间的相关性、异质性测量尺度、缺失值、技术可变性和其他背景噪声（Hasin 等人，2017 年）。Rahnenführer 等人（2023 年）为高维组学和电子健康记录数据分析提供了通用指南，并讨论了生存建模的统计挑战和机遇。Bøvelstad 等人 （2007） 和 Bøvelstad 等人 （2009） 比较了各种机器学习方法，并表明结合多个数据源的系数收缩方法提高了生存预测性能，在他们的案例中是临床和基因表达数据。Bommert 等人 （2022） 最近的另一项研究对高维基因表达数据中的特征选择进行了过滤方法的基准测试，用于生存预测，并建议在拟合 $l2$ 正则化 Cox 回归模型之前使用简单的方差过滤器 （Simon et al. 2011） 进行准确的生存预测和稳定的特征选择。Vinga （2021） 回顾了用于分析高维组学数据以进行生存预测或评估的结构化惩罚回归。 多项研究表明，通过考虑更复杂的结构（例如生物途径），或通过识别组学特征之间的功能关系中的显着特征，可以进一步提高预测准确性和特征选择（Chekouo 等人，2015 年，Kundu 等人，2018 年，Wang 和 Liu 2020 年，Madjar 等人，2021 年）。

In this tutorial, we describe a general workflow for survival analysis with omics data, as well as review the commonly used statistical methods for feature selection and survival prediction; importantly, we provide a step-by-step R tutorial using publicly available omics data from The Cancer Genome Atlas (TCGA) project (http://cancergenome.nih.gov). In this example dataset, the overall survival time, demographic and gene expression data from primary invasive breast cancer patients in TCGA (TCGA-BRCA) were retrieved from the Genomic Data Commons Data Portal data release 32.0–36.0. Compared to the previous reviews and benchmarks of survival models in bioinformatic applications, this tutorial provides a complete workflow ranging from data preparation to final calibrated models, with a particular focus on building survival models using high-dimensional omics data, and as such covers both the commonly used penalized regressions and Bayesian models for survival analysis with high-dimensional and generally noisy datasets. For this tutorial, we assume that readers have the knowledge of basic statistical methods. Terms beyond basic statistics are explained in corresponding text.
在本教程中，我们描述了使用组学数据进行生存分析的一般工作流程，并回顾了用于特征选择和生存预测的常用统计方法;重要的是，我们使用癌症基因组图谱 （TCGA） 项目 （http://cancergenome.nih.gov） 中公开可用的组学数据提供了分步 R 教程。在此示例数据集中，从 Genomic Data Commons Data Portal 数据版本 32.0-36.0 中检索了 TCGA （TCGA-BRCA） 原发性浸润性乳腺癌患者的总生存时间、人口统计学和基因表达数据。与之前生物信息学应用中生存模型的综述和基准测试相比，本教程提供了从数据准备到最终校准模型的完整工作流程，特别侧重于使用高维组学数据构建生存模型，因此涵盖了常用的惩罚回归和贝叶斯模型，用于高维和通常噪声数据集的生存分析。在本教程中，我们假设读者具有基本统计方法的知识。基本统计以外的术语在相应的文本中进行了解释。

2 Data categories 2 数据类别

We use omics data and overall survival outcomes from cancer patients as an example in this tutorial, but the methods are applicable also to other diseases with similar data types in personalized medicine applications. The ultimate goal of personalized medicine is to identify patient-specific risk factors to guide disease prevention, diagnostic and treatment strategies. The identification of potential risk factors for cancer patients often considers clinical, demographic, genomic and molecular information, and their associations with the patients’ time-to-event data (i.e. survival).
在本教程中，我们使用癌症患者的组学数据和总生存期结果作为示例，但这些方法也适用于个性化医疗应用中具有类似数据类型的其他疾病。个性化医疗的最终目标是确定患者特异性的风险因素，以指导疾病预防、诊断和治疗策略。确定癌症患者的潜在风险因素通常考虑临床、人口统计学、基因组和分子信息，以及它们与患者事件发生时间数据（即生存率）的关联。

2.1 Time-to-event data 2.1 事件发生时间数据

Figure 1.

Pan-cancer survival data and omics signatures. (a) Multiple omics layers [modified from Haukaas et al. (2017) and Jendoubi (2021)]. The network illustrates intra- and cross-layered biological features or molecules (e.g. DNA methylation, mRNA, proteins, metabolites). (b) Different origins of tumors represented by the TCGA pan-cancer patient cohorts. (c) Overall survival times of cancer patients from TCGA. Novel methods are needed to model the high-dimensional multi-omics data and leverage information from heterogeneous cancer cohorts for improved survival prognosis and biomarker identification.

图 1.

泛癌生存数据和组学特征。（a） 多组学层 [修改自 Haukaas 等人 （2017） 和 Jendoubi （2021）]。该网络说明了层内和跨层的生物学特征或分子（例如 DNA 甲基化、mRNA、蛋白质、代谢物）。（b） TCGA 泛癌患者队列代表的肿瘤的不同来源。（c） TCGA 癌症患者的总生存时间。需要新的方法来模拟高维多组学数据，并利用来自异质性癌症队列的信息来改善生存预后和生物标志物鉴定。

Open in new tabDownload slide
在新选项卡中打开下载幻灯片

2.2 Clinical and demographic data
2.2 临床和人口统计数据

There are multiple sources of patient-level information that can be explored to identify risk factors for cancer patients. One can start with routinely collected and commonly used patient data, such as clinical and demographic variables. For example, an older male patient with a low body mass index (BMI) has a relatively high risk of gastric cancer (Nam et al. 2022). Table 2 illustrates selected clinical and demographic variables often available for cancer patients. Clinical and demographic variables are considered important sources of information for predicting survival and are often used to build reference models for omics-based prognostic models (Herrmann et al. 2021).
可以探索多种患者层面的信息来源，以确定癌症患者的风险因素。可以从常规收集和常用的患者数据开始，例如临床和人口统计变量。例如，体重指数 （BMI） 低的老年男性患者患胃癌的风险相对较高（Nam 等人，2022 年）。表 2 说明了癌症患者通常可用的选定临床和人口统计变量。临床和人口统计变量被认为是预测生存的重要信息来源，通常用于为基于组学的预后模型构建参考模型（Herrmann 等人，2021 年）。

2.3 Omics data 2.3 组学数据

Thanks to the rapid development of modern molecular biotechnology, large amounts of human genomic and molecular data have become available from many patient profiling projects. These projects often collect multiple levels of molecular information such as genomics data for DNA variation, transcriptomics data for mRNA expression, proteomics data for protein abundance and metabolomics data for metabolite processes, as illustrated in Fig. 1a. Among the multiple omics levels, metabolomics is the closest to observable phenotypes, such as tumor growth and proliferation (Cairns et al. 2011). To deeply understand the molecular biology of tumor development, multiple levels of omics data may deliver novel insights into the circuits of molecular interactions that underlie the disease initiation and progression (Tarazona et al. 2021).
由于现代分子生物技术的快速发展，许多患者分析项目已经提供了大量的人类基因组和分子数据。这些项目通常收集多个层次的分子信息，例如 DNA 变异的基因组学数据、mRNA 表达的转录组学数据、蛋白质丰度的蛋白质组学数据和代谢物过程的代谢组学数据，如图 1a 所示 。在多个组学水平中，代谢组学是最接近可观察的表型的，例如肿瘤生长和增殖（Cairns 等人，2011 年）。为了深入了解肿瘤发展的分子生物学，多层次的组学数据可能会为构成疾病发生和进展基础的分子相互作用回路提供新的见解（Tarazona 等人，2021 年）。

DNA-level omics data often include single nucleotide polymorphisms (SNP), DNA methylation, and somatic copy number variation, see illustration of these data types in Table 3. Each SNP feature can be coded as ${0,1,2}$⁠, according to the number of minor alleles at the given locus in each individual, i.e. AA, AG, GC, or TT. DNA methylation reveals methyl groups added to the DNA molecule, which can be quantified as a $β$-value $β=MM+U+a∈[0,1]$⁠, where M and U are the fluorescence intensities of the methylated and unmethylated DNA at a CpG locus, and the offset a is often set to 100 to stabilize the $β$-value when M and U are small. Somatic copy number variation measures the number of repeated sections of the tumor genome. Table 3 shows the assumed distributions for the downstream statistical modeling. For example, DNA methylation $β$-value is a proportion, often either close to 1 or 0, which can be characterized by a mixture of two beta distributions.
DNA 水平组学数据通常包括单核苷酸多态性 （SNP）、DNA 甲基化和体细胞拷贝数变异，参见表 3 中这些数据类型的插图。每个 SNP 特征都可以编码为 ${0,1,2}$，具体取决于每个个体中给定基因座的次要等位基因数量，即 AA、AG、GC 或 TT。DNA 甲基化揭示了添加到 DNA 分子中的甲基，可以量化为 -值 ，其中 M 和 U是甲基化和未甲基化 DNA 在 CpG 位点的荧光强度，当 M 和 U 较小时，偏移量 a 通常设置为 100 以稳定 -值。体细胞拷贝数变异测量肿瘤基因组重复切片的数量。表 3 显示了下游统计建模的假设分布。例如，DNA 甲基化 $β$ 值是一个比例，通常接近 1 或 0，这可以通过两个 β 分布的混合来表征。

RNA-level omics data usually include messenger RNA (mRNA) expression and microRNA (miRNA) expression. Traditionally, DNA microarrays were used to measure the expression levels of DNA sequences called probes, which acted as a proxy for the amount of reads representing a genomic feature of interest. The reads of microarray expression level can be characterized by a negative binomial distribution directly, or a Gaussian distribution after log₂ transformation (see Table 3). Nowadays, RNA sequencing (RNA-seq) has replaced DNA microarrays, since it allows for the sequencing of the whole transcriptome, while DNA microarrays only profile predefined transcripts or gene probes. miRNAs are small noncoding regulatory RNAs that play an important role in regulating gene expression and are highly evolutionary conserved in mammals (Bartel 2018).
RNA 水平组学数据通常包括信使 RNA （mRNA） 表达和 microRNA （miRNA） 表达。传统上，DNA 微阵列用于测量称为探针的 DNA 序列的表达水平，探针充当代表目标基因组特征的读数量的代理。微阵列表达水平的读数可以直接用负二项分布或 log₂ 转换后的高斯分布来表征（见表 3）。如今，RNA 测序 （RNA-seq） 已经取代了 DNA 微阵列，因为它允许对整个转录组进行测序，而 DNA 微阵列仅分析预定义的转录本或基因探针。miRNA 是小的非编码调节 RNA，在调节基因表达中起重要作用，并且在哺乳动物中高度进化保守 （Bartel 2018）。

Protein-level omics data usually originates from mass spectrometry (MS)-based proteomics profiling of a particular cell, tissue, organ, which can detect and measure abundance levels of entire or phosphorylated proteins. In contrast to global proteomics with MS, TCGA consortium has produced more targeted proteomics profiles, involving a set of 180–250 protein features using reverse-phase protein arrays (RPPA) (Akbani et al. 2014). In contrast, The Clinical Proteomic Tumor Analysis Consortium (CPTAC) profiling has produced >10 000 protein features using MS technology (Edwards et al. 2015). Protein expression data can often be considered approximately Gaussian distributed after log2 transformation, depending on the data-generating platform. In global proteomics, there are often sample-specific missing values due to detection limits of protein quantification.
蛋白质水平组学数据通常来源于基于质谱 （MS） 的特定细胞、组织、器官的蛋白质组学分析，可以检测和测量完整或磷酸化蛋白质的丰度水平。与使用 MS 的全球蛋白质组学相比，TCGA 联盟产生了更具针对性的蛋白质组学概况，涉及使用反相蛋白质阵列 （RPPA） 的一组 180-250 个蛋白质特征（Akbani 等人，2014 年）。相比之下，临床蛋白质组学肿瘤分析联盟 （CPTAC） 分析使用 MS 技术产生了 >10 000 个蛋白质特征（Edwards 等人，2015 年）。根据数据生成平台，在 log2 转换后，蛋白质表达数据通常可以被认为是近似高斯分布的。在全球蛋白质组学中，由于蛋白质定量的检测限，经常存在样品特异性缺失值。

Metabolite-level omics data has similar statistical properties to proteomics data, and metabolomics profiling is usually done with MS-based technologies, enabling the detection, and quantification of many thousands of metabolic features simultaneously. Nuclear magnetic resonance (NMR) spectroscopy is the other main analytical technology to profile metabolic processes. An important limitation of NMR spectroscopy is its relatively low sensitivity, which may lead to relatively few detected metabolites. Metabolic concentrations often follow logarithmic Gaussian distribution. Similar to large-scale proteomics, missing values due to detection limits should be treated differently than missing values due to measurement artefacts, which are more frequent in metabolite-level omics data (Sun and Xia 2023).
代谢物水平的组学数据具有与蛋白质组学数据相似的统计特性，代谢组学分析通常使用基于 MS 的技术完成，能够同时检测和定量数千个代谢特征。核磁共振 （NMR） 波谱是描述代谢过程的另一种主要分析技术。NMR 波谱的一个重要局限性是其灵敏度相对较低，这可能导致检测到的代谢物相对较少。代谢浓度通常服从对数高斯分布。与大规模蛋白质组学类似，应区别对待由于检测限导致的缺失值与由于测量伪影导致的缺失值，后者在代谢物水平的组学数据中更为常见（Sun 和 Xia 2023）。

Single-cell sequencing is becoming increasingly more prevalent in many profiling studies. The modern single-cell omics technologies can produce multiple levels of omics data derived from the same samples, such as transcriptomics and chromatin accessibility. The single-cell data types are similar to the illustration in Table 3, but each measurement originates from the level of individual cells. This article focuses mainly on survival analysis with single bulk omics data.
单细胞测序在许多分析研究中变得越来越普遍。现代单细胞组学技术可以产生来自相同样本的多层次组学数据，例如转录组学和染色质可及性。单单元数据类型与表 3 中的插图类似，但每个测量值都源自单个单元的级别。本文主要关注使用单个批量组学数据的生存分析。

2.4 Missing data 2.4 缺失数据

Missing values are often observed in many types of high-dimensional omics data due to various experimental reasons (Aittokallio 2009, Kong et al. 2022). For example, mRNA transcriptomics data from microarrays have 1%–10% missing values affecting up to 95% of genes due to corruption of image, scratches on the slides, poor hybridization, inadequate resolution, fabrication errors (de Brevern et al. 2004, Tuikkala et al. 2005); MS-based metabolomics data have 10%–20% missing values affecting up to 80% of variables due to lack of peak identification by chromatogram, limitation of computational detection, measurement error, and deconvolution errors (Hrydziuszko and Viant 2012, Sun and Xia 2023). The aforementioned technical reasons can lead to missing data that is either missing at random (MAR) or missing not at random (MNAR). When dealing with MNAR data, traditional imputation methods like multiple imputation may introduce bias. It is thus recommended to remove omics features which have large proportion (e.g. 50%) missingness over patients, and then apply imputation methods (e.g. k-nearest neighbor imputation) for the rest of the features with missing values before doing any statistical analysis or modeling. Alternatively, imputation-free methods (e.g. mixture models) that can deal with missing values can be applied directly (Taylor et al. 2022). Single-cell RNA-seq (scRNA-seq) data has a vast number of zeros, so-called gene dropout events, leading highly scarce data matrices. Jiang et al. (2022) discussed the sources of biological and nonbiological zeros in scRNA-seq data and the existing approaches for handling them.
由于各种实验原因，在许多类型的高维组学数据中经常观察到缺失值（Aittokallio 2009，Kong et al. 2022）。例如，由于图像损坏、载玻片上的划痕、杂交不良、分辨率不足、制造错误，来自微阵列的 mRNA 转录组学数据有 1%-10% 的缺失值，影响高达 95% 的基因（de Brevern 等人，2004 年，Tuikkala 等人，2005 年）;由于无法通过色谱图鉴定峰、计算检测的限制、测量误差和反卷积误差，基于 MS 的代谢组学数据存在 10%-20% 的缺失值，影响高达 80% 的变量（Hrydziuszko 和 Viant 2012，Sun 和 Xia 2023）。上述技术原因可能导致随机缺失 （MAR） 或非随机缺失 （MNAR） 缺失数据。在处理 MNAR 数据时，传统的插补方法（如多重插补）可能会引入偏倚。因此，建议删除患者缺失比例较大（例如 50%）的组学特征，然后在进行任何统计分析或建模之前对其余缺失值的特征应用插补方法（例如 k 最近邻插补）。或者，可以直接应用可以处理缺失值的无插补方法（例如 混合模型）（Taylor 等人，2022 年）。单细胞 RNA-seq （scRNA-seq） 数据有大量的零，即所谓的基因缺失事件，导致数据矩阵高度稀缺。江 et al. （2022） 讨论了 scRNA-seq 数据中生物和非生物零的来源以及处理它们的现有方法。

3 Survival analysis with low-dimensional input data
3 使用低维输入数据进行生存分析

Let us assume we have data $D={(T˜i,δi,Xi):i=1,…,n}$ for n patients, where $T˜i$ is the observed survival time, $δi$ the censoring indicator and $Xi$ contains p covariates including clinical, demographic and omics features. To estimate a survival function $S(t)$ given the data $D$⁠, one needs to keep track both of the number of patients at risk and those who left the study at any time point (here we only consider the case of right censoring and assume no delayed entry). At time t there are $Y(t)=∑i=1nYi(t)$ patients at risk, where $Yi(t)=1{T˜i≥t}$ is the indicator that patient i is at risk. The nonparametric Kaplan–Meier (KM) estimator (Kaplan and Meier 1958) uses the multiplication rule for conditional probabilities to obtain an estimation of the survival function $S^(t)=∏k=1Tk≤tK{1−dkY(Tk)},$ where all events happen (e.g. patients die) at K distinct times, $T1<T2<⋯<TK (K≤n)$ and there are $dk≥1$ events happened at the time $Tk$⁠. If no two events happen at the same time point, $dk=1$ and $k=1,…,n$⁠. The KM estimator gives an estimate of the marginal survival function, i.e. when you disregard the information from the covariates.
假设我们有数据 表示 n 名患者，其中 $T ̃i$ 是观察到的生存时间，$δi$ 是删失指标，$Xi$包含 p 个协变量，包括临床、人口统计学和组学特征。为了在给定数据 $D$ 的情况下估计生存函数 $S（t），$需要跟踪处于风险中的患者数量和在任何时间点退出研究的患者数量（这里我们只考虑右删失的情况，并假设没有延迟进入）。 在时间 t 有 $Y（t）=∑i=1nYi（t）$ 患者处于风险中，其中 $Yi（t）=1{T ̃i≥t}$ 是患者 i 处于风险中的指标。 非参数 Kaplan-Meier （KM） 估计器 （Kaplan 和 Meier 1958） 使用条件概率的乘法规则来获得生存函数 $S^（t）=∏k=1 的估计值Tk≤tK{1−dkY（Tk）}，$所有事件都发生的地方（例如 患者死亡）在 K 个不同时间，$T1<T2<⋯<TK （K≤n）$ 和有 $dk≥1$ 事件发生在 $Tk$ 时间。如果同一时间点没有两个事件发生，$dk=1$ 和 $k=1，...，n$.KM 估计器给出了边际生存函数的估计值，即当您忽略协变量的信息时。

Figure 2a shows the KM curve for TCGA-BRCA primary breast cancer patients data. Some basic statistics can be revealed from the survival curve. For example, the estimated median survival time, i.e. the time when the survival probability is 50%, of all the breast cancer patients is 10.8 years (dashed line in Fig. 2a), and 1-, 5-, and 10-year survival probabilities are 0.988, 0.853, and 0.658, respectively. A log-rank test (Peto and Peto 1972) can be used to test whether two groups of patients [e.g. with treatment (pharmaceutical/radiation therapy) or nontreatment] have the same (null hypothesis) or different survival functions (alternative hypothesis), and provide a corresponding P-value (see Fig. 2b). The log-rank test can also be used to compare the survival probabilities of any subgroups of patients based on other categorical variables.
图 2a 显示了 TCGA-BRCA 原发性乳腺癌患者的 KM 曲线数据。从生存曲线中可以揭示一些基本统计数据。例如，所有乳腺癌患者的估计中位生存时间，即生存概率为 50% 的时间为 10.8 年（图 2a 中的虚线），1 年、5 年和 10 年生存概率分别为 0.988、0.853 和 0.658。对数秩检验（Peto 和 Peto 1972）可用于检验两组患者 [例如接受治疗（药物治疗/放射治疗）或未接受治疗] 是否具有相同（零假设）或不同的生存函数（替代假设），并提供相应的 P 值（见图 2b）。对数秩检验还可用于根据其他分类变量比较任何患者亚组的生存概率。

Figure 2.

Kaplan–Meier curves of TCGA-BRCA data. (a) KM curve of all the TCGA-BRCA patients’ survival data. (b) KM curves of the TCGA-BRCA patients’ survival data grouped by treatment (i.e. pharmaceutical/radiation therapy) or nontreatment. The log-rank test is used to compare the two survival distributions corresponding to the two groups of patients.

图 2.

TCGA-BRCA 数据的 Kaplan-Meier 曲线。（a） 所有 TCGA-BRCA 患者生存数据的 KM 曲线。（b） TCGA-BRCA 患者生存数据的 KM 曲线，按治疗（即 药物/放射治疗）或非治疗。log-rank 检验用于比较对应于两组患者的两个生存分布。

Open in new tabDownload slide
在新选项卡中打开下载幻灯片

It is often of interest to gain insights into the multiple factors and their cooperation for the survival outcomes. Multivariable statistical analysis plays an important role in such multi-modal survival modeling. The univariate Cox model (1) can be straightforwardly generalized by including multiple clinical, demographic or omics variables of interest, as long as the total number of covariates is much smaller than the number of samples; this often requires the use of variable or feature selection methods. Heinze et al. (2018) provides pragmatic recommendations for practitioners of multivariable analysis with variable selection methods in low-dimensional modeling problems.
深入了解多种因素及其对生存结果的合作通常是有趣的。多变量统计分析在这种多模态生存建模中起着重要作用。单变量 Cox 模型 （1） 可以通过包括多个感兴趣的临床、人口统计学或组学变量来直接推广，只要协变量的总数远小于样本的数量;这通常需要使用变量或特征选择方法。Heinze et al. （2018） 为在低维建模问题中使用变量选择方法进行多变量分析的从业者提供了务实的建议。

4 Survival analysis with high-dimensional input data
4 使用高维输入数据的生存分析

Since some omics data contain high-number of variables (e.g. RNA-seq with ca. 60 000 transcriptomic features and DNA methylation with ca. 450 000 features), there is a need to reduce the computational and modeling burden in multivariable analyses. One heuristics approach to pre-select a subset of features is to include only omics features at a specific statistical significance level when fitting a univariate Cox model (1); the pre-specified significance level would usually be higher than the commonly used 0.05 threshold, e.g. 0.1 or 0.2, to avoid losing important features. However, this univariate approach focuses on features that are independently associated with the outcome and might miss variables that are important in combination with other features (Okser et al. 2013).
由于一些组学数据包含大量变量（例如，具有约 60 000 个转录组特征的 RNA-seq 和具有约 450 000 个特征的 DNA 甲基化），因此需要减少多变量分析中的计算和建模负担。一种预先选择特征子集的启发式方法是在拟合单变量 Cox 模型时仅包括特定统计显着性水平的组学特征 （1）;预先设定的显著性水平通常会高于常用的 0.05 阈值，例如 0.1 或 0.2，以避免丢失重要特征。然而，这种单变量方法侧重于与结果独立相关的特征，并且可能会错过与其他特征相结合的重要变量（Okser 等人，2013 年）。

Another simple approach is to pre-select omics features by variance, since larger variability across patients usually implies higher biological information, or at least predictive signal. One can for example pre-select omics features explaining 50% of the total cumulative variance (Zhao and Zucknick 2020). Such unsupervised feature pre-selection is a recommended method to reduce dimensionality when dealing with hundreds or thousands of omics features, with the aim to improve the stability of final feature selection. For example, Bommert et al. (2022) showed that the simple variance filter was the best method among all considered filter methods in terms of the predictive accuracy, run time and the feature selection stability in their benchmark study. Zhao et al. (2024) showed that the stability of final feature subset critically depends on the pre-selected feature set when using a standard Bayesian stochastic search variable selection method (see Section 4.3), and their proposed method that used known biological relationships between omics features lead to a more stable feature selection and slightly improved outcome prediction.
另一种简单的方法是通过方差预先选择组学特征，因为患者之间更大的变异性通常意味着更高的生物学信息，或者至少是预测信号。例如，可以预先选择组学特征来解释总累积方差的 50%（Zhao 和 Zucknick 2020）。这种无监督特征预选是在处理成百上千个组学特征时降低维度的推荐方法，目的是提高最终特征选择的稳定性。例如，Bommert 等人（2022 年）表明，在他们的基准研究中，就预测准确性、运行时间和特征选择稳定性而言，简单方差滤波器是所有考虑的滤波器方法中最好的方法。Zhao et al. （2024） 表明，当使用标准贝叶斯随机搜索变量选择方法时，最终特征子集的稳定性在很大程度上取决于预先选择的特征集（参见第 4.3 节），他们提出的方法使用组学特征之间已知的生物学关系导致更稳定的特征选择和略微改善的结果预测。

When drawing conclusions on survival differences solely from the univariate Cox model (1), it is important to adjust the P-values of risk features for multiple comparisons to control false positives globally, e.g. by controlling the family-wise error rate (FWER) or false discovery rate (FDR). The problem of multiple testing is beyond the scope of this tutorial, but the interested readers are referred to a recent review article (Korthauer et al. 2019). We note also that univariate analysis does not consider any relationships between multiple omics features, e.g. potential confounders for survival outcomes (Clark et al. 2003). To avoid making seriously misleading conclusions in such cases, it is necessary to perform multivariable survival analysis (Bradburn et al. 2003).
当仅从单变量 Cox 模型 （1） 得出关于生存差异的结论时，重要的是要调整风险特征的 P 值以进行多重比较以全局控制假阳性，例如通过控制家庭错误率 （FWER） 或错误发现率 （FDR）。多重测试的问题超出了本教程的范围，但感兴趣的读者可以参考最近的一篇评论文章（Korthauer et al. 2019）。我们还注意到，单变量分析不考虑多个组学特征之间的任何关系，例如 生存结果的潜在混杂因素 （Clark et al. 2003）。为避免在这种情况下做出严重误导性的结论，有必要进行多变量生存分析（Bradburn 等人，2003 年）。

Omics data can have hundreds of thousands of variables measured at various molecular levels, which greatly challenges the classical multivariable regression models for time-to-event endpoints, since the number of variables is often much larger than the number of patients (i.e. $p≫n$⁠). To proceed with survival analysis, one option is to reduce the dimensionality of the omics features via unsupervised learning, and then investigate the association of the learned low-dimensional variables with survival outcomes (see next section). An alternative approach is to directly use supervised learning methods, such as penalized regressions or sparse Bayesian models (Table 4), which enable the modeling of high-dimensional omics features, and the selection of key important features associated with the survival outcomes.
组学数据可以有数十万个在不同分子水平上测量的变量，这极大地挑战了用于事件发生时间终点的经典多变量回归模型，因为变量的数量通常比患者的数量大得多（即 $p≫n$）。要进行生存分析，一种选择是通过无监督学习来降低组学特征的维数，然后研究学习到的低维变量与生存结果的关联（见下一节）。另一种方法是直接使用监督学习方法，例如惩罚回归或稀疏贝叶斯模型（表 4），它们能够对高维组学特征进行建模，并选择与生存结果相关的关键重要特征。

4.1 Unsupervised learning
4.1 无监督学习

Unsupervised learning methods aim to identify hidden patterns or data groupings, and are for example useful when a phenotype (e.g. a cancer type) is to be divided into several subtypes (e.g. to explain heterogeneity among patients). For example, breast cancer has been traditionally categorized into five conceptual molecular classes, originally using pairwise average-linkage cluster analysis of DNA microarray data, to better understand tumor biology and guide clinical decision making (Perou et al. 2000). Unsupervised methods learn underlying patterns from unlabeled data by transforming high-dimensional omics features into a lower dimensional space. Principal component analysis (PCA) is a classical multivariate technique that represents high-dimensional features in a low-dimensional space by building orthogonal (uncorrelated) linear combinations of the features that best capture the variance in the high-dimensional data.
无监督学习方法旨在识别隐藏的模式或数据分组，例如，当表型（例如癌症类型）要分为几个亚型时（例如，解释患者之间的异质性）时，这种方法非常有用。例如，乳腺癌传统上被分为五个概念分子类别，最初使用 DNA 微阵列数据的成对平均连锁聚类分析，以更好地了解肿瘤生物学并指导临床决策（Perou 等人，2000 年）。无监督方法通过将高维组学特征转换为低维空间，从未标记的数据中学习基本模式。主成分分析 （PCA） 是一种经典的多元技术，它通过构建最能捕获高维数据中方差的特征的正交（不相关）线性组合来表示低维空间中的高维特征。

Different from the distance-based PCA with linear transformation, nonlinear techniques have recently emerged, such as t-stochastic neighbor embedding (t-SNE). t-SNE uses pairwise similarities of individuals based on Kullback–Leibler divergence to project similarities into a lower dimensional space, ensuring that individuals with similar omics features are close in the generated embedding (van der Maaten and Hinton 2008). The focus of t-SNE is on preserving local distances between neighboring data points. Another nonlinear alternative to PCA is UMAP (Uniform Manifold Approximation and Projection), which is a general dimension reduction method built upon Riemannian geometry (McInnes et al. 2018). Compared to other nonlinear methods of dimension reduction such as t-SNE, UMAP can sometimes provide better visualization quality in a shorter amount of time, while preserving the global structure of the omics data better.
与线性变换的基于距离的 PCA 不同，最近出现了非线性技术，例如 t 随机邻域嵌入 （t-SNE）。t-SNE 使用基于 Kullback-Leibler 散度的个体成对相似性，将相似性投影到较低维空间中，确保具有相似组学特征的个体在生成的嵌入中接近（van der Maaten 和 Hinton 2008）。t-SNE 的重点是保持相邻数据点之间的本地距离。PCA 的另一种非线性替代方案是 UMAP（均匀流形近似和投影），这是一种建立在黎曼几何基础上的通用降维方法（McInnes 等人，2018 年）。与其他非线性降维方法（如 t-SNE）相比，UMAP 有时可以在更短的时间内提供更好的可视化质量，同时更好地保留组学数据的全局结构。

Unsupervised methods only make use of the input data matrix $X$ and are agnostic to the survival information. To find out whether a given omics profile (i.e. the full set of omics features, not individual features) is associated with survival outcomes, a straightforward approach is to use a few representative components from PCA, t-SNE, or UMAP as covariates in a multivariable Cox model. An alternative is to use semi-supervised methods (Bair and Tibshirani 2004) that combine the clustering procedure and survival modeling together. However, such principal component regression and semi-supervised methods lose interpretability of the individual omics features, since each component is a linear or nonlinear combination of all omics features.
无监督方法仅使用输入数据矩阵，并且与生存信息无关。为了找出给定的组学概况（即完整的组学特征集，而不是单个特征）是否与生存结果相关，一种简单的方法是使用 PCA、t-SNE 或 UMAP 中的一些代表性成分作为多变量 Cox 模型中的协变量。另一种方法是使用半监督方法 （Bair 和 Tibshirani 2004），将聚类过程和生存建模结合在一起。然而，这种主成分回归和半监督方法失去了单个组学特征的可解释性，因为每个成分都是所有组学特征的线性或非线性组合。

4.2 Supervised learning via penalized regressions
4.2 通过惩罚回归进行监督学习

Figure 3.

Coefficient trace plot of an Elastic Net Cox model for overall survival prognosis of breast cancer patients from TCGA based on transcriptomic data and mandatory demographic variables age and ethnicity. The y-axis shows the magnitude of each feature’s coefficient given the strength of penalization displayed on the x-axis (from left to right the penalization decreases). The vertical line indicates the optimal $λ=0.032$ (maximizes the partial likelihood via cross-validation) and its corresponding selected feature names are listed on the right side. Note that the demographic variables age and ethnicity were not penalized, so that their coefficient paths did not start from zero in the figure.

图 3.

基于转录组数据和强制性人口统计变量年龄和种族的 TCGA 乳腺癌患者总生存预后的弹性网考克斯模型的系数轨迹图。y 轴显示给定 x 轴上显示的惩罚强度的每个特征系数的大小（从左到右，惩罚减少）。垂直线表示最佳 λ = 0.032 （通过交叉验证最大化部分似然）及其相应的所选特征名称列在右侧。请注意，人口统计变量 age 和 ethnicity 没有受到惩罚，因此它们的系数路径在图中不是从零开始的。

Open in new tabDownload slide
在新选项卡中打开下载幻灯片

There are many alternative penalties that will achieve feature selection which have been applied to Cox proportional hazards regression for survival outcomes, such as the Adaptive Lasso Cox model that incorporates different penalties for different coefficients to retain important variables (Zhang and Lu 2007), Group-Lasso that performs group selection on (predefined) groups of variables (Kim et al. 2012), Sparse Group-Lasso that introduces sparse effects both on a group and within group level (Simon et al. 2013), the smoothly clipped absolute deviation (SCAD) Cox model that overcomes substantially biased estimates for large coefficients in ultra-sparse models (Fan and Li 2002), sure independence screening (SIS) procedure in combination with Cox model that speeds-up the feature selection dramatically and can also improve the accuracy of estimation when dimensionality becomes ultra-high, i.e. $log(p)=O(nξ)$ for some $ξ>0$ (Fan et al. 2010). However, all these penalized Cox models do not directly provide uncertainty of feature selection or survival prediction. One empirical way for uncertainty quantification is through additional resampling-based methods, see Section 5 for more details.
有许多替代惩罚可以实现特征选择，这些惩罚已应用于生存结果的 Cox 比例风险回归，例如自适应套索 Cox 模型，该模型对不同系数进行了不同的惩罚以保留重要变量（Zhang 和 Lu 2007），对（预定义）变量组执行组选择的 Group-Lasso（Kim 等人，2012 年）、在组级别和组内引入稀疏效应的稀疏组-套索（Simon et al. 2013），平滑裁剪绝对偏差 （SCAD） Cox 模型克服了超稀疏模型中大系数的实质性偏差估计（Fan 和 Li 2002），确定独立性筛选 （SIS） 程序与 Cox 模型相结合，大大加快了特征选择的速度，并且还可以提高维数变得超高，即 $log（p）=O（nξ）$ 对于某些 $ξ>0$ （Fan et al. 2010）。然而，所有这些被惩罚的 Cox 模型并不直接提供特征选择或生存预测的不确定性。不确定性量化的一种实证方法是通过其他基于重采样的方法，有关详细信息，请参见第 5 节。

4.3 Supervised learning via Bayesian priors
4.3 通过贝叶斯先验进行监督学习

Bayesian inference is an appealing approach for survival analysis due to its ability to provide straightforward uncertainty quantification (e.g. credible intervals) of parameters and survival probabilities. For instance, Lee et al. (2011) proposed a Bayesian version of the Lasso Cox (Bayesian Lasso Cox) model that provides credible intervals of coefficients fairly straightforward (see Fig. 4), but it is not easy to derive confidence intervals of coefficients in a Lasso-type model.
贝叶斯推理是一种有吸引力的生存分析方法，因为它能够提供参数和生存概率的直接不确定性量化（例如可信的区间）。例如，Lee 等人（2011 年）提出了套索考克斯（Bayesian Lasso Cox）模型的贝叶斯版本，该模型提供了相当简单的可信系数区间（见图 4），但在套索型模型中推导出系数的置信区间并不容易。

Figure 4.

Point estimation and uncertainty quantification of regression coefficients using Bayesian Lasso Cox model with Laplace prior for overall survival prognosis of breast cancer patients from TCGA based on transcriptomic data and mandatory demographic variables age and ethnicity. Solid dots indicate the posterior mean over 20000 MCMC iterations excluding burn-in period, and horizontal lines show the corresponding 95% credible intervals.

图 4.

使用贝叶斯 Lasso Cox 模型和拉普拉斯先验根据转录组数据和强制性人口统计变量年龄和种族对 TCGA 乳腺癌患者总生存期预后的回归系数进行点估计和不确定性量化。实线表示 20000 次 MCMC 迭代的后验平均值，不包括老化期，水平线表示相应的 95% 可信区间。

Open in new tabDownload slide
在新选项卡中打开下载幻灯片

Figure 5.

Density of shrinkage priors in Bayesian survival modeling. (a) Density of the double exponential (Laplace) prior. (b) Density of the mixture spike-and-slab prior. The spike component $δ0(βj)$ induces $βj=0$ and the slab component $N(0,τj2)$ induces $βj≠0$⁠. (c) Density of the horseshoe prior. The shrinkage weight $κj$ close to 0 shrinks $βj$ toward zero, and $κj$ close to 1 allows $βj$ to escape the shrinkage.

图 5.

贝叶斯生存模型中的收缩先验密度。（a） 双指数 （Laplace） 先验的密度。（b） 混合物 spike-slab 的密度。尖峰分量 $δ0（βj）$ 导致$βj=0$，板分量 $N（0，τ j2）$ 诱导 $βJ≠0$。（c） 马蹄形先验的密度。接近 0 的收缩权重 $κj$ 会$βj$ 向零收缩，而接近 1 $的 κj$ 会允许 $βj$ 逃避收缩。

Open in new tabDownload slide
在新选项卡中打开下载幻灯片

Note that the Bayesian version of the Lasso with Laplace-type priors in practice do not result in automatic feature selection, because only the posterior modes of the coefficients are equivalent to the frequentist Lasso solution, while in Bayesian inference one usually focuses on the posterior means as point estimates. As an alternative, a particular omics feature can be excluded if the estimated credible interval of the corresponding coefficient covers zero.
请注意，在实践中，具有拉普拉斯类型先验的 Lasso 的贝叶斯版本不会导致自动特征选择，因为只有系数的后验模式等价于频率拉索解，而在贝叶斯推理中，人们通常将后验均值作为点估计。作为替代方案，如果相应系数的估计可信区间覆盖零，则可以排除特定的组学特征。

Although Bayesian models can quantify uncertainty of estimators more straightforward than penalized regressions, most Bayesian Cox-type models for high-dimensional covariates do not have user-friendly and standalone R packages on CRAN or GitHub. The main reason is the high computational cost of running a high-dimensional Bayesian Cox model. Advanced users with programming capabilities can contact the corresponding authors for original scripts. Since Bayesian priors are more flexible than frequentist Lasso-type penalties, it can be easier to extend Bayesian models by changing shrinkage priors while keeping almost the same algorithm framework. This means that the Bayesian framework can be more suitable if one is interested in tailoring the shrinkage effects, e.g. to include prior knowledge about the importance of omics features, e.g. features corresponding to a molecular pathway that is known to be affected in the disease under study. For more information on different Bayesian priors in cancer prognosis, we suggest a recent review by Chu et al. (2022) which summarized other different shrinkage priors on regression coefficients, such as Gaussian-gamma, Gaussian, Cauchy, pMOM (product moment distribution), piMOM (product inverse moment distribution) and peNMIG (parameter-expanded normal-mixture-of-inverse-gamma distribution) priors.
尽管贝叶斯模型可以比惩罚回归更直接地量化估计量的不确定性，但大多数用于高维协变量的贝叶斯 Cox 型模型在 CRAN 或 GitHub 上没有用户友好的独立 R 包。主要原因是运行高维 Bayesian Cox 模型的高计算成本。具有编程能力的高级用户可以联系相应的作者获取原始脚本。由于贝叶斯先验比频率 Lasso 类型的惩罚更灵活，因此在保持几乎相同的算法框架的同时，通过更改收缩先验来扩展贝叶斯模型会更容易。这意味着如果有人对定制收缩效应感兴趣，例如包括关于组学特征重要性的先验知识，例如对应于已知在所研究疾病中受影响的分子途径的特征，贝叶斯框架可能更合适。有关癌症预后中不同贝叶斯先验的更多信息，我们建议 Chu 等人（2022 年）最近的一篇综述，该综述总结了回归系数的其他不同收缩先验，例如高斯伽马、高斯、柯西、pMOM（乘积矩分布）、piMOM（乘积逆矩分布）和 peNMIG（参数扩展的正态混合逆伽马分布）先验。

5 Survival model validation
5 生存模型验证

Model validation plays an important role in identifying potential issues, such as model misspecification or overfitting. This is achieved by revisiting the model’s specifications and assumptions following model estimation. For example, the Cox model (1) requires proportional hazards and the logarithm of the hazard to be linear with respect to the model covariates. The former assumption can for example be checked by the cumulative Schoenfeld residuals (Grambsch and Therneau 1994), and the latter assumption can be checked by plotting a nonlinear functional form (e.g. spline) for the effect of a covariate. If the Cox model assumptions are not satisfied, one can try certain transformations of covariates (e.g. Box-Cox power transformations, Box and Cox 1964), allow time-varying coefficients or model interactions among covariates (Ng et al. 2023), or investigate patient stratification using unsupervised approaches (Cristescu et al. 2015). However, the assumption checks are usually suitable only for low-dimensional models, i.e. for a few clinical variables or a few factors projected from the high-dimensional omics feature space. Novel approaches for assumption checks in general high-dimensional settings require further methodological developments. Johnson and Long (2011) used heuristic methods to investigate the Cox model assumptions by separately fitting univariate Cox models one feature at one time, and check P-values for the score tests of individual features and P-values for testing the proportional hazards assumption of univariate Cox models. But the univariate Cox models do not take into account confounding variables. An alternative approach is to loosen the model assumptions for a more robust modeling approach. One example developed specifically for feature selection under possibly nonproportional hazards in a high-dimensional space is concordance regression (Dunkler et al. 2010).
模型验证在识别潜在问题（例如模型错误指定或过度拟合）方面起着重要作用。这是通过在模型估计后重新审视模型的规范和假设来实现的。例如，Cox 模型 （1） 要求比例风险和风险的对数相对于模型协变量呈线性。例如，前一个假设可以通过累积 Schoenfeld 残差 （Grambsch and Therneau 1994） 来检查，后一个假设可以通过绘制协变量效应的非线性函数形式（例如样条）来检查。如果不满足 Cox 模型假设，可以尝试协变量的某些变换（例如 Box-Cox 幂变换、Box 和 Cox 1964），允许协变量之间的时变系数或模型交互（Ng 等人，2023 年），或使用无监督方法研究患者分层（Cristescu 等人，2015 年).然而，假设检查通常仅适用于低维模型，即从高维组学特征空间投射的少数临床变量或少数因素。在一般高维设置中进行假设检查的新方法需要进一步的方法论发展。Johnson 和 Long （2011） 使用启发式方法通过一次单独拟合单变量 Cox 模型一个特征来研究 Cox 模型假设，并检查单个特征的分数测试的 P 值和用于测试单变量 Cox 模型的比例风险假设的 P 值。但是单变量 Cox 模型没有考虑混杂变量。 另一种方法是放宽模型假设，以获得更稳健的建模方法。专门为高维空间中可能的非比例风险下的特征选择开发的一个例子是一致性回归（Dunkler 等人，2010 年）。

5.1 Feature stability analysis
5.1 特征稳定性分析

One important aspect in model validation when using omics or other high-dimensional data is the potential instability of feature selection (Kalousis et al. 2007). Feature selection using penalized regressions as described in Section 4.2 heavily depends on the values of the penalty parameters [e.g. for the $λ$ parameter in Lasso Cox model (2)]. The penalty parameters are often optimized by cross-validation (CV) or other resampling methods, and the uncertainty associated with the random selection of subsets may result in uncertainty in the feature selection, e.g. different CV folds will typically result in different selected features. A straightforward way to identify the most stable features is to find the overlap of identified omics features among different data subsets (e.g. CV folds or resamples) to avoid high false discovery rate (Zucknick et al. 2008). One can also perform stability selection (Meinshausen and Bühlmann 2010), which allows to select the most stable features at a given Type-I error level for a Lasso or Elastic Net Cox model (Sill et al. 2014).
使用组学或其他高维数据时，模型验证的一个重要方面是特征选择的潜在不稳定性（Kalousis 等人，2007 年）。如第 4.2 节所述，使用惩罚回归的特征选择在很大程度上取决于惩罚参数的值 [例如，对于 Lasso Cox 模型 （2） 中的参数]。惩罚参数通常通过交叉验证 （CV） 或其他重采样方法进行优化，与随机选择子集相关的不确定性可能会导致特征选择的不确定性，例如，不同的 CV 折叠通常会导致不同的选定特征。识别最稳定特征的一种直接方法是找到不同数据子集（例如 CV 折叠或重新采样）之间已识别的组学特征的重叠，以避免高错误发现率（Zucknick 等人，2008 年）。还可以进行稳定性选择（Meinshausen 和 Bühlmann 2010），这允许在给定的 I 类误差级别为 Lasso 或 Elastic Net Cox 模型选择最稳定的特征（Sill 等人，2014 年）。

For the Bayesian models in Section 4.3, feature selection stability is naturally assessed by the uncertainty of coefficients’ estimators, as reflected in the posterior variances of $βj$ or the posterior selection probabilities $p(γj|D)$ (in SSVS), which is a natural benefit of utilizing full Bayesian inference. Although the uncertainty in feature selection introduces increased variability in the predicted survival probabilities, in the Bayesian framework, this can be addressed quite naturally by averaging the survival predictions over all models using Bayesian model averaging (Volinsky et al. 1997). If one is interested in a single model, rather than model averaging, the median probability model (Barbieri and Berger 2004) can be used for uncertainty analyses in survival and high-dimensional omics data (Madjar et al. 2021).
对于第 4.3 节中的贝叶斯模型，特征选择稳定性自然是通过系数估计量的不确定性来评估的，这反映在 $βj$ 的后验方差或后验选择概率 $p（γj|D）$ （在 SSVS 中），这是利用完整贝叶斯推理的自然好处。尽管特征选择的不确定性增加了预测生存概率的可变性，但在贝叶斯框架中，这可以通过使用贝叶斯模型平均对所有模型的生存预测进行平均来自然地解决（Volinsky 等人，1997 年）。如果对单个模型感兴趣，而不是模型平均，中位概率模型（Barbieri 和 Berger 2004）可用于生存和高维组学数据的不确定性分析（Madjar 等人，2021 年）。

5.2 Survival prediction and calibration
5.2 生存预测和校准

The fundamental goal of any statistical prediction model is to achieve a better prediction performance than an existing statistical model which we could call the ‘conventional model’ (Gerds and Kattan 2021). The special nature of the combination of clinical and/or other known prognostic factors (typically low-dimensional, with established effect) and novel omics features (high-dimensional, with unknown effect) should especially be taken into account. Therefore, a survival model consisting of only established clinical and/or other known prognostic factors should serve as the benchmark (i.e. conventional model) for the upcoming modeling. The inclusion of new covariates (e.g. omics features) into a prognostic model only makes sense if the new covariates have added prognostic value over the established clinical prognostic factors (De Bin et al. 2014), i.e. the new prognostic model consisting of the new covariates plus benchmark covariates improves the prediction performance over the conventional model.
任何统计预测模型的基本目标都是实现比现有统计模型更好的预测性能，我们可以称之为 “传统模型”（Gerds 和 Kattan 2021）。应特别考虑临床和/或其他已知预后因素（通常是低维的，具有既定效果）和新的组学特征（高维，效果未知）组合的特殊性质。因此，仅由已建立的临床和/或其他已知预后因素组成的生存模型应作为即将进行的建模的基准（即 常规模型）。只有当新的协变量在已建立的临床预后因素上增加了预后价值时，将新的协变量（例如组学特征）纳入预后模型才有意义（De Bin 等人，2014 年），即由新协变量和基准协变量组成的新预后模型提高了预测性能优于传统模型。

To confirm that the identified clinical and omics features have prognostic power with respect to the prediction of patients’ survival outcomes, a model should be both accurate (low prediction error) and precise (low prediction uncertainty). The simplest way to demonstrate the prognostic power is to dichotomize the prognostic scores [i.e. linear predictor $Xiβ$ in the Cox model, Cox (1957)] by its median value, and then use a log-rank test to compare the survival probabilities of the patients in the two groups, see Fig. 2b. Similarly, one can categorize the prognostic scores by multiple quantiles (e.g. 25%, 50% median, and 75%) into multiple groups of patients and perform a log-rank test.
为了确认已确定的临床和组学特征在预测患者生存结果方面具有预后能力，模型应同时准确（低预测误差）和精确（低预测不确定性）。证明预后能力的最简单方法是将预后评分 [即 Cox 模型中的线性预测因子 $Xiβ，Cox$ （1957）] 按其中位数，然后使用对数秩检验来比较两组中患者的生存概率，见图 2b.同样，可以按多个分位数（例如 25%、50% 中位数和 75%）将预后评分分类为多组患者，并执行对数秩检验。