Characteristics of B lymphocyte infiltration in HPV+ head and neck squamous cell carcinoma.

Human papillomavirus (HPV) is an important etiological factor of head and neck squamous cell carcinoma (HNSCC). HPV+ HNSCC patients usually have a better prognosis, which probably results from the higher infiltration of B lymphocytes. This study was purposed to detect the infiltration of B lymphocyte subsets and the correlation between B lymphocyte subsets and the prognosis in HPV-related HNSCC. In this study, 124 HPV+ and 513 HPV- HNSCC samples were obtained from the Gene Expression Omnibus (GEO) database and The Cancer Genome Atlas (TCGA) database for transcriptomic analysis. Infiltration of B lymphocytes subsets was detected with 7 HPV+ HNSCC and 13 HPV- HNSCC tissues through immunohistochemistry and immunofluorescence. One HPV- HNSCC sample was detected with single-cell sequencing for chemokine analysis. In the results, the infiltration of plasma cells (CD19+ CD38+ ) and memory B cells (MS4A1+ CD27+ ) was higher in HPV+ HNSCC samples. High infiltration of plasma cells and memory B cells was related to a better prognosis. High density of B lymphocytes was positively correlated with high CXCL13 production mainly from CD4+ T lymphocytes in HNSCC. These results indicated that a high density of plasma cells and memory B cells could predict excellent prognosis. CD4+ T lymphocytes might affect B lymphocytes and their subsets through the CXCL13/CXCR5 axis in HNSCC.

INTRODUCTION

Head and neck squamous cell carcinoma (HNSCC) is the sixth most common cancer worldwide. In addition to tobacco and alcohol, human papillomavirus (HPV) has been shown to be another important risk factor for subsets of HNSCC. Interestingly, HPV+ HNSCC carries a more favorable prognosis and better curative effect following radiotherapy compared with HPV− HNSCC. These clinical differences may be related to the HPV effect on the immune system.

The immune system plays an important role in HNSCC. The infiltration of T lymphocytes and their subsets in HNSCC has been extensively studied over recent years. However, B lymphocytes and their subsets in HNSCC have been far less investigated. B lymphocytes have been demonstrated to promote tumor progression or tumor suppression balanced by different B‐cell subsets. Recent studies have shown the B lymphocytes and their subsets are associated with clinical outcome in many cancers. The infiltration of MS4A1+ B lymphocytes in the tumor microenvironment (TME) is related to a favorable outcome in breast cancer, cervical cancer, and colorectal cancer. , , The high infiltration of CD138+ B lymphocytes had a positive effect on the prognosis of colorectal cancer, but it was correlated with a poor outcome in melanoma patients. These findings further suggested that the prognostic value of B lymphocyte subsets may vary in HNSCC, however the infiltration of B lymphocytes subsets and the effect on the prognosis in HNSCC are still unknown factors. The mechanism behind the different infiltration of B lymphocytes and their subsets also requires further study.

At this time, large numbers of studies have confirmed that chemokines and their receptors play an important role for T and B lymphocytes. Different lymphocyte subset interactions through chemokines and chemokine receptors constitute a complex network. CD4+ T and CD8+ T cells lymphocytes can migrate into tumors in response to CXCL9 and CXCL10 stimuli, and B lymphocytes may be recruited by CXCL12 into the cancer microenvironment. The chemokine CXCL13 and its corresponding receptor CXCR5 play a central role in B‐cell trafficking in rheumatoid arthritis. Tumor‐infiltrating CXCL13 producing (CXCR5−) Tfh cells induce B cells and are tightly linked to B‐cell maturation in breast cancer. In HNSCC, a high density of activated B lymphocytes interacts with high levels of tumor‐infiltrating CD8+ T cells through CXCL9 production. This suggests that chemokines are key factors in HNSCC, however the mechanism of chemokines on B lymphocyte infiltration in HNSCC is still unknown. The chemokine network in HNSCC also needs further exploration.

This study, therefore for the first time assessed the infiltration characterization of plasma cells (CD19+CD38+), memory B cells (MS4A1+CD27+), and regulatory B lymphocytes (CD19+IL10+) between HPV+ and HPV− HNSCCs using immunohistochemistry (IHC) and immunofluorescence (IF). The correlation between B lymphocytes subsets and the prognosis of HNSCC was also estimated. Furthermore, differential infiltration mechanisms for B lymphocytes and their subsets were analyzed using the TCGA and GEO databases and single‐cell sequencing data. Finally, a network diagram of the interaction between chemokines and immune cells was constructed to better understand the role of chemokines in the HNSCC microenvironment. These findings extend the current knowledge of the immune landscape in HPV‐related HNSCC and provide new biological markers for prognosis evaluation.

MATERIALS AND METHODS

Patient characteristics

In total, 20 clinical samples from head and neck surgery from the Third Affiliated Hospital of Harbin Medical University in Heilongjiang Province from 2008 to 2015 were chosen. All participants provided written informed consent regarding this study, and ethical approval for the study was obtained from the Ethics Committee of Harbin Medical University (HMUIRB20180026). Participant information was fully protected.

Immunohistochemistry

Immunohistochemical (IHC) staining was performed as previously described. Briefly, serial 3 µm‐thick HNSCC tissues sections were prepared. After deparaffinization and rehydration, slides were covered with Tris‐EDTA (pH 9.0) buffer for antigen retrieval. Slides were blocked with 3% H2O2 solution for 10 min and then the sections were blocked for 1 h with 5% goat serum and incubated with p16, CD4, CD8, CD19, MS4A1, or CXCL13 primary antibodies overnight at 4°C. Next, the slides were incubated with Two‐Step IHC reagents and 3,3‐diaminobenzidine (DAB) solution following the manufacturer’s instruction (antibodies and their optimal dilution ratios are list in Table S1). A series of randomly selected 10 high‐power fields was observed under an Olympus BX51 microscope with a ×20 objective. The cytoplasm of CD4+, CD8+, CD19+ and MS4A1+ cells was stained dark brown. CXCL13+ tissues were stained brown with diffuse distribution. Evaluation criteria of p16+ tissues are described in a previous study. Negative controls were treated identically but without the primary antibodies. IHC staining was assessed using the IHC Profile plugin in ImageJ software.

Double immunofluorescence staining

Immunofluorescence staining was performed as described. Briefly, serial 3 µm‐thick formalin‐fixed, paraffin‐embedded (FFPE) slides were prepared. After blocking with 5% goat serum, slides were incubated with a mixture of 2 different primary antibodies (anti‐CD19 and anti‐CD38, anti‐CD19 and anti‐IL10, anti‐MS4A1 and anti‐CD27 primary antibodies) overnight at 4°C. Then, Coralite594‐conjugated anti‐mouse IgG (H + L) secondary antibody was mixed with Coralite488‐conjugated anti‐rabbit IgG (H + L) secondary antibody (Table S1). Slides were incubated with the mixture for 1 h at 37°C. Vector® TrueVIEW™ Autofluorescence Quenching Kit with DAPI (Table S1) was performed as instructed to reduce the influence of background spontaneous fluorescence. At least 3 randomly high‐power fields were acquired on the Olympus BX51 microscope at the appropriate excitation wavelength for the fluorophore.

The Cancer Genome Atlas and Gene Expression Omnibus data

HNSCC patient sequencing data used in this study were obtained from TCGA (https://cancergenome.nih.gov) with 42 HPV+ and 383 HPV− samples. An additional 212 HNSCC samples with definitive HPV status were obtain from the GEO database (https://www.ncbi.nlm.nih.gov/gds/) including 82 HPV+ and 130 HPV− HNSCC sequencing data. GEO serial numbers are GSE3292, GSE6791, and GSE40074.

Single‐cell RNA‐seq data and analysis

HPV− HNSCC tissue was dissociated into single cells and sequenced using the BD Rhapsody™ analysis system. The gene expression matrix was converted to Seurat objects using the Seurat R package (Version 3.0.2). HNSCC cells were robustly grouped by principal component analysis (PCA) and a t‐distributed statistical neighbor embedding (t‐SNE) algorithm. Then immune cells and tumor cells were separately grouped and annotated based on the cell signature marker. Data processing and analysis were performed as previously described.

Differentially expressed gene analysis

RNA‐seq and array data were analyzed using the R packages DESeq2 and limma respectively. Gene expression landscape associated with key mRNAs and co‐expressed genes was identified using weighted gene co‐expression network analysis (WGCNA). Function annotations were performed using DAVID (Version 6.8) (https://david.ncifcrf.gov) and KEGG (Kyoto Encyclopedia of Genes and Genomes) databases. mRNA‐mRNA interaction networks were analyzed using STRING (https://string‐db.org) and Cytoscape (Version 3.6.0) software.

Immunocyte infiltration and prognosis analysis

Immune infiltrates and prognosis analysis were estimated using the Tumor IMmune Estimation Resource (TIMER; https://cistrome.shinyapps.io/timer). The CIBERSORT (http://cibersort.standford.edu/) tool was used to estimate immune cell types in HNSCC. Immune cells fractions were calculated using The Cancer Immune Atlas (TCIA; https://tcia.at/home).

Statistical analysis

Statistical analyses comparing 2 groups were carried out using the nonparametric Mann‐Whitney U‐test or two‐sided Student t test with SPSS and GraphPad Prism 6 software. Correlation analysis was tested with R package corrplot.

RESULTS

Differential infiltration of immune cells in HPV+/− HNSCC microenvironment

First, 516 HNSCC samples were downloaded from TCGA database, and included mRNA expression profile data and clinical information. Based on the clinical information, 425 samples with complete clinical information were selected, including smoking history, drinking history, tumor recurrence, lymph node metastasis, radiotherapy dose, and HPV status. Based on the HPV status, the HNSCC samples were divided into 2 groups, of which 383 cases were HPV− and 42 cases were HPV+. By analyzing the mRNA expression profile of HPV+ and HPV− HNSCC, 3688 differential genes with P < .05 were obtained. Next, 3688 genes RNA sequence data from 425 samples were analyzed using WGCNA for functional enrichment. Based on the correlation of gene expression, 20 co‐expressed genes modules were constructed (Figure 1A). The correlation between these gene expression profiles and specific clinical traits was analyzed and 1353 genes that were highly correlated with HPV infection status in 5 modules were selected (Figure 1B). These genes were then analyzed for functional enrichment using DAVID and KEGG. From different biological pathways, 6 pathways were immune‐related (P < .05). These pathways were involved in activation and differentiation of B lymphocytes, somatic mutation and gene recombination of immunoglobulin, and function of natural killer T lymphocytes. In addition, these genes were also related to p53 pathway and cell cycle (Figure S1A), therefore HPV could affect the immune response in the microenvironment of HNSCC.

FIGURE 1

Differential infiltration of immune cells around HPV+ and HPV− HNSCC. A, Module functional enrichment analysis of the relationships between differential mRNA with clinical traits. B, Consensus module and HNSCC clinical trait relationships. Each row in the table corresponds to a consensus module. The table is color‐coded to indicate correlation (as shown in the color legend). Red module shows the highest positive correlation with the HPV state. C, Unsupervised hierarchical clustering analysis of immune‐related genes in HPV+ / − HNSCC based on TCGA database. The red horizontal line indicates HPV+ HNSCC clustering, while the blue horizontal line represents HPV− HNSCC clustering. D, Average fractions of immune cells in HPV+ / − HNSCC in TCGA datasets evaluated by TCIA. E, Analysis of infiltrating immune cell of HPV+ and HPV− HNSCC in TCGA database calculated by TCIA. ****P < .0001

To compare the difference between infiltrating immune cells in the microenvironment of HPV+ and HPV− HNSCC, 59 signature genes of different immune cells were screened. Interestingly, the upregulated genes of different immune cells in HPV+ and HPV− HNSCC were closely related to B lymphocytes, T lymphocytes and macrophages (Figure 1C). Next, the different infiltration of 10 types of immune cells in HPV+ and HPV− HNSCC were analyzed using the TCIA database. Relative quantitative analysis showed that B lymphocytes, natural killer cells, CD8+ T lymphocytes and regulatory T lymphocytes in HPV+ HNSCC were significantly increased (P < .0001) (Figure 1D,E).

The effect of different immune cells on the prognosis of HPV+ and HPV− HNSCC was evaluated through TIMER, including B lymphocytes, dendritic cells, neutrophils, macrophages, CD4+ T lymphocytes and CD8+ T lymphocyte. Kaplan‐Meier method was used to show survival curve and log rank test was used to determine the correlation between immune cells and prognosis. The results showed that the high infiltration of B lymphocytes was positively correlated with the prognosis of HNSCC patients (P < .05). Interesting, infiltrations of B lymphocytes had a more significant effect on the prognosis of HPV+ HNSCC patients (Figure 2).

FIGURE 2

Prognostic significance of immune cells assessed in HPV+ and HPV− HNSCC using TIMER

In short, HPV had a great impact on the immune microenvironment of HNSCC, and there were huge differences in immune cell infiltration between HPV+ and HPV− HNSCC samples. Among them, B lymphocytes had a great significant difference that was related to the better survival of patients. The high infiltration of B lymphocytes might be an immune reason for better prognosis of HPV+ HNSCC patients.

Difference in infiltration of B lymphocyte subpopulation in microenvironment between HPV+ and HPV− HNSCC

Further classification of B lymphocytes helped us to understand the role of immune cells in the microenvironment in the progression of HNSCC. Based on a literature review and arrangement, the main surface markers of these 4 B lymphocyte subsets were summarized, and subsequently 2 markers were selected for experimental verification (Figure 3A).

FIGURE 3

The frequency of B lymphocytes and major subsets in HPV+ HNSCC increased. A, The major subsets (plasma cell, memory B cell, and Breg) of B lymphocytes and the representative genes of different B lymphocytes summarized from recent literature. The markers in red indicate populations selected for experiments. B, Analysis of relative fractions of total B cells and B cells subtypes evaluated by TCIA. C, Representative IHC images of p16+ and p16− HNSCC tissues. D, Statistic analysis of IHC scores of CD19+/MS4A1+ B lymphocytes in HNSCC with HPV status. E, Representative IHC images of CD19+B lymphocytes and MS4A1+ B lymphocytes in HPV+ and HPV− HNSCC. Tumor structure shown by H&E staining. **P < .01, ****P < .0001

The infiltration of different B lymphocyte subsets in HNSCC was analyzed using RNA‐seq date from a GEO dataset (GSE6792) using CIBESORT. The results showed that the infiltration proportion of 22 immune cell subsets including memory B cells and plasma cells increased in HPV+ HNSCC (Figure S1B,C). TCIA was used to analyze the infiltration of B lymphocyte subpopulations with the HNSCC data using TCGA database. These results showed that the infiltration of memory B cells and plasma cells was increased in HPV+ HNSCC samples (Figure 3B). Unsupervised hierarchical clustering analysis also showed the differential expression of B lymphocyte subset characteristic genes between HPV+ and HPV− HNSCC using TCGA dataset (Figure S1D).

Immunohistochemistry (IHC) further confirmed the results of bioinformatics analysis. Firstly, 20 HNSCC tissues could be divided into HPV+ and HPV− groups based on the expression of p16 protein (Figure 3C). Clinical information of these patients is listed in Table 1. The infiltration of B lymphocytes into tissues was estimated and the results showed that the infiltration of MS4A1+ and CD19+ B lymphocytes increased significantly in HPV+ HNSCC. Most lymphocytes were concentrated in the tumor margin, however the infiltration of MS4A1+ and CD19+ B lymphocytes was slightly different (Figure 3D,E). To detect the infiltration and location of major B‐cell subsets, memory B cells (MS4A1+CD27+), plasma cells (CD19+CD38+) and Breg cells (CD19+IL10+), double immunofluorescence was detected in HNSCC samples. The results showed that the infiltration of MS4A1+CD27+ B cells/CD19+CD38+ B lymphocytes was increased in HPV+ HNSCC samples, while CD19+IL10+ B lymphocytes showed no significant difference (Figure 4). The correlation between B lymphocyte subset‐related genes and patient survival was evaluated using the Kaplan‐Meier method with TIMER. The results showed that the high expression levels of plasma cell and memory B‐cell related genes was more significantly correlated with the prognosis of patients in HPV+ HNSCC (P < .05) (Figure 5).

TABLE 1

Clinical characteristics of HNSCC patients

Characteristics	HPV−(N = 13)	HPV+(N = 7)	P‐value
Gender, n (%)
Male	10 (72.4)	4 (27.6)	.357
Female	3 (50.0)	3 (50.0)
Age (y)
Mean	58	57	\
Range	45‐77	38‐74
Anatomic site, n (%)
Larynx	9 (64.3)	5 (35.7)	.919
Oropharynx	4 (66.7)	2 (33.3)
Tobacco, n (%)
Yes	9 (60.0)	6 (40.0)	.417
Never	4 (80.0)	1 (20.0)
Alcohol, n (%)
Yes	4 (57.1)	3(42.9)	.589
Never	9 (69.2)	4 (30.8)
T stage, n (%)
T1‐T2	10 (66.7)	5 (33.3)	.787
T3‐T4	3 (60.0)	2 (40.0)
N stage, n (%)
N0	7 (53.8)	6 (46.2)	.154
N1‐N3	6 (76.7)	1 (23.3)
M stage, n (%)
M0	12 (63.2)	7 (36.8)	.452
M1	1 (100)	0 (0.00)

FIGURE 4

Infiltration of MS4A1+CD27+ B lymphocytes and CD19+CD38+ B lymphocytes in HPV+ HNSCC. A, The infiltration and location of major B‐cell subsets, memory B cells (MS4A1+ CD27+), plasma cells (CD19+CD38+) and Breg cell (CD19+IL10+) in HPV+ HNSCC detected by immunofluorescence. B, Statistic analysis of immunofluorescence results in HNSCC with HPV status. **P < .01. ****P < .0001

FIGURE 5

Infiltration of MS4A1+ CD27+ B lymphocytes and CD19+CD38+ B lymphocytes in HPV+ HNSCC related to better prognosis

Infiltration of B lymphocytes in HPV+ HNSCC affected by CXCL13

Through analysis of the interaction network of the differential genes in GEO data for HNSCC (GSE3292, GSE6791, GSE40074), important gene modules were selected, among which the most significant gene module networks were chemokine interaction networks (Figure 6A,B).

FIGURE 6

Infiltration of B lymphocyte affected by CXCL13 in HPV+ HNSCC. A, Volcano plots of differential mRNA expression in HPV+ and HPV− HNSCC with GSE3292, GSE6791, GSE40074 datasets, from left to right respectively. Red and green circles indicate high and low mRNA expression, respectively. Black circles indicate mRNAs with P‐value > .05. B, The mRNA‐mRNA functional network with differentially expression genes and the hub network of differentially expressed chemokines in HPV+ and HPV− HNSCC. Red and blue squares indicate high and low mRNA expression, respectively. C, The interaction between B lymphocyte markers (CD19, MS4A1) with the chemokines selected in the module. D, Expression of CXCL13 and CXCR5 was related to B lymphocytes (CD19, MS4A1). Spearman's correlation was calculated using TIMER

Chemokines are secreted proteins responsible for immune cell transport. Further analysis of the relationship between B lymphocyte markers and the chemokine module showed that only 9 genes were directly related to CD19 and MS4A1 (Figure 6C). Interestingly, CXCL13 and its receptor CXCR5 played an important role in the module. CXCL13 is a chemokine ligand originally termed B‐cell attracting chemokine 1 (BCA‐1), the CXCL13/CXCR5 axis had been shown to be involved in regulating B lymphocyte migration and promoting inflammation. Furthermore, TCGA dataset analysis showed that the expression levels of CXCL13, and CXCR5 were significant correlated to CD19 and MS4A1 in HNSCC (Figure 6D) and were higher in HPV+ HNSCC with statistical significance (P < .05) (Figure 7A). Immunohistochemical results of CXCL13 showed the same trend (Figure 7B,C). Correlation analysis between CXCL13 and B‐cell subsets showed that CXCL13 was mainly related to CD19+CD38+ B lymphocytes (Figure 7D). Furthermore, the high expression levels of CXCL13 and CXCR5 were statistically positively correlated with the prognosis of HPV+ HNSCC patients (Figure S2A).

FIGURE 7

Expression of CXCL13 in HPV+ HNSCC increased. A, Differential expression of chemokines CXCL13 and CXCR5 in HPV+ and HPV− HNSCC using TCGA datasets. B, Statistical analysis of CXCL13 IHC scores in HNSCC with HPV status. C, Representative images of CXCL13 expression in HPV+ and HPV− HNSCC analyzed by IHC. D, Correlation analysis of CXCL13 and B lymphocytes subsets. **P < .01, ****P < .0001

In brief, these results suggested that overexpression of CXCL13 regulated the infiltration of B lymphocytes in HPV+ HNSCC.

CD4+ T lymphocytes interact with B lymphocytes through the CXCL13 / CXCR5 axis in an HNSCC microenvironment

To further explore the source of CXCL13, an HPV− HNSCC sample was analyzed using single‐cell sequencing and the infiltration of 5 types of cell was analyzed (cancer cells, immune cells, endothelial, epithelial, fibroblasts). Interestingly, cancer cells and immune cells secreted CXCL13 at the same time, but the effect of immune cells was more obvious (Figure 8A,B). Single‐cell sequencing data were further classified into 8 group of immune cell subsets based on their characteristic gene expression (Figure 8C). CD4+ and CD8+ T lymphocytes were the main CXCL13‐secreting cells among them (Figure 8D). Correlation analysis results showed that CD4+ T lymphocytes may be related to CD19+CD38+ B lymphocyte (Figure S2B), however the different infiltration of CD4+ and CD8+ T lymphocytes also affected the expression of CXCL13. In this sample, 454 CD4+ T lymphocytes and 785 CD8+ T lymphocytes were identified. In total, 72% of CD4+ T lymphocytes could secrete CXCL13, while only 57% of CD8+ T lymphocytes were CXCL13‐secreting cells (Figure 8E,F). Considering the number of T lymphocytes and the expression level of CXCL13, the ability of CD4+ T lymphocytes to secrete CXCL13 was higher compared with that of CD8+ T lymphocytes in HNSCC samples (Figure 8G,H). IHC analysis of serial sections also revealed that the infiltration of CD4+ T lymphocytes was often associated with a high infiltration of B lymphocytes in the same location of HPV+/− HNSCC samples (Figure 9). The other 2 groups showed that this phenomenon was not accidental (Figure S2C). CD4+ T lymphocytes were the core effector cells in anti‐tumor immunity. The infiltration frequency of CD4+ T lymphocytes was positively correlated with the prognosis of HPV+ HNSCC. These results suggested that CD4+ T cells mainly affected the infiltration of B lymphocytes by secreting chemokine CXCL13 in HNSCC (Figure 10).

FIGURE 8

Analysis of CXCL13 secretory cells by single‐cell sequencing. A, Cell clustering and gene annotation for single‐cell sequencing of head and neck tissue with a t‐distributed statistical neighbor embedding (t‐SNE) algorithm. B, Expression of CXCL13 in different cells shown with a violin chart. C, The expression of CXCL13 in immune cells with t‐SNE algorithm. D, The expression of CXCL13 in different immune cells shown with a violin chart. E, The different proportion of CXCL13‐secreting cells in CD4+ and CD8+ T lymphocytes. F, Number of CD4+ and CD8+ T lymphocytes in single‐cell sequencing samples. G, The relative expression levels of CXCL13 in CD4+ and CD8+ T lymphocytes. H, The ability of CD4+and CD8+ T lymphocytes to secrete CXCL13. ****P < .0001

FIGURE 9

The infiltration of B lymphocytes, T lymphocytes and CXCL13 in HPV ± HNSCC tissues. A, Representative images of MS4A1+ B lymphocytes, CD4+ T lymphocytes, CD8+ T lymphocytes and CXCL13 in HPV+ and HPV− HNSCC analyzed by IHC. Tumor structure shown by H&E staining. B, Correlation analysis of MS4A1+ B lymphocytes, CD4+ T lymphocytes and CD8+ T lymphocytes

FIGURE 10

Diagram of CD4+ T lymphocytes and B lymphocytes interaction through the CXCL13/CXCR5 axis. The HPV gene could integrate with the host gene at different sites that makes HPV+ HNSCC express more neoantigens and be more immunogenic, leading to the different infiltration of T and B lymphocytes.

DISCUSSION

This study for the first time clarified the unique frequency and distribution characteristics of memory B cells (MS4A1+CD27+), plasma cells (CD19+CD38+) and Breg cells (CD19+IL10+) in the HPV+ HNSCC microenvironment, and concluded that CD4+ T lymphocytes pass through the CXCL13/CXCR5 axis to affect the infiltration of B lymphocytes and its subsets.

Generally, after the immune response, typical memory B cells will recirculate, while special memory B cells resident in tissue are responsible for immune monitoring and a rapid immune response. In our results, the infiltration of plasma cells (CD19+CD38+) and memory B (MS4A1+CD27+) cells increased in HPV+ HNSCC samples.Compared with MS4A1+CD27+ B cells, the infiltration frequency of CD19+CD38+ B cells was higher in the same area. It suggested there were mass memory B cells resident in the head and neck position for immune surveillance of HPV infection, and that HPV could activate more B lymphocytes for antibody producing. Trust analysis of HPV+/− HNSCC sequencing data also found that B lymphocytes could produce more abundant types of immunoglobulins in HPV+ HNSCC, and indicated that the better prognosis of HPV+ HNSCC patients might be related to the secretion of HPV‐reactive antibodies (Figure S2D). Furthermore, whether there remains a mutual transformation between resident memory B cells and plasma cells still needed further clarification.

Among the analysis of the reasons for the higher infiltration of B lymphocytes, chemokines were found to be an interesting factor. By constructing the chemokine network, a group of chemokines CXCL13/CXCR5 was found to be closely related to B lymphocytes, and the higher expression of CXCL13 was positively related to the HPV+ HNSCC prognosis. CXCL13 is commonly secreted by follicular dendritic cells, macrophages, and germinal center T cells. Its receptor CXCR5 is often distributed on peripheral mature B lymphocytes and a small subset of T lymphocyte. Previously, CXCL13 has been shown to be highly expressed in a variety of cancers and is associated with tumor metastasis and prognosis, but the mechanism of CXCL13 in cancer remains controversial. Zhenyu Zh et al showed that the CXCL13/CXCR5 axis could promote the growth and invasion of colon cancer cells through the PI3K/AKT pathway. Ohandjo et al suggested that the high expression of the CXCL13/CXCR5 axis in prostate tumors enhanced tertiary lymphoid structure (TLS) formation. In HPV+ HNSCC, the existing TLS that could explain better prognoses, but more experiments are needed to verify this proposal. Furthermore, the molecular mechanism of the CXCL13/CXCR5 axis and the role of other chemokines in HNSCC remain to be understood.

In this study, we first demonstrated the source of CXCL13 in HNSCC. With single‐cell sequencing analysis, results showed that immune cells secreted more CXCL13 compared with tumor cells, especially CD4+T lymphocytes. CD4+T lymphocytes serve a variety of biological functions related to their subpopulations. Th1 plays an active role in anti‐tumor immunity, while Tregs play the opposite role. For a productive immune response, activated B cells must interact with activated CD4+ T lymphocytes that are specific for the same antigen. Nordfors et al found a significantly higher number of CD4+ TILs in HPV+ OPSCCs and Chen et al also showed that the infiltration of CD4+ T cells increased in HPV+ HNSCC. This suggested that CD4+ T lymphocytes could activate B lymphocytes for immune response in HPV+ HNSCC. In the tumor microenvironment, chemokines can be expressed by tumor cells and other cells, including immune cells and stromal cells. Single‐cell sequencing analysis showed that tumor cells and CD8+ T cells also secreted CXCL13 as well as CD4+ T lymphocytes in HNSCC. Gu‐Trantien et al demonstrated that CXCL13 producing (CXCR5−) TFH cells promoting local memory B‐cell differentiation in breast cancer. Although the number of CD8+ T cells is large, its ability to secrete CXCL13 is not as good as that of CD4+ T cells, based on single‐cell sequencing data. This may be due to the different interaction between CD8+ T cells and B lymphocytes. Kamila et al observed in oropharyngeal squamous cell carcinoma that CD8+ T cells could also interact directly with B lymphocytes. The complex regulatory network between T lymphocytes and B lymphocytes remains to be explored. The absence of sufficient HPV+/− HNSCC samples limited the single‐cell sequencing analysis results and more HNSCC samples are needed to prove the main source of CXCL13. In addition, single‐cell sequencing analysis could not divide the data into more specific immune cell subsets as the depth and the method of sequencing was not perfect. Whether there are other unclassified CXCL13 secreting immune cells still needs further study.

Compelling reports have suggested that IFN‐I signaling in cancer cells and tumor microenvironment critically modulates the recruitment and activation of effector T lymphocytes through chemokines such as CXCL9 and CXCL10, which promote the recruitment of T lymphocytes to the tumor area. , Interestingly, these 2 chemokines were also highly expressed in HPV+ HNSCC, which suggested that tumor cells could recruited CD4+ T cells through CXCL10. In addition, Kamila and colleagues suggested that tumor‐infiltrating B lymphocytes might recruit CD8+ T cells through CXCL9 and, due to a highly activated phenotype, contribute by secondary costimulation for maintenance of CD8+ T lymphocytes in the oropharyngeal squamous cell carcinoma tumor microenvironment. The interaction between other immune cells and chemokines in the HNSCC microenvironment remains to be further studied.

It has been stated that the infiltration of plasma cells (CD19+CD38+) and memory B cells (MS4A1+CD27+) was higher in HPV+ HNSCC compared with HPV− HNSCC, whereas there was no significant difference in the infiltration of regulatory B lymphocytes (CD19+IL10+) between HPV+ and HPV− HNSCC using immunofluorescence. These B lymphocyte subtypes in the HPV+ HNSCC microenvironment were related to a better prognosis and, for the first time, these results demonstrated that CD4+ T cells could affect B lymphocyte infiltration through the CXCL13/CXCR5 axis in HPV+HNSCC. These findings suggested that CXCL13 may be a key target and a new prognostic indicator for regulating HNSCC in the future. At the same time, these results further deepened the understanding of the role of different subtypes of B lymphocytes in the HPV+ HNSCC microenvironment; this would be more helpful to study the impact of HPV on the HNSCC microenvironment and expand the current understanding of the HPV‐related HNSCC immune landscape. It also provided an experimental basis and a theoretical basis for optimizing the HNSCC immunotherapy scheme through B lymphocytes and chemokines, and provided new biological markers for prognostic evaluation, which is of great significance for optimizing HNSCC treatment and improving patient prognosis.

DISCLOSURE

The authors declare that they have no conflicts of interest.

Fig S1

Click here for additional data file.

Fig S2

Click here for additional data file.

Table S1

Click here for additional data file.