The sex dichotomy in the efficacy of immunotherapy

After carefully screening and selection, four trials were included for the final analysis based on our search strategy. The main characteristics of these studies were summarized in Suppl. Table 1. All four trials were international, multi-center, phase 3 RCTs. Due to the success of these eligible studies [3], FDA approved their applications in clinical practice in 2019 (KEYNOTE-181) [26], 2020 (ATTRACTION-3) [27], 2021 (KEYNOTE-590) [28], and 2022 (CHECKMATE-648) [29]. OS was the primary endpoint for all trials. The method qualities of these RCTs were generally good as evaluated by Jadad scores, the main issue affecting quality was lack of blinding since only KEYNOTE-590 was double-blind [28].

The analysis for OS included 2360 patients, most of them are men (n = 1977, 84%), and 1336 (57%) subjects were treated with ICIs. As expected, compared with conventional chemotherapy, ICIs decreased the risk of death by 25% (HR, 0.75; 95% CI, 0.69–0.82; P < 0.001; Fig. 1). However, further analysis revealed immunotherapy was associated with favorable outcomes only in men (HR, 0.71; 95% CI, 0.65–0.79; P < 0.001), but not in women (HR, 0.98; 95% CI, 0.78–1.23; P = 0.84). There was a significant difference in the efficacy of immunotherapy between male and female (Pinteraction=0.02). It should be noted that immunotherapy fails to show superior over chemotherapy in women in every single comparison. In contrast, in all the eligible trials, men can benefit from the application of ICIs, suggesting the dichotomy between male and female was robust and conclusive. No substantial heterogeneities were discovered in male (Q = 0.9; I2 = 0.0%; P = 0.92), female (Q = 3.3; I2 = 0.0%; P = 0.51), and overall population (Q = 10.2; I2 = 11.8%; P = 0.34).

Fig. 1

Comparison of the overall survival between immunotherapy and conventional chemotherapy in male, female, and overall population. C, chemotherapy; CI, confidence interval; HR, hazard ratio; I, ipilimumab; N, nivolumab

The genomic landscape of esophageal cancer

To explore the genomic mechanisms underlying the sex dichotomy, we selected, processed, and integrated WES/WGS information from 1425 EC patients in 13 datasets (male, n = 1100; female, n = 325). The key features of the eligible studies were presented in Suppl. Table 2. The potential impact of the heterogeneities in data collection, sequencing method, and analysis approach among various datasets were minimized with tremendous efforts in data verification. An overview of the pooled EC cohort was illustrated in Suppl. Figure 1. The frequencies of non-silent mutations remained relatively constant across various datasets (Suppl. Figure 1 A). The distributions of t-SNE clusters were characterized mainly by mutant genes, no obvious batch effects could be identified (Suppl. Figure 1B). Totally, we identified 134,049 non-silent mutations occurred in 17,572 genes. As shown in Suppl. Figure 1 C, the most common mutant genes were TP53 (78%), TTN (34%), and MUC16 (16%). Their frequencies in single datasets and overall dataset were similar to the pooled mutational frequencies.

In female, 91.08% (n = 296) tumors harbored gene non-silent mutations, and the median numbers in every patient was 80 (interquartile range, 48–115; Fig. 2A). 83.70% of the mutations were missense mutation (Fig. 2B). The most common mutant genes were TP53 (71%), TTN (34%), and MUC16 (16%) (Fig. 2C), and the most common SNV class was C > T (Fig. 2D). We further examined the gene network affected by the most common 50 mutant genes in female by conducting GO analysis (Fig. 2E) and KEGG analysis (Fig. 2F). In male, non-silent mutations were discovered in 94.89% of the samples (n = 1040). Medially, 78 non-silent mutations were identified in every patient (interquartile range, 50–108; Fig. 2G). Missense mutations accounted for 82.25% (Fig. 2H). The highest mutant frequencies were found in TP53 (80%), TTN (35%), and MUC16 (15%) (Fig. 2I), and the most common SNV class was C > T (Fig. 2J). Similarly, we also investigated the gene network affected by the most common 50 mutant genes in men (Fig. 2K and L).

Fig. 2

Genomic mutation landscape in female (upper panel) and male (lower panel) EC patients. The non-silent mutation burden (A and G), mutation subtype (B and H), top 10 most common mutant genes (C and I), SNV class (D and J), GO pathways (E and K), and KEGG pathways (F and L) in male and female EC patients. SNV, single nucleotide variant; GO, gene ontology; KEGG, Kyoto encyclopedia of genes and genomes

For genes whose mutant frequencies over 5%, only 8 genes showed significantly different between men and women (P < 0.05). The mutation of seven genes, namely DMD, FBWX7, ZNF750, OBSCN, MUC4, FAT2, and DNAH11, were enriched in female, while more TP53 mutations were discovered in male. The gene network affected by these 8 mutant genes were illustrated in Suppl. Figure 2.

The sex dichotomy of mutation signatures and reactome pathways

Previous study of mutational processes revealed that both endogenous processes and exogenous exposures resulted in distinctive patterns of mutations, known as mutational signatures [30]. Here, we conducted non-negative matrix factorization analysis of mutational signatures with deconstructSigs [31], then the extracted mutation patterns were compared with Catalogue of Somatic Mutations in Cancer (COSMIC) reference signatures to estimate the mutation burden for each COSMIC signature. Totally, the frequencies of 78 mutational signatures (Suppl. Figure 3) were compared between male and female patients. As shown in Fig. 3, the frequencies of 12 signatures were statistically different between sexes. Among them, the frequencies of SBS1 (known etiology, spontaneous deamination of 5-methylcytosine), SBS2 (activity of APOBEC family of cytidine deaminases), SBS18 (damage by reactive oxygen species), SBS33 (unknown), SBS37 (unknown), and SBS40 (unknown) were increased significantly in women. However, SBS16 (unknown), SBS24 (Aflatoxin exposure), SBS42 (haloalkane exposure), SBS86 (unknown chemotherapy treatment), SBS87 (Thiopurine chemotherapy treatment), and SBS92 (tobacco smoking) were more common in men. Interestingly, C > T were enriched in all these SBS signatures.

Fig. 3

The frequencies of 12 COSMIC reference signatures were significantly different between male and female in EC. Bold black, SBS signature and its known etiologies. Red number, SBS frequency in female patients. Blue number, SBS frequency in male patients. P valued, the difference of SBS frequencies between male and female

The reactome pathway knowledgebase, summarizing the molecular details of signal transduction, DNA replication, and other cellular processes, was a useful tool for discovering the functional associations from somatic mutation profiles in cancer [32]. With this knowledgebase, we investigated 2022 reactome pathways in 1421 eligible patients. As illustrated in Suppl. Figure 4, the most common altered reactome pathways occurred in EC were signal transduction (98.00%), immune system pathway (97.04%), and metabolism of proteins (97.04%). Most of the identified pathway alterations were missense gene mutations (85.90%). Moreover, we examined the alteration of 14 druggable genes, the biomarkers for potential targeted therapy. 14.39% (n = 205) patients harbored these druggable genes in EC, and the most common gene were BRCA2 (3%), BRCA1 (2%), ROS1 (2%), ALK (2%), and EGFR (2%). Among these pathways, compared with female, the frequencies of 84 reactome pathways increased and 47 pathways decreased significantly in male (P < 0.05). In Fig. 4, we showed 5 typical pathways whose altered frequencies were significant higher in female (Fig. 4A), and 5 pathways that were higher in male (Fig. 4B). The full list of 131 altered reactome pathways was illustrated in Suppl. Table 3.

Fig. 4

Ten typical altered reactome pathways occurred in EC patients. (A) Five pathways whose altered frequencies were significant higher in female. (B) Five pathways whose altered frequencies were significant higher in male

Development and validation of a novel sex-related signature (SRS) to predict the efficacy of immunotherapy

Although the pooled analysis revealed that 8 mutations, 12 SBS, and 131 pathways were significantly different between male and female, none of these features demonstrated comparatively dispersive associations with the efficacy of immunotherapy in MSK cohort with 60 EC patients [33], suggesting single feature was insufficient to impact the whole landscape of anti-cancer immune response. Accordingly, we constructed a risk model to develop a comprehensive molecular signature that can predict the efficacy of immunotherapy.

A multivariable Cox regression analysis of the above candidate features was conducted for the OS in the MSK cohort [33], and 21 potential markers related to the efficacy of immunotherapy emerged. After careful evaluation, 6 features including 1 gene mutation (FBWX7) and 5 reactome pathways (signaling by NTRKs, regulation of RAS by GAPs, TP53 regulates transcription of DNA repair genes, DNA double-strand break repair, and FBXW7 mutants and NOTCH1 in cancer) were selected to construct as a risk model defined as SRS. The details of these 5 pathways and their frequencies in male and female were illustrated in Suppl. Figure 5. This model was calculated for every patient with the following formula derived from the alteration status (0 or 1) of the selected six features weighted by their regression coefficient:

SRS score = (0.81×TP53 regulates transcription of DNA repair genes) – (1.14×signaling by NTRKs) – (0.42×regulation of RAS by GAPs) – (1.90×FBXW7) – (0.38×DNA double-strand break repair) + (2.27×FBXW7 mutants and NOTCH1 in cancer).

Based on SRS score and OS, X-tile were applied to determine the optimal cutoff value and categorized patients into high-risk ( > = 1.31) and low-risk (< 1.31) subgroups. With this SRS model, 25 (41.7%) EC patients with low-risk score showed favorable outcomes compared with 35 patients (58.3%) with high-risk score (HR, 0.42; 95% CI, 0.18–0.99; P = 0.03) (Fig. 5A). To prove the generalization of SRS in predicting the efficacy of immunotherapy, we further evaluated the performance of this model in two cohorts enrolled patients with gastric cancer [33, 34]. As expected, low-risk score (n = 29, 50.0%) was associated with longer OS compared with high-risk score (n = 29, 50.0%; HR, 0.16; 95% CI, 0.05–0.50; P < 0.001) (Fig. 5B). In another cohort enrolled 55 patients, patients in low-score subgroup (n = 39) achieved higher objective response rate (ORR; 35.9% vs. 0.0%; P = 0.02; Fig. 5C) compared patients with high-score (n = 11). Tumor mutation burden (TMB) was an FDA-approved biomarker for immunotherapy [3]. As shown in Fig. 5D, E and F, the performances of TMB were not so powerful as our signature in predicting the efficacies of immunotherapy in all three cohorts.

Fig. 5

Development and validation of a novel sex-related signature (SRS) to predict the efficacy of immunotherapy. (A and D) Kaplan-Meier survival analysis in subgroups stratified by SRS (A) or TMB (D) in the training cohort with 60 EC patients. (B and E) Kaplan-Meier survival analysis in subgroups stratified by SRS (B) or TMB (E) in the validation cohort with 58 GC patients. (C and F) Comparison of objective response rates in subgroups stratified by SRS (C) or TMB (F) in another validation cohort with 55 GC patients. GC, gastric cancer; TMB, tumor mutation burden; CI, confidence interval; HR, hazard ratio; CR, complete response; PR, partial response; SD, stable disease; PD, progressive disease

The sex dichotomy of tumor antigens in EC

Due to the differences in XY chromosomes, hormone levels, and genomics, the tumor antigens could induce different immune responses in male and female [35]. With the latest IMGT/HLA database as a reference, here we determined human leukocyte antigen (HLA) alleles with 6-digit precision using HLA-HD [36] in 183 patients with EC. The frequencies of all subtypes of HLA-I antigens (including HLA-A, HLA-B, and HLA-C) were similar between male and female. However, six subtypes of HLA-II antigens were enriched in women. The identified antigens and their prevalence were HLA-DMB*01:03:01 (female vs. male, 57.9% vs. 38.6%; P = 0.03), HLA-DOB*01:04:01 (21.1% vs. 7.6%; P = 0.02), HLA-DQB1*05:02:01 (26.3% vs. 12.4%; P = 0.04), HLA-DQB1*06:09:01 (10.5% vs. 1.4%; P = 0.02), HLA-DRB1*13:02:01 (10.5% vs. 2.1%; P = 0.03), and HLA-DRB1*14:54:01 (13.2% vs. 2.8%; P = 0.02).

Single-cell transcriptome atlas of esophageal cancer

The tumor immune micro-environment played important roles in tumor growth, metastasis, and immunotherapy response [37]. To explore the cell populations within EC, we conducted scRNA-seq and T cell receptor (TCR)-seq analysis on immune cells from 60 patients [38] (Male, n = 44; female, n = 16; Suppl Table 4). After quality controls, we removed the batch effects and integrated the single-cell information with Harmony [39]. Totally, the transcriptome of 105,145 immune cells (CD45+) were included. Clusters obtained from the uniform manifold approximation and projection (UMAP) was annotated using established marker genes. As shown in Fig. 6A, seven major cell populations were identified: CD8+ T cells (n = 35,814), CD4+ T cells (n = 27,224), NK cells (n = 3952), myeloid cells (n = 16,605), plasma cells (n = 7635), B cells (n = 11,997), and Mast cell (n = 1918). The heat-map of marker genes of all seven cell populations were presented in Fig. 6B and Suppl Table 4.

Fig. 6

The single-cell transcriptome atlas of immunity in EC. Major subtypes of immune cells identified by uniform manifold approximation and projection (UMAP). Key marker genes for the classification of seven cell populations. Composition of different subtypes of immune cells in male and female EC patients. The differential expressed genes in male and female. The association between TCR clonotypes and immune cell populations. The proportions of different TCR clonotypes in male and female. TCR, T cell receptor

Previous studies revealed that, compared with men, women tended to accumulate gene mutations which could strongly affect the mutation presented by MHC-II, hence a higher proportion of CD4+ T cell was observed in female than in male [40]. Here, we systematically investigated the composition of EC (Fig. 6C) and also found that there was a higher infiltration of CD4+ T cells in female compare with male patients (31% vs. 23%, P = 0.05). This trend remained relatively consistent across various tissue types and tumor stages. Interestingly, the proportions of other major types of immune were similar between men and women. The differential expressed gene (DEGs) between men and women in all 7 subtypes of immune cells were also examined (Fig. 6D and Suppl. Table 4), many DEGs played key roles in the physiological function of lymphocytes. We then investigated the richness of TCR clonotype in immune cells with scRepertoire [41]. As shown in Fig. 6E, most of the single or small clonotypes were identified in CD4+ T cells, while hyperexpanded, large, and medium clonotypes were mainly discovered in CD8+ T cells. As presented in Fig. 6F, sex dichotomy was also observed in term of proportions of clonotypes. Men harbored more hyperexpanded, large, and medium clonotypes, while single and small clonotypes were enriched in women.

Exhausted CD8+ T cells were highly infiltrated in male patients with EC

Considering the central role of CD8+ T cells in cancer immunotherapy [35], next we investigated the characteristics of this specific cell population in male and female. As shown in Fig. 7A, six major subtypes of CD8+ T cells were identified in UMAP, namely exhausted CD8+ T cells (CD8-Tex, n = 14,170), naïve CD8+ T cells (CD8-Tn, n = 8300), effector CD8+ T cells (CD8-Teff, n = 3128), exhausted terminal CD8+ T cells (CD8-Tex-Term, n = 4520), Effector memory CD8+ T cells (CD8-Tem, n = 1173), and central memory CD8+ T cells (CD8-Tcm, n = 545). The marker genes used in cell classification were presented in Suppl Table 4. Further pseudo-time analysis confirmed the identities of these cell populations. Interestingly, TCR analysis revealed that single clonotypes were enriched in the naïve CD8+ T cells (Fig. 7B).

Fig. 7

The key characteristics of CD8+ T cells in EC. (A) Major subtypes of identified CD8+ T cells. (B) The association between TCR clonotypes and CD8+ T cell populations. (C) Composition of different subtypes of CD8+ T cells in male and female. (D) Pathway activities estimated by GSVA in male and female. (E) The dynamic immune states of CD8+ T cells evaluated by Monocle and their proportions in male and female. (F-I) The proportion of exhausted CD8+ T cells in responders and non-responders in breast cancer (F), basal or squamous cell cancer (G), triple-negative breast cancer (H), and lung cancer (I)

Next, we compared the abundances of different CD8 + subtype cells in male and female EC patients. As shown in Fig. 7C, compared with women, there were significant higher infiltrations of CD8-Tex (49% vs. 31%; P < 0.001), CD8-Tex-Term (16% vs. 9%; P < 0.001), and CD8-Tem (4% vs. 1%; P < 0.001) cells in men, while CD8-Tn (22% vs. 40%; P < 0.001) and CD8-Tcm (0.5% vs. 5%; P < 0.001) were less abundant. The proportion of CD8-Teff were similar between two subgroups. Totally, the exhausted CD8+ T cells accounted for 65% of all identified CD8+ T cells in men. We then conducted GSVA analysis on CD8+ T cells in male and female patients, respectively. As illustrated in Fig. 7D, exhaustion and cytotoxicity signature pathways were highly enriched in men, while naïve signaling pathway obtained the highest score in women. These results further confirmed the robust sex dichotomy in the compositions of CD8+ T cells. We also deciphered the dynamic immune state by inferring the state trajectories with Monocle [24]. As shown in Fig. 7E, seven states along the pseudo-developmental stages were identified in CD8+ T cells. More late-stage of CD8+ T cells were discovered in male, while early-stage of CD8+ T cells were enriched in female.

Based on these results, we hypothesized the proportion of CD8-Tex could be a powerful predictive biomarker in immunotherapy. Indeed, we evaluated the performance of this predictor in four independent cohorts, namely EGAS00001004809 (breast cancer; responders vs. non-responders, 9 vs. 18; Fig. 6F) [42], GSE123814 (basal or squamous cell cancer; responders vs. non-responders, 3 vs. 5; Fig. 6G) [43], GSE169246 (triple-negative breast cancer; responders vs. non-responders, 3 vs. 4; Fig. 6H) [44], and GSE179994 (lung cancer; responders vs. non-responders, 6 vs. 3; Fig. 6I) [45]. As expected, in all four cohorts, the proportion of CD8-Tex in responders were higher than those in non-responders. We also assessed the predictive value of naïve CD8+ T cells in these cohorts. As demonstrated in Suppl. Figure 6, the performance of CD8-Tn as a biomarker was found to be only marginally effective.

Identification and characterization of other cell populations in male and female

Five subtypes of B cells were identified in EC (Suppl. Figure 7 A). The marker genes for classification were presented in Suppl. Figure 7B and Suppl Table 4. The proportions of these five subtypes of B cells showed no significant difference between female and male (Suppl. Figure 7 C). Similar analysis was also conducted on CD4+ T cells (Suppl. Figure 8), and macrophage and DC cells (Suppl. Figure 9). No substantial sex dichotomies were discovered in these cell populations. CopyKAT [46], a computational tool to separate normal cells from malignant cells, was applied to study the epithelial cells. The proportion of ten subtypes of epithelial cells were also similar between male and female (Suppl. Figure 10).

Interestingly, we discovered the abundance of one specific subgroup of NK cells were significantly increased in male patients (Suppl. Figure 11). The featured maker gene in these NK cells was BAG3, which could mediate the CD8+ T cell recruitment [47]. Further investigations were needed to disclose the underlying mechanisms between this specific cell population and the efficacy of cancer immunotherapy.