Feasibility of Inferring Spatial Transcriptomics from Single-Cell Histological Patterns for Studying Colon Cancer Tumor Heterogeneity

Overview of Cells2RNA Framework: Bridging Histological Patterns with Single-Cell Expression

Cells2RNA was crafted to infer single-cell expression insights from discernible histological patterns in instances where spatial transcriptomics and single-cell data might be lacking (Figure 1). The crux of the challenge lies in deducing single-cell nuances solely from histological patterns surrounding pinpointed cells (Figure 2A). Historically, prior research has been limited to interpreting aggregated spot-level data. Yet, when disaggregated to the individual cell level, a richer tapestry of heterogeneity emerges, which becomes our focal point for inference. This innovation holds the potential to augment findings across expansive cohorts without the constraints of traditional spatial analyses that heavily rely on Visium data. The goal of this study is to derive molecular insights paralleling the depth of Visium-based investigations, but strictly from histological imaging.

• Download figure
• Open in new tab

Figure 1: Overview of Cell2RNA’s Co-Mapping Approach

(A) Low-dimensional visualization of single-cell RNA profiles, clusters indicating cell-type. (B) Spatial layout of identified cells across the tissue slide (assignment to spots represented by hexagons), color-coded by distinct gene expression patterns mapped from single cell profiles featured in (A). (C) In-depth view of cells located within a specific Visium spot, illustrating connectivity and cell relationships. Expression-related histological features, represented by grey rectangles, are shared among neighboring cells through red curves via a graph neural network. (D) A side-by-side low-dimensional comparison of scRNASeq profiles for a representative Visium spot: actual expression (red), model-predicted expression using the co-mapping training approach (green), and expression prediction without co-mapping training (blue).

• Download figure
• Open in new tab

Figure 2: Schematic Representation of the Neural Network Workflow for Single-Cell Analysis

During the training phase, (A) a pre-trained Mask R-CNN model is applied to raw histology images to detect individual cells, after which a 6-nearest neighbors graph is constructed for the detected cells. (B) Features for each cell are extracted using a ResNet-50 neural network, and the aggregation of neighboring cell information is modeled using a Graph Attention Network (GAT). For each Visium spot, the node features are aggregated using sum pooling. (C) Pre-pooled node values are jointly optimized against single-cell RNA-sequencing (scRNA-seq) data, and (D) pooled Visium spot predictions are optimized against the corresponding ground truth data, using a mean-squared-error loss computed across log-transformed counts. (E)-(G) Visual description of neural network architectures and penalizations employed: (E) a two-stage neural network comprising a feature extraction stage and a prediction stage, this was not used in this work, (F) an end-to-end neural network encompassing the entire process from cell detection to feature extraction, graph convolutions and prediction, utilized in this study, and (G) the incorporation of single-cell-level penalties into the loss function to enforce consistent predictions with scRNA-seq data.

Central to our approach is a co-mapping methodology. Here, histological patterns detected at the cellular level are intricately aligned with single-cell expression data (Figure 1A). Spatial transcriptomics serves as a vital intermediary in this process: during training, it not only maps single-cell RNASeq data to corresponding Visium spots (Figure 1B) where cells are located but also acts as a crucial inference target for the expression-centric histological attributes derived from these located cells.

To accomplish this, Cell Graph Neural Networks (CGNNs) are employed (Figure 2B,E,F). They efficiently extract the essence of histological features from located cells within Visium spots (Figure 1C). These features are updated during model training through simultaneous harmonization with scRNASeq data via optimal transport methodologies (Figures 1D,2G), forging a cohesive link between cellular histological characteristics and their intricate gene expressions.

The culmination of this methodological synergy is a framework where the expression profiles predicted from detected cells align with the broader single-cell expression landscape. Although this alignment might not be perfect, it closely mirrors genuine single-cell expression dynamics within each Visium spot. Thus, our anticipation is twofold: the accurate convergence of aggregated inferred single-cell expression to spot-level benchmarks and a holistic reflection of single-cell expression diversity within the spot when applied to external, held-out slides (Figure 2C,D).

Using Visium and paired 40X resolution whole slide imaging from a cohort of nine stage pT3 colorectal patients (see section “Data Collection and Preprocessing”), the co-mapping technique was benchmarked against patch-level models (Inceptionv3) and other CGNNs that utilize alternative information extraction methods. We assessed their performance based on predicting spot-level expression, capturing cellular heterogeneity within spots (using Wasserstein distance), maintaining tissue architectural relationships, and conducting a pathway analysis. A detailed examination of these methodologies and comparisons is provided in the Methods section.

Model Comparison

Overall, models have strong performance– selecting the top CGNN model per gene resulted in an AUROC of 0.8138 ± 0.0069 and Spearman’s statistic of 0.5724 ± 0.0133 (Table 1). However, across all experiments, model performances did not appear significantly different from each other, though we noticed several important trends (Figure 3,4). CGNN models performed on par with the Inception model (AUC=0.8204 ± 0.0073). The most predictive cell-based model had an AUROC of 0.8093 ± 0.0083, similar to the InceptionV3 model’s AUROC interval of 0.8204 ± 0.0073, which leveraged additional information beyond the cell’s immediate neighborhood and may have also benefited from the built-in structural feature extraction of CNNs. There was high agreement in top-performing genes between CGNN methods using graph contrastive learning or single-cell penalization as compared to a CGNN with no penalization/pretraining (Supplementary Figure 1, Supplementary Table 1). Notable genes with high predictive performance include TMSB4X, which encodes an actin-sequestering protein vital for actin polymerization, cell proliferation, migration, differentiation, and bypasses X inactivation with a homolog on chromosome Y; and ELF3, which regulates the inflammatory response (Supplementary Table 1).

• Download figure
• Open in new tab

Figure 3: Performance comparison between methods.

A) Boxplot of AUROC scores from each method; B) comparison of AUROC for best CGNN and CNN for each gene

• Download figure
• Open in new tab

Figure 4: Predicted expression for various genes

CNN, CGNN, compared to ground truth for genes CDX1, COL3A1, CDH1 and EPCAM across sections from all nine patients

Single-Cell Attribution Maps Point to Spatial Cellular Heterogeneity

Single-cell regularization was able to significantly improve the alignment of cellular information extracted from located cells within held-out slides with their corresponding single-cell profiles, as measured by the Earth Mover’s (Wasserstein) distance between cells assigned to spots using Tangram and their closestdetected matches (EMD=0.1415 ± 0.0018 with penalization, 0.2113 ± 0.0018 without penalization). This improvement does not negatively impact AUROC. Cells were embedded using UMAP based on the ground truth and predicted expression, with and without penalization with scRNASeq. Visual inspection of these UMAP embeddings confirmed the quantitative results of differences in EMD (Supplementary Figure 2), that single-cell penalization causes node-level predicted expression from cellular histomorphology for genes to more closely resemble the distribution of single-cell data assigned to the Visium spot with Tangram.

An optimal transport approach matched single-cell profiles between the predicted and ground truth. We then assessed the correlation between individual cells for specific genes based on their performance on the Visium spot level. Overall, more than 80% of the genes exhibited a positive correlation between ground truth and predicted single-cell expression when single-cell regularization was employed, compared to around 20-30% of the genes without such regularization was not used (Supplementary Figure 3). Examining inferred cell-level expression before the final pooling layer makes the source of our network’s predicted gene expression evident. As illustrated in Figure 5E-G, we juxtapose the predicted level of EPCAM expression for each cell against ground truth data from a Visium assay. Notably, the Visium assay offers aggregate expression measurements spanning a broader region. Our model’s predictions and the ground truth at cellular resolution are visually consistent (Figure 5A-D), further corroborating with both the high accuracy results reported from the previous section as well as the lower EMD reported through single-cell penalization.

• Download figure
• Open in new tab

Figure 5: Alignment of True and Predicted Single-Cell and Visium-Spot Level Expression on a Histological Section.

The figure illustrates the relationship between true and predicted single-cell expression on a histological section for genes CD24, KRT8, and S100A6. Ground truth expression was assigned to each cell using optimal transport based on comparing predicted and actual expressions rather than their spatial arrangement. A) and C) display the ground truth of single-cell expression with and without single-cell regularization, respectively. B) and D) visualize the respective predicted single-cell expressions. Progressing from individual cellular predictions to a broader view, D)-G) detail the transition through EPCAM expression: from predicted cell-level expression in D) to an overarching cell graph across multiple Visium spots in E) and concluding with spot-level Visium expression in G).

Cell-level attribution (i.e., inferred expression) maps mirror tissue architecture. A closer inspection of a tissue section reveals patterns, such as cells exhibiting high expression juxtaposed with a cluster of low-expression cells. Such patterns potentially highlight activities conforming to the underlying tissue structure. The structural coherence of granular cell-level information likely stems from our regularization strategies, compelling the neural networks to internalize a general representation of cell graphs and the sheer volume of available data. In instances where Visium spots exhibit diverse expression levels yet have similar cell graphs, the graph neural network may pinpoint certain cells that act as primary influencers or representatives of the expression within that spot (i.e., have markedly similar/different expression for specific genes). These specific cells could either epitomize the overall expression profile or significantly sway the expression heterogeneity at that location. Further validation is necessary to corroborate these findings, and this will be a central focus of our future research.

Topological Consistency of Inferred Expression Patterns

Across all capture areas, predicted spot level expression clustered similarly to the true expression as the relative positioning between the true and predicted clusters in the UMAP plots was preserved (Figure 6). However, there were notable differences. Overlaying the clusters assigned to ground truth embeddings over the predicted expression embeddings, we found that clusters were less separated, more connected, clumped, and fuzzier than the ground truth. Nonetheless, overlaying cluster assignments across the whole slide image demonstrates the ability of these models to derive expression signatures that can delineate key histological architectures that will be the subject of inquiry in future work.

• Download figure
• Open in new tab

Figure 6: UMAP embeddings of tissue slides from selected capture areas, color-coded by HDBSCAN clusters.

Comparisons include CGNN, CGNN with single-cell penalization, and patch-based methods against the ground truth. Clusters derived from the ground truth are overlaid on the slides for context.

Pathway Analysis

To compare performance across several potential prediction targets, we selected pathways from MSigDB’s Hallmark Gene Sets [35,36] and reported the average AUC for genes from these sets. Across all modeling approaches, genes involved in DNA repair and E2F targets were predicted with higher performance as compared to other molecular pathways (Supplementary Figure 4). Dysregulation of DNA repair can impact the acceleration of tumor progression [37], and therefore accurately detecting the presence of relevant genes may be useful in the prognostication of the speed of progression. We did notice that for some pathways, e.g., Epithelial to Mesenchymal Transition, penalizing by single-cell expression did lead to some loss of performance in distinguishing these molecular signatures (Supplementary Figure 4).

We performed a pathway analysis by subsetting the top 10% of genes (100 out of 1000) per modeling approach for further analysis using the Enrichr software/database [38–41]. The top ten statistically-significant pathways (from the Elsevier pathway collection), divided by each modeling approach, are shown in Supplementary Table 2. We also found that the WNT in Epithelial to Mesenchymal Transition in Cancer pathway, a chief contributor to the migration and metastasis of cancer cells, and several pathways associated with desmosome assembly (which regulate intercellular adhesion between metastasizing cells) were among the top ten most statistically significant gene sets detected in all four techniques, and EPCAM in Cancer Cell Motility and Proliferation is a statistically significant gene set in all four techniques. The WNT in Epithelial to Mesenchymal Transition in Cancer pathway has an AUROC of 0.8686 ± 0.0273 for the Inception model and 0.8638 ± 0.0238 for the “vanilla” cell graph model.

Comparison of Cell-Level Approach to Patch-Based Methods

The performance of the CNN model does not surpass that of the cell-based approaches. Interestingly, our basic cell model, devoid of any augmentations or pretraining, demonstrates a bootstrapped AUROC confidence interval overlapping with that of the Inception model. This indicates that, even when operating with potentially less diverse information like the extracellular matrix and connective tissue, the cell-based model remains competitive against its CNN counterpart. We posit that patch-based CNNs have an inherent advantage due to their richer input dataset and innate spatial reasoning granted by 2D convolutional kernels. In contrast, graph neural networks (GNNs) rely exclusively on adjacency matrices for spatial information. The crux here is the trade-off between robustness and interpretability. Although CNN may show a slight performance advantage, its insights are limited to single-pixel attributions, neglecting the broader scope of cell-cell interactions. Conversely, the GNN model offers superior explainability, permitting direct visualization of pivotal cell-cell interactions for particular genes, with tools like Captum, GNNExplainer, and topological methods that can be used to decipher important structural motifs.

Impact of Single-Cell Penalization

Several of our experiments, including single-cell penalization and contrastive pretraining, showed minimal influence on the final outcome. This lack of influence indicates that employing single-cell penalization can shed light on the spatial nuances of cellular disparities without compromising performance, emphasizing the significance of cellular phenomena and connectivity. We believe this is due to the large dataset size (more than 60,000 Visium spots), which may mitigate the need or potential benefit of pretraining. Additionally, although we hoped that single-cell penalization would improve the model’s robustness (by grounding predictions in real single-cell RNA quantification), the penalization provided modest performance gains over methods that did not employ this penalization. This modest gain suggests that models may produce the same optimum regardless of the intermediate feature values (i.e., cell-level predictions). However, we do note that single-cell penalization causes cell-level values to be more aligned with true single-cell data, which may indicate potential future applications in deconvolution at single-cell resolution pending further validation. Cell regularization seems to have a slight negative impact across all pathways, but this might be expected since it more directly impacts node predictions rather than graph predictions. Without this regularization, the model is more directly trained towards Visium alignment, on which it is validated. However, with the regularization, we notice that agreement between node-level predictions and single-cell profiles improves (as measured by a decrease in mean-squared error). Therefore, we propose that this method would be useful in exploring the potential for graph neural networks in cell deconvolution.

Revisiting Topological Consistency and Intermediate Histologically-Associated Molecular States

We discovered that although the predicted expression patterns mirrored the essential topological relationships tied to specific histological structures, they were more intertwined compared to the true expression, resulting in less pronounced clustering. Such mixed clustering might suggest that these clusters signify different degrees of cellular activity for various phenomena. It seems easier for machine learning models to distinguish between low and high activity levels, but interpolating intermediate levels of activity poses a challenge from a visual standpoint. In some instances, while the ground truth UMAP expression plot showed nuanced expression, the predicted expression indicated centralized expression around specific profiles. This implies a potential loss of intricate genetic interrelations, which is expected when deducing genomic data mainly from visual information in histology. This nuanced loss in the distinctiveness of the models’ topological embeddings could indicate the limitation of extracting molecular data from histology alone, suggesting richer molecular details exist beyond just histology. A potential enhancement might involve using single-cell expression profiles to better guide the model toward a more accurate gene distribution representation. Nevertheless, overall, the model’s predictions are topologically in line with the ground truth. Areas of tissue with similar ground truth measurements also exhibit similar predicted expressions.

Reflections on Pathway Analysis

The WNT in Epithelial to Mesenchymal Transition in Cancer and EPCAM in Cancer Cell Motility and Proliferation were notable pathways from the results section. Wnt/β-catenin signaling is implicated in cell differentiation and proliferation and has been implicated in increasing the number of “stem-like” cells in a tumor [42]. EPCAM is responsible for modulating epithelial cell adhesion, and - while having conflicting trends in recent research - can result in adhesive and migratory cell activity, potentially impacting the potential for metastasis [43].

Immunological Considerations

Our approach to unveil single-cell heterogeneity from whole slide images through alignment with single-cell expression bears several important immunological implications. First, the spatial arrangement of immune cells not only influences processes governing the anti-tumoral response but may offer insights as to the efficacy of immunotherapies including checkpoint inhibitors which has been a timely subject of inquiry [44,45]. Deciphering the spatial make-up may also further reveal how tumors can establish immunosuppresive environments or contribute to an immune exhaustion phenotype [46–48]. These topics underscore work being done to study how tumors can alter their immunogenicity and immune evasion tactics, potentially informing CAR T-cell therapies or selection of specific antibodies which can be applied in a personalized manner [49–51]. Revealing additional heterogeneity may refine the selection of adjuvant therapy choices outside of existing prognostic measures (e.g., pTNM staging) depending on the kind of immune reactivity.

Limitations and Future Directions

Our study, while promising, has several limitations that offer avenues for future research. First, the generalizability of our detection neural network models may be constrained due to the limited sample size of the Lizard dataset used to train the cell detection models, which may not be representative of the broader colorectal tumor population. Inaccurate cell detection may have hampered the predictive performance of the cell-graph neural networks. Additionally, while our cohort of nine samples is large for a Visium study, given its cost, plans are underway to amass a larger, more diverse cohort to bolster the robustness of our findings by accounting for further tumor heterogeneity. As our cohort was restricted to pT3 patients, future work will examine the predictiveness of these algorithms at additional tumor sites and levels of invasiveness, particularly in relation to the tumor-immune microenvironment. Another limitation is that our cell graphs solely relied on local connectivity for information flow, lacking positional embeddings or integration with patch-level information. The slight performance advantage of convolutional neural networks (CNNs) over graph neural networks (GNNs) could be attributed to CNNs’ inherent structural benefits. Future work could explore the benefits of deeper or more extensive cell graphs and the incorporation of single-cell RNA priors for more accurate predictions. While cell regularization and similar strategies might not enhance the ultimate predictive capability for pooled Visium expression, they do bolster performance concerning optimal transport between single-cell profiles assigned to Visium spots using deconvolution methods (like Tangram). Direct model predictions can display notable correlations in gene expression levels across the cell distribution within each spot, even without an input single-cell profile. In the future, these methods could serve as adjustable enhancements or alternatives to Tangram. Inaccurate mapping of single-cell profiles to Visium spots may have also impacted the validity of the associations between single-cell inferred expression and could improve with the adoption of other single-cell spatial mapping methods. Potentially this approach could be used to map single cell profiles, where scRNASeq information has been collected, to their precise locations within slides without the use of Visium as an intermediary. Our approach could also extend to identify other spatially-resolved molecular features, such as protein interactions [52,53]. Validation of our findings using pathological examination of underlying cellular phenomena, through immunohistochemistry and other single-cell level spatial analysis platforms like Xenium/merScope is essential [54,55], as this study serves as a proof of concept. We anticipate applying our models to a larger cohort to identify metastasis and recurrence predictors at a population scale. Furthermore, the trade-off in model accuracy is mitigated by the potential for scalability, allowing for the assessment of a larger number of genes with reasonable predictive performance and statistical power. Overall, our study signifies a crucial step towards improving cancer diagnostics and prognosis by incorporating spatial transcriptomics into histological images, and future efforts will focus on refining these techniques.

Ultimately our findings could shed light on the molecular alterations occurring in these immune cell subsets, which may help identify potential immune evasion mechanisms employed by tumor cells, such as the upregulation of immune checkpoint molecules or the recruitment of immunosuppressive cell populations. By integrating this spatial and molecular information, we can potentially uncover the functional “hotspots” within the tumor microenvironment, where effective immune responses are either promoted or hindered. These functional “hotspots” could ultimately aid in designing immunotherapeutic strategies that aim to reinvigorate the anti-tumor immune response, either by modulating the activity of immune cells in these hotspots or by altering the spatial organization of immune cells within the tumor microenvironment.

Data Collection and Preprocessing

The dataset used in this study comprised nine patients with pathologic T Stage-III (pT3) colorectal cancer. Following IRB approval, these patients were selected through a retrospective review of pathology reports from 2016 to 2019. Notably, this cohort is distinct from our initial cohort of four slides, profiled using the original Visium assay without using CytAssist.

Patients were matched based on various criteria such as age, sex, tumor grade, tissue size, mismatch repair/microsatellite instability (MMR/MSI/MSS) status, and tumor site. MSI status was determined through the loss of expression of MLH1 and PMS2 as assessed through immunohistochemistry (there was no loss of expression for MSH2 and MSH6).

The dataset included both patients with concurrent tumor metastasis and patients without metastasis, with equal representation from each group. Patients were carefully selected to ensure a balanced representation across various factors, including tumor site, grade, node status, MSI status, and sex. Tissue blocks were sectioned into 5-micron thick layers. Specific regions of interest within these sections, including epithelium, tumor-invasive front, intratumoral areas, and lymphatics, were annotated by a board-certified GI pathologist. Following annotation, these regions were dissected from the tissue, and subjected to H&E staining, imaging, and Visium profiling at the Pathology Shared Resource at Dartmouth Cancer Center and Single Cell Genomics Core in the Center for Quantitative Biology.

Preprocessing and Augmentation

Deep Learning to Integrate Information from Localized Cells to Predict Spatial Gene Expression

Cell graph neural networks (CGNN) facilitate the exchange of messages between adjacent cells, enabling the exchange/incorporation of contextual information [62–67]. This approach effectively captures the relationships between different cell populations within the tissue, including tumor cells and surrounding immune and other cell subpopulations. The important difference between this approach and a more general convolutional neural network model is that cells and the relationships between them are modeled explicitly. Leveraging these relationships can enhance the predictive performance of our spatial RNA inference algorithms while providing additional information as to relevant cells for these predictions.

We implemented an end-to-end training strategy that integrates the simultaneous training of a Convolutional Neural Network (CNN) and a Graph Neural Network (GNN). The CNN is designed to extract cell-level features from histological images, while the GNN contextualizes these features by incorporating information from neighboring cells during the model-fitting process. This stands in contrast to previous cell-level neural network models, which typically employ a two-stage approach. In these two-stage models, features are initially extracted using a CNN, and these extracted features serve as input to train the GNN in a separate, subsequent step. Such a two-stage paradigm can potentially lead to suboptimal results, as the CNN’s parameters are fixed before the GNN training begins, thereby limiting the CNN’s ability to adapt based on the contextual information that the GNN could provide. Our end-to-end approach aims to harmonize the feature extraction and contextualization processes, enabling the CNN to learn cell-level features that are more effectively contextualized through iterative, integrated training with the GNN (Figure 1,2).

The backbone of the model is a four-layer graph attention network (GAT) [68,69], which uses self-attention mechanisms to update the representation of each cell with the information of its neighbors. We extract nodal attributes from detected cells using a ResNet-50 model, which is trained jointly with the graph attention layers. The Euclidean distances between the spatial locations of detected cells are used to form k-nearest-neighbor cell graphs (k=6, determined through a sensitivity analysis), connecting cells based on their spatial proximity. Patients within the same Visium capture area were divided according to pathologist-annotated segmentation masks– we ensured that cell graphs did not overlap across patients. The model updates cell-level embeddings through information sharing up to four neighbors away, facilitated by message passing between nodes across four graph attention layers, each mapping cells to 512-dimensional numerical vectors. Final node embeddings pass through a linear layer, producing a vector representing each gene’s relative pseudocount-transformed expression for each cell. Cells corresponding to the same Visium spot are then aggregated through global sum pooling to predict expression at the Visium spot. This is then compared to the pseudocount-transformed ground-truth Visium data with mean squared error.

Comparison of Cell-Graph Neural Network Regularization Strategies

In addition to evaluating the congruence between ground truth and predicted expression at the spot level, we explored the following methodological variations:

1. Vanilla Supervised Learning Objective: This baseline approach focuses solely on the supervised learning objective, serving as a reference for evaluating the potential gains from additional regularization strategies.

2. Incorporating Graph Contrastive Learning: This approach introduces a self-supervised regularization term that encourages the model to learn embeddings through the comparison of augmented viewpoints of the same cell graph / Visium spot to different cell-graphs / Visium spot. This can enhance the model’s sensitivity to spatial patterns in the data, potentially improving its predictive accuracy for spatial transcriptomics patterns.

3. Incorporating Single-Cell RNA-Seq Penalization through Optimal Transport: This strategy introduces a penalty term that encourages the model to align cell-level histological features more closely with corresponding single-cell RNA-Seq data. By leveraging optimal transport theory, this term effectively “guides” the model towards a solution where the spatial patterns inferred from histology are maximally consistent with independent single-cell RNA-Seq measurements, thereby enhancing the biological validity of the model’s predictions.

4. Combining Graph Contrastive Learning and Single-Cell Penalization: This approach synergistically combines both the graph contrastive learning and the single-cell RNA-Seq penalization strategies, aiming to leverage the benefits of both spatial context awareness and alignment with single-cell RNA-Seq data. This dual-regularization strategy is designed to promote a model that is both sensitive to spatial patterns and tightly aligned with independent molecular measurements, potentially offering a balance between spatial sensitivity and biological validity.

Graph Contrastive Learning

Using the PyGCL package, graph contrastive learning was implemented through augmentations to random cell positions in the nearest neighbor graph construction, dropping edges with a probability of 0.1, and masking out features with a probability of 0.3. Graph contrastive learning is a form of self-supervised learning that can improve the generalizability and robustness of graphs [70–72]. By intentionally adding noise to the training cell graphs and comparing these representations at different Visium spots, we aimed to improve the model’s generalizability when tested on held-out data by forcing it to learn a stable, robust, foundational representation.

Incorporating Single-Cell Expression

For two of our patients in the study, we obtained corresponding single-cell RNA-Seq data, which were selected to be representative across sex, microsatellite instability, and metastasis status. Due to budget constraints, we were unable to acquire matching single-cell data for all the slides, a limitation we aim to address in future work. Our approach involves the integration of single-cell RNA-Seq data as a form of model regularization. By encouraging the predictions derived from histological images of individual cells to align closely with the corresponding true single-cell expression profiles, we aim to enhance the interpretability of our models through more consistent and biologically meaningful cellular information. This, in turn, is expected to improve the model’s generalizability. The objective is to develop a nuanced understanding of cell type populations at specific spatial locations while maintaining fidelity to the ground-truth expression data sampled through the Visium platform. By fine-tuning the alignment between these profiles, we strive to increase the likelihood that our predictions accurately reflect the true cellular composition at each spatial location. In contrast to previous attempts to integrate single-cell information with histology, where spatial transcriptomics data were not employed to guide the spatial mapping of single-cell information, our approach uniquely capitalizes on the spatial context offered by the spatial transcriptomics data. Specifically, we align single-cell profiles with individual Visium spots, thereby avoiding spurious cell assignments that can arise when attempting to map single-cell information across entire whole-slide images without the benefit of such spatial guidance.

We initiated our analysis by mapping scRNA profiles to Visium spots using the Tangram tool [73], and we selected the top k most likely cells to be assigned to each spot, where k represents the number of detected cells in that spot. In the context of our CGNN model, which aims to predict log-expression levels for individual cells, it was critical to assess the alignment between these predictions and the assigned single-cell profiles. Originally, we employed a dynamic matching approach wherein the pool of cells assigned to each spot was paired with our detected cells based on the Euclidean distance between their predicted and observed expression profiles. This matching was formulated as a variant of the linear sum assignment problem, solved using the Hungarian algorithm. This algorithm established the optimal one-to-one correspondence between detected cells and their closest-matching assigned scRNA cells, thereby minimizing the Euclidean distance between them. We applied a mean squared error (MSE) loss across the log-transformed expression values to penalize discrepancies between these matched profiles. However, following internal comparisons, we shifted our approach to leverage the Wasserstein loss as a more effective metric for aligning our predictive single-cell expression profiles with the true expression profiles derived from scRNA data. The Wasserstein distance, a loss metric formulated on principles from the theory of optimal transport [74–77], quantifies the minimum cost required to transform one distribution into another, effectively capturing the overall distributional differences between predicted and ground truth single-cell expression profiles. This approach offers more flexible comparisons, as it can accommodate shifts in the distribution and discrepancies in the shape of the distributions without necessitating a direct one-to-one matching between predicted and observed cells.

Our internal comparisons revealed that the Wasserstein loss provided more precise and interpretable mappings between detected and observed cellular profiles, solidifying it as our chosen approach. The Wasserstein loss was computed for each Visium spot and used to compare the effectiveness of different regularization methods, measuring alignment between predicted and ground truth single-cell expression.

Comparison to Convolutional Neural Network Approaches

The CGNN approaches were compared to patch-based convolutional neural network methodologies deemed highly predictive from previous works– namely the InceptionV3 neural network. The Inceptionv3 network can capture features at various receptive fields. This model was trained on images of tissue patches encompassing multiple cells inclusive of surrounding matrixed tissue architecture. We initialize the model with ImageNet weights (with the final layer truncated) on the same dataset. We apply the same visual transformations as for the cell embeddings, (standard scaling, color jitter, random 90° rotations).

Training and Validation

CGNN models were implemented with the torch-geometric Python package [78], which extends the PyTorch machine-learning framework. We use PyGCL [72] to apply graph augmentations. CGNN were trained using the Adam optimizer [79] with a learning rate of 0.0001 on one Nvidia V100 GPU with Dartmouth Research Computing, quickly converging after two epochs. Similarly, the CNN model was trained for around 100000 iterations on a Nvidia V100 GPU.

The final performances of these models were compared using leave-one-patient-out cross-validation. Statistics are reported with the Spearman correlation coefficients. We also sought to assess the performance of predicting dichotomized gene expression (low/high), a binary classification problem– this was accomplished by dichotomizing expression according to [30], used to calculate the area under the receiver operating characteristic curve (AUROC) as another performance measure. Models were trained and then cross-validated using all capture areas except one, reserved for testing, and repeated for all cross-validation folds. Performance statistics were generated for each cross-validation fold, including Spearman’s correlation coefficients and area under the receiver operating characteristic curves (AUROCs) by gene. The results were then averaged across all folds to assess the best-performing model on a gene-specific basis. We calculated 95% confidence intervals for all performance statistics, reported using 1000 sample non-parametric bootstrapping.

Model Interpretation through Single Cell Attributions, Gene Embeddings and Pathway Analysis