Development and Validation of a Machine Learning Model for Detection and Classification of Tertiary Lymphoid Structures in Gastrointestinal Cancers

This diagnostic/prognostic study investigates the use of a machine learning model in evaluating tertiary lymphoid structures and their association with survival in gastrointestinal cancers.


Image preprocessing
All H&E images were reviewed to ensure sufficient image quality.Whenever available, images at the 40× magnification were processed and analyzed; 65 slides were scanned at 20× and the corresponding images were used.To minimize the influence of image artifacts, we used the Openslide software to down-sample the whole-slide images by a factor of 32, and then removed those regions with pen marks, folding and blurring artifacts by using appropriate color filters (https://github.com/histolab/histolab).

Tumor detection
Since only TLS within or around the tumor area are relevant, we first trained a deep learning model for automated tumor detection in histopathology images.For this purpose, we used publicly available and previously annotated H&E-stained whole-slide images of colorectal cancer and stomach cancer as the training dataset 1 (https://doi.org/10.5281/zenodo.2530789).A total of 94 whole-slide images from 81 patients were used to create 11977 image tiles of 512×512 at 0.5 µm/pixel.Each image tile was manually annotated as tumor and non-tumor (including adipose tissue, mucus, stroma or muscle).We used the ResNet18 deep learning model and adopted the same experimental setting described previously 2 .During the training process, we used horizontal/vertical flipping and translation to augment the training dataset.
We employed a cross-entropy loss function and used the Adam optimizer with a learning rate of 5×10 -6 for training, and counteracted overfitting by an L2-regularization of 1×10 -4 .The batch size was 64 and training was run for 25 epochs.Finally, tumor segmentation was expanded via image dilation by 0.5 mm to include the invasive margin.

Nuclei segmentation and lymphocyte classification
For single-cell analysis, we trained a Mask R-CNN deep learning model to segment and classify individual nuclei into tumor cells, lymphocytes, and other nonmalignant cells.For this purpose, we used our previously curated public dataset set, which contains a total of manually annotated 17,582 tumor cells, 22,550 lymphocytes, and 10,675 other non-malignant cells in 1358 image patches from 66 patients in the TCGA-LIHC dataset (https://github.com/zilanjiuwan/Single-Cell-Imaging-Analysis-of-HCC-Data).We adopted the same experimental setting described in our previous study 3 .Specifically, we used horizontal/vertical flipping, scaling, rotation, contrast normalization, affine transformation, and Gaussian blurring to augment the training dataset.We employed a stochastic gradient decent with momentum of 0.9, a weight decay of 1×10 -4 , and a batch size of 4. The network was trained for 20,000 iterations, starting from a learning rate of 0.001, and reducing to 0.0002 at 16,000 and 0.0001 at 18,000 iterations.The Mask R-CNN model was trained in TensorFlow and Keras platform and trained using the NVIDIA Tesla V100 (32GB).
In addition to TLS scoring, we also computed the density of tumor-infiltrating lymphocytes (excluding TLS) per tumor area based on the results of nuclei segmentation and classification.

Molecular evaluation of the imaging-based TLS scores
We used the matched histopathology image and gene expression data in TCGA cohorts to explore the underlying molecular features associated with the TLS scores computed on images.Since TLS consists of multiple cell types including B and T lymphocytes, we assessed the correlation between TLS scores and the abundance of immune cell infiltrate estimated from gene expression data.In addition, because cytokines play a critical role in mediating the formation and maturation of TLS, we also assessed the relations between TLS scores and cytokine gene expression levels.We further developed a gene signature for the imagingbased TLS score using the TCGA-STAD cohort as the training set.Finally, we evaluated the prognostic effect of the TLS gene signature in gastric and colorectal cancers given their high incidence rates among GI cancers.These cohorts include TCGA-STAD, TCGA-COAD/READ, and 6 additional largest datasets for which gene expression and survival data are publicly available (gastric cancer: GSE62254, GSE84437 and GSE15459; colorectal cancer: GSE39582, GSE14333 and GSE37892).

Genes associated with the TLS scores
We used the gene expression data and TLS scores available for 335 patients in the TCGA-STAD cohort for this analysis.The absolute abundance of 10 different cell types in the TME were estimated from bulk gene expression data using the MCPcounter algorithm.The cytokine genes in the Cytokine-cytokine receptor interaction pathway of the Kyoto Encyclopedia of Genes and Genomes (KEGG) database were included in the following analysis.One-hundred and four cytokine genes were measured in all the patients.All the features were centered and scaled.Due to high percentage of patients with no detectable TLS, we applied a univariate tobit model (R package VGAM), which was specifically designed to model zero-inflated data to assess the correlations between cytokines or immune cell abundance with the imaging-based TLS score.

Development of a gene expression signature of TLS score
We used the TCGA-STAD cohort as the training set given that it had the highest levels of TLS.
Considering that different types of TLS (1 vs. 2 vs. 3) are biologically related and may share common molecular regulatory mechanisms, instead of modeling each TLS separately, we employed a multi-task learning scheme to conduct feature selection and model construction to predict the 3 TLS scores simultaneously.First, we used the univariate tobit model to assess the correlations between cytokines with each of the 3 TLS scores.The resulting P values for each gene was adjusted using Benjamini & Hochberg method.We selected cytokines that were significantly associated with two or all the TLS scores (FDR < 0.1) as the candidate genes.Based on these features, we trained a multi-task linear regression model in the TCGA-STAD cohort (R package, RMTL), with 'L21' norm and penalty strength determined by 10-fold cross validation.In the final model, 11 cytokines were included with non-zero weights associated with each of the 3 TLS scores.The gene signature for the overall TLS score was calculated based on the weights determined previously.

Validation of the TLS gene signature in independent cohorts
The prognostic effect of gene signature of TLS score was assessed in 8 large independent cohorts of gastric or colorectal cancers.Overall survival information is available for all the gastric cohorts and was used to assess the prognostic effect of TLS gene signature.Recurrence free survival is available for all the colorectal cohorts and its relationship with gene signature was assessed.The combined prognostic effect in multiple cohorts was estimated using the generic inverse variance method (R package meta).
Clinicopathologic variables such as age, gender, tumor stage, and the estimated abundance of CTLs were included in the multivariate Cox regression analysis.

eResults Different weighting of individual TLS scores
We then investigated different weighting of individual TLS scores.Instead of using TCGA-STAD as the training set, we retrained a linear model using the SMU-STAD, TCGA-PAAD, or combined dataset.The overall TLS scores remained highly stable with Pearson correlation >0.93 and prognostic patterns were similar to the original results (eFig.22).

Molecular correlates and gene signature of TLS score
We investigated the molecular features related to TLS and set out to develop a gene signature for TLS scores.Among the 10 immune and stromal cell types estimated from gene expression data, score for TLS 1 had the strongest correlations with T cell abundance including CD8 T cells, while score for TLS 2 was most correlated with B cell abundance.The score for TLS 3 had strong positive correlations with both T and B cells, consistent with the enrichment of both cell types in mature TLS (eFig.23).All TLS scores had negative correlations with neutrophils and fibroblasts, which often play an immunosuppressive role in the tumor microenvironment.Among 104 cytokine genes, we identified 23 cytokines that were significantly correlated with the imaging-based TLS score (eFig.24).Considering that the TCGA-STAD cohort had the highest level of TLS compared with other cancer types, we used this cohort as the training set to obtain a gene signature for TLS score (eTable 11).The final model consists of 11 cytokine genes: 5 genes including CXCL13, CXCL11, CXCL10, CD40LG, LTA had positive weights and were upregulated in tumors with a high TLS score; while 6 genes including INHBB, INHBA, TGFB2, VEGFB, PDGFA, CLCF1 had negative weights and were downregulated in tumors with a high TLS score (eFig.25A).There was significant association (P < 0.0001 and P = 0.00138) between the 11-gene signature and TLS score in TCGA-STAD and TCGA-COAD/READ cohorts.
We finally assessed the prognostic effect of the 11-gene signature of TLS in two most common GI cancers, i.e., gastric and colorectal cancers (eFig.25B-C).We observed a consistently favorable prognostic effect of the TLS gene signature in 4 independent cohorts of gastric cancer, with overall HR = 0.75 (95% CI: 0.7 to 0.81, P < 0.001).Similarly, the 11-gene signature also had a favorable prognostic effect in 4 independent cohorts of colorectal cancer with overall HR = 0.68 (95% CI: 0.61 to 0.75, P < 0.001).In multivariable analysis, the TLS gene signature remained an independent prognostic factor when adjusting for clinicopathologic factors and abundance of cytotoxic lymphocytes in both gastric cancer (HR = 0.74, 95% CI: 0.68−0.80,P < 0.001) and colorectal cancer (HR = 0.78, 95% CI: 0.70−0.87,P < 0.001), as shown in

eFigures. Supplementary Figures
Proposed Workflow for Automated Tertiary Lymphoid Structure Evaluation on Hematoxylin-Eosin-Stained Whole-Slide Images Computational pipeline for artificial intelligence-based detection, classification, and quantitative evaluation of TLS.Our computational pipeline uses whole-slide images (WSIs) as input and consists of three modules: (1) a ResNet18 model for segmenting tumor areas in WSIs, (2) a Mask RCNN model for segmenting lymphocytes in tumor areas, and (3) the classification and regression trees (CART) algorithm for classifying individual TLSs.Flow Chart of Patient Inclusion and Exclusion 378 WSIs of 353 patients available 4 WSIs of 4 patients with poor image quality 382 WSIs of 357 patients available 18 WSIs of 18 patients without survival information 400 WSIs of 375 patients available 42 WSIs of 41 patients without image resolution (µm per pixel) 442 H&E-stain WSIs of 416 patients collected from TCGA-STAD 415 WSIs of 409 patients available 5 WSIs of 3 patients with poor image quality 420 WSIs of 412 patients available 22 WSIs of 22 patients without survival information 442 WSIs of 434 patients available 17 WSIs of 17 patients without image resolution (µm per pixel) 459 H&E-stain WSIs of 451 patients collected from TCGA-COAD 208 WSIs of 182 patients available 1 WSIs of 1 patients without survival information 209 H&E-stain WSIs of 183 patients collected from TCGA-PAAD 332 WSIs of 332 patients available 110 WSIs of 110 patients without survival information 442 H&E-stain WSIs of 442 patients collected from SMU-STAD 7 WSIs of 7 patients with poor image quality 146 WSIs of 145 patients available 3 WSIs of 3 patients with poor image quality 149 WSIs of 148 patients available 9 WSIs of 9 patients without survival information 158 WSIs of 157 patients available 8 WSIs of 8 patients without image resolution (µm per pixel) 166 H&E-stain WSIs of 165 patients collected from TCGA-READ 363 WSIs of 355 patients available 1 WSIs of 1 patients with poor image quality 364 WSIs of 356 patients available 8 WSIs of 7 patients without survival information 372 WSIs of 363 patients available 7 WSIs of 2 patients without image resolution (µm per pixel) 379 H&E-stain WSIs of 365 patients collected from TCGA-LIHC 157 WSIs of 155 patients available 1 WSIs of 1 patients without survival information 158 H&E-stain WSIs of 156 patients collected from TCGA-ESCA % eFigure 3. Confusion Matrices for Nuclei Classification on Training and Testing Data Set Values are the percentage and number of nuclei correctly and incorrectly classified by the Mask R-CNN model eFigure 4. Example Images for Tertiary Lymphoid Structure Segmentation and Classification Example image patches (A) with cell segmentation and classification (B), lymphocyte density maps (C) and TLS segmentation and classification (D).Red, blue, and green outlines correspond to TLS1, TLS2 and TLS3, respectively.Scale bar: 100 .Distributions of Tertiary Lymphoid Structure Across 6 Cancer Types in 7 CohortsThe percentage of tumors in which TLS1-3 were detected (A) and any TLS were detected (B).
Gene expression surrogate in CRC: HR (95% CI, p−value) Association Between Tertiary Lymphoid Structure Score and Tumor Stage or Grade in 7 Cohorts Univariate and Multivariate Survival Analysis of Individual Tertiary Lymphoid Structure Scores in 7 Cohorts and Combined Data Set Univariate and Multivariate Survival Analysis of Overall Tertiary Lymphoid Structure Scores in The Cancer Genome Atlas Esophageal Carcinoma Univariate and Multivariate Survival Analysis of Overall Tertiary Lymphoid Structure Scores in The Cancer Genome Atlas Stomach Adenocarcinoma Univariate and Multivariate Survival Analysis of Overall Tertiary Lymphoid Structure Scores in Southern Medical University Stomach Adenocarcinoma Univariate and Multivariate Survival Analysis of Overall Tertiary Lymphoid Structure Scores in The Cancer Genome Atlas Colon Adenocarcinoma Univariate and Multivariate Survival Analysis of Overall Tertiary Lymphoid Structure Scores in The Cancer Genome Atlas Rectum Adenocarcinoma Univariate and Multivariate Survival Analysis of Overall Tertiary Lymphoid Structure Scores in The Cancer Genome Atlas Pancreatic Adenocarcinoma The 11 Cytokine Genes and Corresponding Weights in Tertiary Lymphoid Structure Molecular Signature