Multi-Parameter Quantitative Imaging of Tumor Microenvironments Reveals Perivascular Immune Niches Associated With Anti-Tumor Immunity

Tumors are populated by a multitude of immune cell types with varied phenotypic and functional properties, which can either promote or inhibit anti-tumor responses. Appropriate localization and function of these cells within tumors is critical for protective immunity, with CD8 T cell infiltration being a biomarker of disease outcome and therapeutic efficacy. Recent multiplexed imaging approaches have revealed highly complex patterns of localization for these immune cell subsets and the generation of distinct tumor microenvironments (TMEs), which can vary among cancer types, individuals, and within individual tumors. While it is recognized that TMEs play a pivotal role in disease progression, a better understanding of their composition, organization, and heterogeneity, as well as how distinct TMEs are reshaped with immunotherapy, is necessary. Here, we performed spatial analysis using multi-parameter confocal imaging, histocytometry, and CytoMAP to study the microanatomical organization of immune cells in two widely used preclinical cancer models, the MC38 colorectal and KPC pancreatic murine tumors engineered to express human carcinoembryonic antigen (CEA). Immune responses were examined in either unperturbed tumors or after immunotherapy with a CEA T cell bispecific (CEA-TCB) surrogate antibody and anti-PD-L1 treatment. CEA-TCB mono and combination immunotherapy markedly enhanced intra-tumoral cellularity of CD8 T cells, dominantly driven by the expansion of TCF1-PD1+ effector T cells and with more minor increases in TCF1+PD1+ resource CD8 T cells. The majority of infiltrating T cells, particularly resource CD8 T cells, were colocalized with dendritic cells (DCs) or activated MHCII+ macrophages, but largely avoided the deeper tumor nest regions composed of cancer cells and non-activated macrophages. These myeloid cell – T cell aggregates were found in close proximity to tumor blood vessels, generating perivascular immune niches. This perivascular TME was present in untreated samples and markedly increased after CEA-TCB therapy, with its relative abundance positively associated with response to therapy. Together, these studies demonstrate the utility of advanced spatial analysis in cancer research by revealing that blood vessels are key organizational hubs of innate and adaptive immune cells within tumors, and suggesting the likely relevance of the perivascular immune TME in disease outcome.

the corresponding paired .wsp and .fcs files after saving them in FlowJo. Once imported the following functions were used in CytoMAP to analyze the data.
Cluster CEA and PD-L1 Spots CEA spot objects and PD-L1 spot objects were clustered using a self-organizing map. The mean fluorescent intensity of all 13 channels, standardized per sample, were used as the clustering parameters. Clusters were broadly grouped into subcategories based on their prevalence, channel intensity, and location in the image.

Make surface
To find the distance to the tumor border a surface was created around regions found by clustering neighborhoods using only the CEA and PD-L1 spot sub-types. This surface was created in CytoMAP using the Make Surface function. This function uses MATLAB's alphaShape function, with user defined parameters, to wrap the currently plotted points in a surface. This surface defines the CEA + Tumor region border shown in Figure 2B.

Calculate distance
This function was used to find the physical distance to the tumor border as well as various cell types shown in Figures 2C and S5.

Raster scan neighborhoods
This function was used to define neighborhoods for subsequent analysis shown in Figures 2, 3, 4, and 6. For the analysis presented here neighborhood analysis calculated the number of cells within a cylindrical volume in the tissue with a radius of 50 µm. The positions of the neighborhoods are evenly distributed throughout the tissue in a grid pattern with a distance between neighborhood centers of 25 µm.

Classify neighborhoods into regions
This function was used to define tissue regions. For all figures in this manuscript, we used the Global Composition (number of cells divided by the maximum number of cells in all neighborhoods across all samples). In this manuscript, the physical position of the neighborhoods was not used for region definition and the minimum of the Davies-Bouldin function was used to automatically determine the number of regions. If the global minimum of the Davies-Bouldin function was 2, then the next lowest minimum was chosen. The neighborhoods were clustered using the SOM function. For visualization purposes the names and colors of the regions were changed with the Annotate Regions function in CytoMAP. Some regions were combined for clarity. The region heatmaps and prevalence of the color-coded neighborhoods was plotted using the Region Statistics function in CytoMAP. The region heatmap plots the fold change in the number of objects per neighborhood within the indicated regions compared to the number of objects from all neighborhoods. The spatial distribution of the regions was visualized by generating a new figure in CytoMAP, plotting the positions of the neighborhoods, and selecting the regions for the 'c' axis to color-code the neighborhoods by region type.
Reduce dimensions UMAP and t-SNE dimensionality reduction algorithms were used within CytoMAP to visualize tissue structure and complexity, as well as for treatment group comparisons. The Global composition of the neighborhoods was used as the dimensions to be reduced. We used the MATLAB implementation of UMAP, provided by the Herzenberg Lab at Stanford University available for download at: https://www.mathworks.com/matlabcentral/fileexchange/71902uniform-manifold-approximation-and-projection-umap.   (ch11-ch10).*(ch11>ch10) All Myeloid Channel V2 Cd11c + MHC-II + CD103 + Sirpa + CD206