Pathology-and-genomics Multimodal Transformer for Survival Outcome Prediction

Fig 1 illustrates our multimodal transformer framework. Our method includes an unsupervised multimodal data fusion pretraining and a supervised flexible-modal finetuning. From Fig 1a, in the pretraining, our unsupervised data fusion aims to capture the interaction pattern of image and genomics features. Overall, we formulate the objective of multimodal feature learning by converting image patches and tabular genomics data into group-wise embeddings, and then extracting multimodal patient-wise embeddings. More specifically, we construct group-wise representations for both image and genomics modalities. For image feature representation, we randomly divide image patches into groups; Meanwhile, for each type of genomics data, we construct groups of genes depending on their clinical relevance [22]. Next, as seen in Fig 1b and c, our approach enables three types of finetuning modal modes (i.e., multimodal, image-only, and genomics-only) towards prognostic prediction, expanding the downstream data utility from the pretrained model.

We define the group-wise genomics representation by referring to $N=8$ major functional groups obtained from [22]. Each group contains a list of well-defined molecular features related to cancer biology, including transcription factors, tumor suppression, cytokines and growth factors, cell differentiation markers, homeodomain proteins, translocated cancer genes, and protein kinases. The group-wise genomics representation is defined as $G_{n}\in\mathbb{R}^{1\times d_{g}}$, where $n\in N$, $d_{g}$ is the attribute dimension in each group which could be various. To better extract high-dimensional group-wise genomics representation, we use a Self-Normalizing Network (SNN) together with scaled exponential linear units (SeLU) and Alpha Dropout for feature extraction to generate the group-wise embedding $G_{n}\in\mathbb{R}^{1\times 256}$ for each group.

For group-wise WSIs representation, we first cropped all tissue-region image tiles from the entire WSI and extracted CNN-based (e.g., ResNet50) $d_{i}$-dimensional features for each image tile k as $h_{k}\in\mathbb{R}^{1\times d_{i}}$, where $d_{i}=1,024$, $k\in K$ and K is the number of image patches. We construct the group-wise WSIs representation by randomly splitting image tile features into N groups (i.e., the same number as genomics categories). Therefore, group-wise image representation could be defined as $I_{n}\in\mathbb{R}^{k_{n}\times 1024}$, where $n\in N$ and $k_{n}$ represents tile k in group n. Then we apply an attention-based refiner (ABR) [17], which is able to weight the feature embeddings in the group, together with a dimension deduction (e.g., fully-connected layers) to achieve the group-wise embedding. The ABR and the group-wise embedding $I_{n}\in\mathbb{R}^{1\times 256}$ are defined as:

where w,V1 and V2 are the learnable parameters.

To aggregate patient-wise multimodal feature embedding from the group-wise representations, as shown in Fig 1a, we propose a pathology-and-genomics multimodal model containing two model streams, including a pathological image and a genomics data stream. In each stream, we use the same architecture with different weights, which is updated separately in each modality stream. In the pathological image stream, the patient-wise image representation is aggregated by N group representations as $I_{p}\in\mathbb{R}^{N\times 256}$, where $p\in P$ and P is the number of patients. Similarly, the patient-wise genomics representation is aggregated as $G_{p}\in\mathbb{R}^{N\times 256}$. After generating patient-wise representation, we utilize two transformer layers [27] to extract feature embeddings for each modality as follows:

where MSA denotes Multi-head Self-attention [27] (see Appendix 1), l denotes the layer index of the transformer, and $H_{p}$ could either be $I_{p}$ or $G_{p}$. Then, we construct global attention poolings [17] as Eq.1 to adaptively compute a weighted sum of each modality feature embeddings to finally construct patient-wise embedding as $I_{embedding}^{p}\in\mathbb{R}^{1\times 256}$ and $G_{embedding}^{p}\in\mathbb{R}^{1\times 256}$ in each modality.

Due to the domain gap between image and molecular feature heterogeneity, a proper design of multimodal fusion is crucial to advance integrative analysis. In the pretraining stage, we develop an unsupervised data fusion strategy by decreasing the mean square error (MSE) loss to map images and genomics embeddings into the same space. Ideally, the image and genomics embeddings belonging to the same patient should have a higher relevance between each other. MSE measures the average squared difference between multimodal embeddings. In this way, the pretrained model is trained to map the paired image and genomics embeddings to be closer in the latent space, leading to strengthen the interaction between different modalities.

In the single modality finetuning, even if we use image-only data, the model is able to produce genomic-related image feature embedding due to the multimodal knowledge aggregation already obtained from the model pretraining. As a result, our cross-modal information aggregation relaxes the modality requirement in the finetuning stage. As shown in Fig 1b, for multimodal finetuning, we deploy a concatenation layer to obtain the fused multimodal feature representation and implement a risk classifier (FC layer) to achieve the final survival stratification (see Appendix 2). As for single-modality finetuning mode in Fig 1c, we simply feed $I_{embedding}^{p}$ or $G_{embedding}^{p}$ into risk classifier for the final prognosis prediction. During the finetuning, we update the model parameters using a log-likelihood loss for the discrete-time survival model training [6](see Appendix 2).

All image and genomics data are publicly available. We collected WSIs from The Cancer Genome Atlas Colon Adenocarcinoma (TCGA-COAD) dataset (CC-BY-3.0) [21, 8] and Rectum Adenocarcinoma (TCGA-READ) dataset (CC-BY-3.0) [20, 8], which contain 440 and 153 patients. We cropped each WSI into 512 $\times$ 512 non-overlapped patches. We also collected the corresponding tabular genomics data (e.g., mRNA sequence, copy number alteration, and methylation) with overall survival (OS) times and censorship statuses from Cbioportal [2, 14]. We removed the samples without the corresponding genomics data or ground truth of survival outcomes. Finally, we included 426 patients of TCGA-COAD and 145 patients of TCGA-READ.

We implement two types of settings that involve internal and external datasets for model pretraining and finetuning. As shown in Fig 2a, we pretrain and finetune the model on the same dataset (i.e., internal setting). We split TCGA-COAD into training (80%) and holdout testing set (20%). Then, we implement four-fold cross-validation on the training set for pretraining, finetuning, and hyperparameter-tuning. The test set is only used for evaluating the best finetuned models from each cross-validation split. For the external setting, we implement pretraining and finetuning on the different datasets, as shown in Fig 2b; we use TCGA-COAD for pretraining; Then, we only use TCGA-READ for finetuning and final evaluation. We implement a five-fold cross-validation for pretraining, and the best pretrained models are used for finetuning. We split TCGA-READ into finetuning (60%), validation (20%), and evaluation set (20%). For all experiments, we calculate the average performance on the evaluation set across the best models.

The number of epochs for pretraining and finetuning are 25, the batch size is 1, the optimizer is Adam [19], and the learning rate is 1e-4 for pretraining and 5e-5 for finetuning. We used one 32GB Tesla V100 SXM2 GPU and Pytorch. The concordance index (C-index) is used to measure the survival prediction performance. We followed the previous studies [6, 5, 7] to partition the overall survival (OS) months into four non-overlapping intervals by using the quartiles of event times of uncensored patients for discretized-survival C-index calculation (see Appendix 2). For each experiment, we reported the average C-index among three-times repeated experiments. Conceptionally, our method shares a similar idea to multiple instance learning (MIL) [9, 23]. Therefore, we include two types of baseline models, including the MIL-based models (DeepSet [30], AB-MIL [17], and TransMIL [26]) and MIL multimodal-based models (MCAT [6], PORPOISE [7]). We follow the same data split and processing, as well as the identical training hyperparameters and supervised fusion as above. Notably, there is no need for supervised finetuning for the baselines when using TCGA-COAD (Table 1), because the supervised pretraining is already applied to the training set.

In Table 1, our approach shows improved survival prediction performance on both TCGA-COAD and TCGA-READ datasets. Compared with supervised baselines, our unsupervised data fusion is able to extract the phenotype-genotype interaction features, leading to achieving a flexible finetuning for different data settings. With the multimodal pretraining and finetuning, our method outperforms state-of-the-art models by about 2% on TCGA-COAD and 4% TCGA-READ. We recognize that the combination of image and mRNA sequencing data leads to reflecting distinguishing survival outcomes. Remarkably, our model achieved positive results even using a single-modal finetuning when compared with baselines (more results in Appendix 3.1). In the meantime, on the TCGA-READ, our single-modality finetuned model achieves a better performance than multimodal finetuned baseline models (e.g., with model pretraining via image and methylation data, we have only used the image data for finetuning and achieved a C-index of 74.85%, which is about 4% higher than the best baseline models). We show that with a single-modal finetuning strategy, the model could generate meaningful embedding to combine image- and genomic-related patterns. In addition, our model reflects its efficiency on the limited finetuning data (e.g., 75 patients are used for finetuning on TCGA-READ, which are only 22% of TCGA-COAD finetuning data). In Table 1, our method could yield better performance compared with baselines on the small dataset across the combination of images and multiple types of genomics data.

We verify the model efficiency by using fewer amounts of finetuning data in finetuning. For TCGA-COAD dataset, we include 50%, 25%, and 10% of the finetuning data. For the TCGA-READ dataset, as the number of uncensored patients is limited, we use 75%, 50%, and 25% of the finetuning data to allow at least one uncensored patient to be included for finetuning. As shown in Fig 3a, by using 50% of TCGA-COAD finetuning data, our approach achieves the C-index of 64.80%, which is higher than the average performance of baselines in several modalities. Similarly, in Fig 3b, our model retains a good performance by using 50% or 75% of TCGA-READ finetuning data compared with the average of C-index across baselines (e.g., 72.32% versus 64.23%). For evaluating the effect of cross-modality information extraction in the pretraining, we kept supervised model training (i.e., the finetuning stage) while removing the unsupervised pretraining. The performance is lower 2%-10% than ours on multi- and single-modality data. For evaluating the genomics data usage, we designed two settings: (1) combining all types of genomics data and categorizing them by groups; (2) removing category information while keeping using different types of genomics data separately. Our approach outperforms the above ablation studies by 3%-7% on TCGA-READ and performs similarly on TCGA-COAD. In addition, we replaced our unsupervised loss with cosine similarity loss; our approach outperforms the setting of using cosine similarity loss by 3%-6%.

The results of this study are based on the data collected from the public TCGA Research Network: https://www.cancer.gov/tcga.