If you don't remember your password, you can reset it by entering your email address and clicking the Reset Password button. You will then receive an email that contains a secure link for resetting your password
If the address matches a valid account an email will be sent to __email__ with instructions for resetting your password
Department of Immunology, Genetics and Pathology, Uppsala University, Uppsala, SwedenCentre for Cancer Biomarkers CCBIO, Department of Clinical Medicine, University of Bergen, Bergen, Norway
Department of Immunology, Genetics and Pathology, Uppsala University, Uppsala, SwedenMolecular Oncology Group, Vall d'Hebron Institute of Oncology, Barcelona, Spain
Immune cells showed specific infiltration patterns within the tumour microenvironment.
•
Infiltration of lymphocytes and plasmacytoid dendritic cells are associated with good prognosis.
•
Distances of immune cells to each other or to tumour cells reveal an independent prognostic impact.
Abstract
Introduction
Immune cells in the tumour microenvironment are associated with prognosis and response to therapy. We aimed to comprehensively characterise the spatial immune phenotypes in the mutational and clinicopathological background of non–small cell lung cancer (NSCLC).
CD4-Eff cells, CD8-Eff cells and M1 macrophages were the most abundant immune cells invading the tumour cell compartment and indicated a patient group with a favourable prognosis in the cluster analysis. Likewise, single densities of lymphocytic subsets (CD4-Eff, CD4-Treg, CD8-Treg, B-cells and pDCs) were independently associated with longer survival. However, when these immune cells were located close to CD8-Treg cells, the favourable impact was attenuated. In the multivariable Cox regression model, including cell densities and distances, the densities of M1 and CD163 cells and distances between cells (CD8-Treg–B-cells, CD8-Eff–cancer cells and B-cells–CD4-Treg) demonstrated positive prognostic impact, whereas short M2–M1 distances were prognostically unfavourable.
Conclusion
We present a unique spatial profile of the in situ immune cell landscape in NSCLC as a publicly available data set. Cell densities and cell distances contribute independently to prognostic information on clinical outcomes, suggesting that spatial information is crucial for diagnostic use.
Checkpoint inhibitors are routinely used in different lines of lung cancer treatment and, for the first time, provide a chance of long-term survival. However, only a minority of patients develop a durable response, and accurate predictors of therapeutic benefits are still lacking. Concepts that explain these variable responses include features of tumour cells, such as programmed death-ligand 1 expression and high tumour mutational burden, or features of the immune environment with hot and cold immune phenotypes [
]. Unfortunately, the predictive value of these biomarkers is modest. Strategies aimed at introducing new immune modulatory agents, either alone or in combination, have been discouraging so far [
]. This can be explained by the selection of combinations based on a reductionist approach rather than a perception that integrates complex cellular interactions in the tumour microenvironment. Thus, a better understanding of the biological signals that contribute to anti-cancer immunity is urgently needed to optimise clinical diagnostics and suggest new treatment targets.
Limitations of in situ analytical tools have partially contributed to the reductionist perspective placing T-cells at the centre of investigations. Immunohistochemical evaluation of T-cell markers indicates that the amount and location of T-cells imply prognostic information [
Multispectral imaging for quantitative and compartment-specific immune infiltrates reveals distinct immune profiles that classify lung cancer patients.
Chemoradiotherapy efficacy is predicted by intra-tumour CD8+/FoxP3+ double positive T cell density in locally advanced N2 non-small-cell lung carcinoma.
Biomarkers associated with beneficial PD-1 checkpoint blockade in non-small cell lung cancer (NSCLC) identified using high-plex digital spatial profiling.
Tumour-infiltrating lymphocyte density is associated with favourable outcome in patients with advanced non-small cell lung cancer treated with immunotherapy.
]. It can be speculated that the impact is attributed to either specific T-cell subsets or other immune cells co-infiltrating the tumour tissue. The coordinated occurrence of immune cells is an intrinsic feature of the cancer microenvironment and is best exemplified in the inflamed phenotype of tumours [
]. The strong correlation between immune cell densities (e.g. cell numbers of CD8-, CD4- or CD20-positive cells) makes it difficult to determine which immune cells are the drivers of cancer immunity and are of the strongest clinical relevance. Indeed, many studies have found that the abundance of other immune cells, such as plasma cells [
The prognostic relevance of tumour-infiltrating plasma cells and immunoglobulin kappa C indicates an important role of the humoral immune response in non-small cell lung cancer.
Multispectral imaging for quantitative and compartment-specific immune infiltrates reveals distinct immune profiles that classify lung cancer patients.
Spatial architecture and arrangement of tumor-infiltrating lymphocytes for predicting likelihood of recurrence in early-stage non-small cell lung cancer.
]. T-cells directly infiltrating the tumour cell nests are likely of higher relevance than T-cells in the surrounding stroma or that T-cells interacting with B-cells represent a state of stronger immune reaction against cancer. These aspects have not been adequately addressed when immune profiles have been characterised in the in situ environment of cancer. This is in part due to technical limitations with regard to staining techniques, imaging systems, and computational power to spatially analyse complex tissue features.
Recent advances in multiplexed staining methods, in association with novel methods for digital image analysis, can address these limitations and open opportunities for biologically well-informed high-content profiling of tumour tissue [
Multispectral imaging for quantitative and compartment-specific immune infiltrates reveals distinct immune profiles that classify lung cancer patients.
]. These efforts have great potential but remain difficult to execute because of the necessary equipment, technical and histopathological expertise and tissue availability.
In this study, we aim to comprehensively describe the tumour microenvironment of non–small cell lung cancer (NSCLC) by quantification of 13 important immune cell subsets and their spatial relations with the use of a state-of-the-art multiplex immune fluorescence staining technique and a multispectral image analysis pipeline.
2. Material and methods
2.1 Patient material
The study population consisted of 359 NSCLC patients who were surgically treated at Uppsala University Hospital between 2006 and 2010 [
]. The patient characteristics are shown in Table 1. Formalin-fixed paraffin-embedded tumour blocks from resection specimens were used to construct a tissue microarray (TMA) based on duplicate 1 mm cores, as described previously [
]. The study was conducted in compliance with the Declaration of Helsinki and the Swedish Ethical Review Act (approved by the Ethical Review Board in Uppsala, #2012/532).
Mutation analysis has been described previously; for mutation analyses, genomic DNA was extracted from either fresh frozen or paraffin-embedded cancer tissue [
]. Targeted deep sequencing was performed with DNA from 352 of 359 patients using the Haloplex system for target amplification (Agilent Technologies, Santa Clara, USA). The analysis included all coding exons of 82 lung cancer–related genes (Supplementary Table 1). Sequencing was performed using 125 bp paired-end reads on the Illumina HiSeq 2500 platform (Illumina, San Diego, USA). The reads were mapped to the reference genome (hg19), and mutations were identified as described by La Fleur et al. [
]. The tumour mutational burden (TMB) was estimated by dividing the number of non-synonymous mutations in a sample by the size (0.47 Mb) of the sequenced genome.
Corresponding gene expression data were available for 197 patients obtained from RNA sequencing, as described in a previous study [
]. RNA was extracted from fresh frozen tissue and prepared for sequencing using the Illumina TruSeq RNA Sample Prep Kit v2 with polyA selection. The sequencing was performed based on the standard Illumina RNAseq protocol with a read length of 2 × 100 bases. The raw data together with clinical information are available on the gene expression omnibus with accession number GSE81089.
2.3 Multiplex immunofluorescence
Multiplex immunofluorescence staining was performed on 4 μm thick TMA sections that were deparaffinised, rehydrated and rinsed in distilled H2O. Three panels were used: a tumour-infiltrating lymphocyte panel (TIL), including CD4, CD8, CD20, FoxP3 and PanCK; a natural killer cell/macrophage panel (NK/MP), including CD3, NKp46, CD56, CD68, CD163 and PanCK; and an antibody presenting cell panel (APC), including CD1a, CD208, CD123, CD15, CD68 and PanCK. The staining procedure was performed in accordance with an earlier published study [
Multispectral imaging for quantitative and compartment-specific immune infiltrates reveals distinct immune profiles that classify lung cancer patients.
] based on the modified Opal Multiplex IHC assay (Akoya, Marlborough, USA). Detailed reagent references are provided in Supplementary Table 2. For each panel, specific staining protocols were established (TIL, NK/MP, and APC), including five to six cycles of microwave treatment, primary antibody incubation, and fluorophore labelling, with a final cycle of DAPI staining and slide mounting. After staining, the slides were scanned using the VectraPolaris system (Akoya) at a 2 pixels/μm resolution in multispectral mode and analysed in the inForm software, where spectral unmixing was used to generate an oligo-layer image with layers corresponding to the specific staining, 4′,6-diamidino-2-phenylindole (DAPI), and autofluorescence. The inForm software was used to define tumour and stroma compartments within each tissue core. The algorithm was trained on pathologist-annotated samples. Cell segmentation was based on DAPI nuclear staining. Representative subsets of the included markers were annotated as either positive or negative, and the inform software was then trained on these to phenotype all other cells accordingly. The intensity of each marker expression was used to calculate the thresholds for marker positivity. A pathologist reviewed each image and curated it with regard to artefacts, staining defects and necrosis. Tissue segmentation was also controlled.
Each patient was represented by one to two TMA cores (1 mm in diameter). When two cores were available, total cell number and tissue area were used to compute density scores. When considering distances, the mean intercellular distances for the separate cores were used.
Thresholds for marker intensity were calculated in the R programming environment by GeneVia Technologies (Tampere, Finland). The distribution of intensities for representative subsets of the included markers, annotated earlier by a pathologist as either positive or negative, was used to define the thresholds. Probability density distributions were estimated for each marker by smoothing the intensity values using Gaussian kernel estimation, with automatic bandwidth detection using the density function of the ‘stats’ package in R.
If the intensities of the positive and negative cells did not overlap, the mean value of the highest intensity of the negative cells and the lowest intensity of the positive cells was used as the intensity threshold. If the intensities did overlap, the intensity value that minimised the overall classification error based on the probability density distribution was used instead. For each established threshold, true/false positive/negative rates and the overall classification error were calculated and controlled. The thresholds were then applied to the raw output data and used to classify the cell phenotype (Fig. 1A and B). The cell phenotypes, their tissue coordinates (location) and the main clinical parameters connected to each case were uploaded at https://doi.org/10.5281/zenodo.6997517.
Fig. 1A multiplex-multispectral imaging pipeline. (A) The cohort was compiled on TMAs and stained using multiplex immunofluorescence (IF) with the tyramide signal amplification (TSA)-based detection method (Opal) and three different antibody panels. Slides were then scanned using the multispectral VectraPolaris system, and singleplex stains were used to calculate a deconvolution model for optimal identification of specific markers in the inForm software. Scanned images were annotated and curated by a specialist pathologist, and machine learning–based tissue and cell segmentation was performed. The data were then exported and further analysed in R. Cut-offs for positive marker expression were calculated for each individual marker, and cells were annotated based on the combination of marker expression. (B) Immune phenotyping: Immune cells were annotated based on the markers in three different antibody panels, resulting in 13 different immune phenotypes. The PanCKsingle phenotype represents tumour cells. TMA, tissue microarray.
], version 4.2.0. Correlation analyses between immune markers and TMB were performed using Spearman's rank correlation. Fisher's exact test was used to analyse the association between immune markers and mutations. Analysis of clinical parameters in relation to immune infiltrates and differences in immune cell densities between histological subtypes was conducted using the Wilcoxon signed-rank test. Survival analysis was undertaken using univariable and multivariable Cox regressions. The ‘maxstat’ package was used to determine maximally selected rank statistics for the cut-off used to dichotomise distance metrics and densities. Hierarchical cluster analysis was performed with Euclidean distance and the complete method using the ‘ComplexHeatmap’ package (version 2.9.3) for R [
] based on normalised values of immune cell densities. The stepwise multivariable Cox regression was performed using the ‘My.stepwise’ package for R with ‘entry’ and ‘stay’ significance levels at 0.15. Diagnostic tests, including the Schoenfeld individual and global test, dfbetas, deviance residuals and bootstrapping simulations were performed, finally indicating a stable model and meeting the proportional hazards assumption for all but two covariates (density metric ‘CD163’ and distance metric ‘PanCKsingle – M2’).
Distance analyses were carried out using the FNN package (version 1.1.3) for R. A kd-tree fast k-nearest neighbour searching algorithm was used to analyse the distance to the nearest neighbouring cells, and a mean value of the nearest neighbours was used as a case-specific metric in further analysis. Notably, this approach means that for a pair of two different cell phenotypes, the mean distance is not commutative.
P-values <0.05 were considered significant, and when performed, the Benjamini-Hochberg procedure was used for multiple testing adjustment within each histological group (all NSCLC, adenocarcinoma and squamous cell cancer).
3. Results
3.1 A pipeline to spatially phenotype immune cells in NSCLC
We established a high throughput workflow that included (1) mIF staining for T- and B-lymphocytes, natural killer (NK) cells, macrophages and subsets of dendritic cells; (2) multispectral scanning; (3) image curation by a specialist pathologist and machine learning–based tissue and cell segmentation; (4) careful data curation and thresholding of marker signals with subsequent immune cell phenotyping; and (5) spatial quantification based on marker expression and cell segmentation (Fig. 1).
Altogether, cancer tissue from 359 NSCLC patients was stained and analysed in this pipeline, resulting in core-specific multilayered images consisting of pixels. Tissue areas with artefacts, necrosis and benign tissue regions were excluded. Finally, 300 cases with high-quality staining for all three marker panels remained, including spatial phenotypes of over 1,300,000 individual immune cells. Cytokeratin staining with a trainable deep learning image analysis algorithm was used for tissue segmentation into tumour and stroma compartments. The expression of different immune markers in each of the immune panels was combined to allocate each cell to one of 13 immune cell subtypes, including T- and B-lymphocytes, macrophages and myeloid cells, NK-cells, natural killer T-cells (NKT) and dendritic cells (DCs) (Fig. 1B). An example of a multiplex immunofluorescent image with corresponding tissue segmentation is shown in Fig. 1A. Each annotated cell was connected to spatial coordinates; thus, cells can be assigned to the tumour or the stroma compartment within the cancer tissue and be related to each other (distances) spatially (Fig. 1A).
3.2 The immune landscape of NSCLC
Quantification of immune cells revealed cell distributions, with the majority of immune cells being located in the tumour stroma and CD4 Eff cells as the most abundant cell type (Fig. 2A, Supplementary Table 3). Only iDCs, mDCs and NK-cells were more frequent in the tumour compartment compared with their stromal numbers, although with generally low absolute numbers. There were only minor differences in the immune environment of adenocarcinoma and squamous cell cancer, with no significant differences after adjustment for multiple testing (Supplementary Table 3).
Correlative analysis revealed that most lymphocyte subtypes occurred together (Fig. 2B; Supplementary Table 4a–c), indicating that immune cells infiltrated co-ordinately. Notably, a generally high lymphocytic infiltration was accompanied by the abundance of M1 macrophages, whereas M2 macrophages were only associated with CD4 Eff and CD8 Eff cells in the stromal area. Furthermore, we observed an unexpected connection between the abundance of plasmacytoid dendritic cells and the density of B-cells in the tumour stroma (Fig. 2B). When correlation analysis was conducted within the main histologic subtypes of NSCLC, specific associations of NK- and NKT-cells with regulatory and effector T-lymphocytes were present in squamous cell carcinoma but not in adenocarcinoma (Supplementary Table 4b–c). The unsupervised hierarchical cluster analysis, including all evaluated immune cell types, resulted in highly segmented immune clusters (Supplementary Fig. 1A). A separated group (n = 37) was seen, represented by high immune infiltration, associated with high numbers of M1 macrophages. This group of patients also showed a prolonged survival in the multimarker context (log-rank test, P = 0.0002, Supplementary Fig. 1B).
Fig. 2Cell densities and correlation matrix. (A) Boxplots showing compartment-specific cell densities (stroma in blue, tumour in red) in all NSCLC, adenocarcinoma (AC) and squamous cell cancer (SCC) histologies, with the upper and lower hinges corresponding to the first and third quartiles and whiskers to the largest and smallest value within 1.5 times the interquartile range. Points outside the whiskers are outliers. (B) Cross-correlation between the densities of all immune cell phenotypes in the complete cohort of NSCLC patients. The stars indicate significance (∗P-value <0.05, ∗∗<0.01, ∗∗∗<0.001). The colour bar at the bottom indicates the strength of correlation, corresponding to the pie charts for each immune cell pair. NSCLC, non–small cell lung cancer.
When analysing the molecular status, we did not detect any association between the estimated tumour mutational burden or any specific mutation and immune cell infiltration after rigorous adjustment for multiple testing (Supplementary Table 5a–d).
3.3 Immune cell infiltration, clinical correlates and survival
Associations between immune cell infiltrations in stroma and tumour compartments and dichotomised clinical parameters—age, performance status, gender, smoking status, and stage—were analysed (Supplementary Table 6). No statistically significant associations were found after adjusting for multiple testing.
In the next step, we analysed immune cell densities in relation to overall survival. This was done for all NSCLC cases, for both main histological subtypes and separately for immune densities in the stroma or cancer compartment, as well as the total cancer area (= cancer stroma area + tumour cell area). High densities of the entirety of all annotated immune cells were associated with longer survival (Supplementary Fig. 2, P < 0.01). The prognostic impact (hazard ratio [HR]: 0.53; P < 0.003) was independent of stage, age and performance status in the multivariable Cox regression analysis (Supplementary Fig. 3). Subsequently, we analysed the prognostic impact for each annotated immune cell subset. The highest prognostic impact was mainly demonstrated when densities in the total tumour area were analysed and not separated in stroma and tumour compartment (Supplementary Fig. 4). In the complete NSCLC cohort, a positive association with longer survival was observed for CD4 Eff, CD4 Treg, CD8 Eff, CD8 Treg, B-cells, M1 and plasmacytoid dendritic cells (pDCs; Supplementary Fig. 4), with HRs ranging from 0.44 to 0.65. Most associations were also observed in the adenocarcinoma and squamous cell cancer subgroups, although with weaker statistical reliability.
After adjustment for the three main clinical features (stage, age and performance status), the positive influence on survival was generally maintained (Supplementary Fig. 5). Besides stage and performance status, the stepwise Cox regression model identified CD4 Eff and CD8 Treg as the most important positive prognostic parameters (Supplementary Fig. 6).
The results indicate that many immune cell subsets have a favourable prognostic impact independent of traditional clinical parameters and that this impact is not restricted to cytotoxic T-lymphocytes.
3.4 Spatial analysis of tumour and immune cells
It is reasonable to predict that not only the quantity but also the location and interrelation of cells are important for the efficacy of cancer immunity. Therefore, we measured the distances between tumour and immune cells through nearest neighbour analysis. NK- and NKT-cells could not be analysed because of their low abundance.
The distribution of the different immune cell subsets in relation to tumour cells showed two main patterns (Fig. 3). In particular, CD4 Eff, CD8 Eff, M1 macrophages and pDCs were aggregated in close vicinity to the tumour cells. In contrast, the other immune cells were more evenly distributed over the area and had a longer mean distance to the tumour cells. Mature dendritic cells revealed an opposite pattern with a predominant distant location.
Fig. 3Distances between tumour cells and immune cells in association with patient survival, unadjusted for clinical parameters. Left: The density plots show the distribution of the mean nearest neighbour distances of tumour cells (PanCKsingle) to immune cell phenotypes. Right: survival forest plot with boxes indicating hazard ratio and hinged whiskers denoting confidence intervals. P-values were adjusted for multiple testing.
In the next step, we analysed the impact of the mean distances between tumour cells and immune cells on overall survival. The same immune cells that demonstrated an association with survival based on densities also showed an association when the mean distances were analysed (Fig. 3). This can be partially explained by the fact that distances correlated with densities. Mature dendritic cells were an exception; a longer distance of these cells to tumour cells was associated with longer survival. After adjustment for clinical parameters, only the proximity of CD8 Eff to tumour cells remained significant in the survival analysis (Supplementary Fig. 7).
When the distances of immune cells to each other were analysed, we found that closer vicinity of adaptive lymphocytic subsets to each other (CD4 Eff, CD4 Treg, CD8 Eff, CD8 Treg and B-cells) were associated with longer survival (Supplementary Fig. 8). For the macrophages and dendritic cell subsets, this relationship was not as clear. When the distances were adjusted to the main clinical parameters (stage, performance status and age; Supplementary Fig. 8), only the proximity between CD8 Treg and B-cells demonstrated an association with longer overall survival.
In the next step, we combined distances and densities in the survival analysis. Predictably, the ratio of densities to distance increased impact on survival in the unadjusted (Fig. 4) and adjusted analysis (stage, performance status and age; Supplementary Fig. 9). However, a remarkable exception was the relation of the lymphocytic subsets to CD8 Treg cells. Here, the favourable prognostic impact of lymphocytes was diminished when the distance to CD8 Treg cells was accounted for.
Fig. 4Survival based on the density/distance metric, unadjusted for clinical parameters. Survival forest plots showing the density through distance metric for all panels and all markers in each histology (all NSCLC, AC, and SCC). Boxes indicate hazard ratios and hinged whiskers indicating confidence intervals. The nearest neighbour analysis does result in different distances depending on which cell is used as the starting point; therefore, the survival associations of the same cell pairs differ (e.g. CD8 Eff–CD8 Treg versus CD8 Treg–CD8 Eff). NSCLC, non–small cell lung cancer.
To evaluate the interdependency of distances and densities with regard to survival, we performed a multivariable stepwise Cox regression analysis, including all immune cell densities and distances of all lymphocytic cell subsets together with age, stage and performance status (Fig. 5). Again, low stage and good performance status showed a strong impact on survival. Independently, high densities of CD163 and M1 macrophages were prognostically favourable. However, not only were pure densities important, but the distances between immune cells also showed an independent prognostic impact. Longer patient survival was observed when CD8 Treg and B-cells, as well as B-cells and CD4 Treg cells, were closer to each other. Also, the closer CD8 Eff cells were located to tumour cells, the longer the survival experienced by the patient, notably with the lowest HR of all selected parameters. Finally, a shorter distance between M2 and M1 macrophages was associated with poor prognosis.
Fig. 5Stepwise multivariable Cox regression including densities and distances. A stepwise multivariable Cox regression model was performed, including all densities (simple cell phenotype) and distances (hyphenated pairs of cell phenotype names; PanCKsingle represents tumour cells) as well as age, stage and performance status. The distances for NK-cells, NKT-cells and dendritic cells were excluded because of high numbers of missing data points due to general low cell counts. The model was iterated with inclusion and stay criteria of P <0.15 until reaching a stable model. The results are ordered by z-score. Boxes indicate hazard ratio and hinged whiskers denote confidence intervals. NK, natural killer; NKT, natural killer T-cells.
Taken together, these results indicate that both the densities and the location of the immune cells bear prognostic information, and both features should be considered in biomarker studies.
4. Discussion
Based on a unique pipeline, combining multiplex immunohistochemical staining, spectral scanning, and advanced image analysis, we were able to quantify and localise major immune cell classes and their subtypes in the in situ environment of diagnostic lung cancer tissue. This analysis not only included traditionally analysed adaptive lymphocyte populations but also their subtypes, as well as the NK-cell family, macrophage subtypes, and the major representatives of dendritic cells. The ‘flowcytometry-like’ purification of immune cell subsets by marker expression allows a more balanced interpretation of the cellular immune reactions, and correlation analysis gives a better picture of inter-relations between cells. Furthermore, the spatial mapping indicates a distinct distribution and neighbourhoods of immune cell subsets, indicating possible mechanistic interplay. We suggest unique cell parameters, including novel immune subsets and spatial interactions that provide prognostic information, independent from the strongest clinical parameters. Finally, we provide a complete data set, including all cell types with their spatial coordinates, and the clinical and molecular background data, as publicly available resources for further, more advanced, or more focused biomarker discovery studies. This first-of-its-kind data set will supplement existing bulk or single-cell RNA sequencing that does not provide the topographical contexture.
There is strong evidence that the type and pattern of immune infiltrates in the local cancer environment have a major impact on the individual prognosis and are associated with response to checkpoint inhibitor therapy. These studies [
], have used immunohistochemistry analysing single or few immune cell markers based on a semiquantitative assessment by bright field microscopy. This unidimensional approach does not allow for subclassification of immune cell types or for a spatial description of cell relations in their in situ cancer environment. These limitations were addressed with our automated pipeline, providing several advantages that extend the description of immune cell patterns in cancer tissue considerably. With higher numeric accuracy, as provided in absolute cell numbers, we confirmed the typical distribution of lymphocytes and macrophages infiltrating either the tumour cell compartment (e.g. CD8 Eff cells) or the surrounding stroma (typically CD4 Eff cells). Here, we provide data on the distribution of NK-cells and NKT-cells as well as dendritic cells, all rarely in situ characterised cell lineages. Obviously, low in absolute numbers, NK-cells predominantly invade tumour nests and are thus in direct contact with cancer cells. The typing of dendritic cell subsets was based on three markers, resulting in three subsets. Notably, while the mature and immature cells were predominantly located in the tumour cell area, the pDC lineage was more evenly distributed over the tumour cell area and stroma. Tumour-infiltrating pDCs are usually linked to a worse prognosis in cancer patients; however, with some exceptions, pDCs are indeed the major type I interferon (IFN) producers of the immune system, a function that is crucially linked to enhanced NK-cell cytotoxicity, Th1 and CD8+ cytotoxic T-cell responses [
]. In line with this, we demonstrated that the presence of pDCs in the tumour cell area was associated with longer patient survival.
In this respect, it should be emphasised that our study was not designed as a biomarker study. By testing multiple markers and combinations, the risk of overfitting is high, despite our stringent adjustment for multiple testing. With the lack of an independent validation cohort, the analyses must be regarded as descriptive. Still, we believe that many interesting observations were made.
Most previous studies have shown that in lung cancer, the number of cytotoxic T-cells, B-cells, plasma cells and M1 macrophages are associated with longer survival [
The prognostic relevance of tumour-infiltrating plasma cells and immunoglobulin kappa C indicates an important role of the humoral immune response in non-small cell lung cancer.
]. Considering that the infiltration of immune cells occurs in a coordinated fashion, that is, most immune cell types infiltrate together and cell numbers correlate strongly, it is difficult to determine the main cellular players. Indeed, our study confirms that the bulk number of immune cells is already a strong independent factor for survival. Also, when each immune cell subset was evaluated separately, most were linked to better prognosis, including CD4 Eff, CD4 Treg, CD8 Eff, CD8 Treg, B-cells, M1 macrophages and NK-cells as well as pDCs. Thus, even the ‘immune suppressive cells’, such as CD8 Treg revealed a paradoxically favourable clinical impact in this unbiased approach. Not until distances were included in the model, it became clear that not only the amount but also the spatial location was important. While the pure density of CD8 Treg cells was associated with an overall good prognosis, the positive clinical impact was eliminated when they were located close to cancer cells or to other immune cells.
This regulatory effect may either be a direct negative effect caused by the CD8 Tregs on the immune environment [
]. Although these findings fit in the perception of immunology, our study provides, to the best of our knowledge, the first in situ evidence of this paradoxical relation.
Another interesting observation is that the CD8 Effcell densities in the adjusted analysis and in the stepwise regression models are not significantly associated with survival. However, when distances either to tumor cells, or to other immune cells were incorporated in the survival analysis, the favourable impact on prognosis was strong and independent (Fig. 4, Fig. S9 and S9). It can be speculated that not all tumour-infiltrating immune cells are antigen specific and the intimacy to tumour cells indicate a certain anti-cancer immune reactivity.
The independent influence of the spatial context is further illustrated in our stepwise Cox regression model. In this unconventional approach, we included all densities and distances and found that distinct distances had a strong and independent impact on survival. For instance, the M1 and M2 phenotype and function were obviously not only dependent on their cell densities or their ratios [
An immune score reflecting pro- and anti-tumoural balance of tumour microenvironment has major prognostic impact and predicts immunotherapy response in solid cancers.
] but also regulated by the topographical relation of both cell types to each other. Many unexpected findings, inconsistent with in vitro and in vivo models in previous studies, might be explained by the fact that the spatial tissue context was not adequately incorporated.
Although our study provides a wealth of data, our analysis is rather conventional with unidimensional metrics (density and distance), alone or in combinations. We did not include more advanced models or methods that could capture complex patterns, such as cell clusters or cell networks. For this purpose, advanced computational approaches, including deep learning models, are necessary [
]. This would ultimately create more powerful classifications of the immune microenvironment, with presumably higher prognostic impact. However, we refrained from this option at this stage to remain in an explainable and interpretable environment where simple distance metrics already convincingly demonstrate the importance of the spatial dimension in the immune context. We are aware that this is just the beginning of a new computational era of immune profiling and explainable artificial intelligence approaches have been developed [
]. To apply such methods or address focused questions, we provide this data set with all cell phenotypes and cell coordinates as a publicly available resource for other researchers.
It should be noted that our immune phenotyping is based on a limited set of markers and might deviate from the definitions used influorescence-activated cell sorting (FACS) analysis. Because cells were located in clusters, and signals from different cell types were colocalized and overlapping, we had to adapt the marker-based annotation. For instance, immune cells in the tumour compartment were often found to overlap with the ‘PanCK’ marker signal. We are aware of this compromise but hope to describe this phenotyping transparently for replication.
Finally, we used TMAs representing only a small area of the whole tumour. A whole section would be the preferable tissue resource for this study because tumours are notoriously heterogenous, including the stroma contexture with different immune cell patterns. This implies the risk that we missed immune phenomena with the approach to focus on centrally located, representative tumor areas. Nevertheless, many biomarker studies demonstrated a surprisingly good agreement between TMA cores and the corresponding whole tissue section [
]. Therefore, we believe that we captured major characteristics of immune cell infiltration in NSCLC.
In conclusion, the combination of state-of-the-art in situ analysis methods with advanced tissue imaging techniques introduces a new dimension to describe the tumour immune microenvironment in cancer that translates to more sensitive and accurate information of clinical relevance. This is not only important for biomarker studies but also essential for understanding the biological complexity of cancer immunity to develop the next generation of immunotherapeutic agents.
Authors' contributions
M.B. contributed to conceptualisation, methodology, software, formal analysis, investigation, data curation and writing, reviewing and editing the article. C.S. contributed to conceptualisation, methodology, validation, formal analysis and reviewing and editing the article. A.L. contributed to methodology, formal analysis and reviewing and editing the article. J.S.M.M., H.E. and M.B. contributed to methodology and reviewing and editing the article. H.B. contributed to methodology, curation, investigation, and resources. A.O.R. contributed to conceptualisation and reviewing and editing the article. M.G. contributed to resources and reviewing and editing the article. J.I. and J.B. contributed to conceptualisation and resources. K.K. and F.P. contributed to resources and writing, reviewing and editing the article. K.L. contributed to conceptualisation and resources. K.L. contributed to conceptualisation, resources and reviewing and editing the article. A.M. contributed to conceptualisation, methodology, software, validation, formal analysis, data curation and writing, reviewing and editing the article. P.M. contributed to conceptualisation, methodology, validation, resources and writing, reviewing and editing the article.
Financial support
This study was partly supported by Swedish Cancer Society, The Lions Cancer Foundation Uppsala, Sweden, Selanders Foundation and The Sjöberg Foundation, Sweden.
Conflict of interest statement
The authors declare the following financial interests/personal relationships which may be considered as potential competing interests:
K.L. is a board member of Cantargia AB, a company developing IL1RAP inhibitors. This does not alter the author's adherence to all guidelines for publication and does not present a conflict of interest. All authors otherwise declare no conflicts of interest.
Appendix A. Supplementary data
The following are the Supplementary data to this article.
Fig. S1Unsupervised hierarchical cluster analysis. A) Unsupervised hierarchical cluster analysis of the cell densities of all immune phenotypes in the total area of NSCLC cases (n = 278). The upper rows specify sex, smoking status, and 2- and 5-year survival, as well as mutation status, estimated TMB, and main histologic subtype. The lowest row indicates the cluster group. B) Kaplan–Meier analysis showing survival of the M1-macrophage associated and non-M1 cluster groups, with better overall survival for the M1 group (log-rank test P = 0.0002).
Fig. S2The general impact of immune cell infiltration on survival. The densities of all immune cells were combined and dichotomised into a high immune cell and low immune cell group in the Kaplan–Meier analysis, including all NSCLC histologies.
Fig. S3Survival association of total immune cell densities. Uni- and multivariable Cox regression with the immune cell densities (total number of immune cells/stroma and tumour area) were analysed together with stage, age, and performance status.
Fig. S4Survival association of individual immune cell densities, unadjusted for clinical parameters. Forest plots showing the association between dichotomised density and overall survival for each immune cell phenotype in each histology (all NSCLC, AC, and SCC) in A) the total area, B) stroma compartment, and C) tumour compartment. Boxes indicate hazard ratios, and hinged whiskers denote confidence intervals.
Fig. S5Survival based on density, adjusted for clinical parameters. Forest plots showing the association between dichotomised density and overall survival for each immune cell phenotype, with adjustment for stage, age and performance status in each histology (all NSCLC, AC and SCC) in A) the total area, B) stroma area and C) tumour area. Boxes indicate hazard ratios and hinged whiskers denote confidence intervals.
Fig. S6Stepwise multivariable Cox regression model based on cell densities. A stepwise multivariable Cox regression was performed that included all cell densities and the most important clinical parameters: age, stage, and performance status. Iteration was done with inclusion and stay criteria of P <0.15 until reaching a stable model. Boxes indicate the hazard ratio, and hinged whiskers denote confidence intervals.
Fig. S7Distances between tumour cells and immune cells in association with patient survival, unadjusted for clinical parameters. Left: The density plots show the distribution of the mean nearest neighbour distances of tumour cells (PanCKsingle) to immune cell phenotypes. Right: survival forest plot with boxes indicating hazard ratio and hinged whiskers denoting confidence intervals, with adjustment for stage, age and performance status. P-values were adjusted for multiple testing.
Fig. S8Distances between immune cell phenotype pairs. Left: density plot showing distribution of mean distance to nearest neighbour between immune cells. Middle: survival forest plot with boxes indicating hazard ratio and hinged whiskers indicating confidence intervals, without and with adjustment for clinical parameters (stage, age and performance status). Right: example images plotted with distinct immune cell phenotypes and distances shown as lines.
Fig. S9Survival based on the density/distance metric, adjusted for clinical parameters. Survival forest plots showing the density through distance metric for all panels and all markers in each histology (all NSCLC, AC, and SCC), with adjustment for stage, age and performance status. Boxes indicating hazard ratio and hinged whiskers indicating confidence intervals. The nearest neighbour analysis does result in different distances depending on which cells are used as the starting point; therefore, the survival associations of the same cell pairs differ (e.g. CD8 Eff–CD8 Treg versus CD8 Treg–CD8 Eff). The PanCKsingle phenotype represents tumour cells.
Multispectral imaging for quantitative and compartment-specific immune infiltrates reveals distinct immune profiles that classify lung cancer patients.
Chemoradiotherapy efficacy is predicted by intra-tumour CD8+/FoxP3+ double positive T cell density in locally advanced N2 non-small-cell lung carcinoma.
Biomarkers associated with beneficial PD-1 checkpoint blockade in non-small cell lung cancer (NSCLC) identified using high-plex digital spatial profiling.
Tumour-infiltrating lymphocyte density is associated with favourable outcome in patients with advanced non-small cell lung cancer treated with immunotherapy.
The prognostic relevance of tumour-infiltrating plasma cells and immunoglobulin kappa C indicates an important role of the humoral immune response in non-small cell lung cancer.
Spatial architecture and arrangement of tumor-infiltrating lymphocytes for predicting likelihood of recurrence in early-stage non-small cell lung cancer.
An immune score reflecting pro- and anti-tumoural balance of tumour microenvironment has major prognostic impact and predicts immunotherapy response in solid cancers.
00100-4&id=gr1.jpg){kind=link}
00100-4&id=gr2.jpg){kind=link}
00100-4&id=gr3.jpg){kind=link}
00100-4&id=gr4.jpg){kind=link}
00100-4&id=gr5.jpg){kind=link}
00100-4&id=figs1.jpg){kind=link}
00100-4&id=figs2.jpg){kind=link}
00100-4&id=figs3.jpg){kind=link}
00100-4&id=figs4.jpg){kind=link}
00100-4&id=figs5.jpg){kind=link}
00100-4&id=figs6.jpg){kind=link}
00100-4&id=figs7.jpg){kind=link}
00100-4&id=figs8.jpg){kind=link}
00100-4&id=figs9.jpg){kind=link}