- 1Jinshan Hospital Center for Tumor Diagnosis & Therapy, Jinshan Hospital, Fudan University, Shanghai, China
- 2Shanghai Key Laboratory of Bioactive Small Molecules, Department of Physiology and Pathophysiology, School of Basic Medical Sciences, Fudan University, Shanghai, China
- 3Shanghai Key Laboratory of Lung Inflammation and Injury, Institute of Clinical Science, Zhongshan Hospital, Fudan University, Shanghai, China
- 4Department of Liver Surgery and Transplantation, Liver Cancer Institute, Zhongshan Hospital, Fudan University, Shanghai, China
- 5Key Laboratory of Carcinogenesis and Cancer Invasion of Ministry of Education, Chinese Academy of Medical Sciences, Shanghai, China
- 6Research Unit of Bench and Clinic Research for Liver cancer Recurrence and Metastasis, Chinese Academy of Medical Sciences, Shanghai, China
- 7Center for Bioinformatics and Computational Biology, Shanghai Key Laboratory of Regulatory Biology, Institute of Biomedical Sciences, School of Life Sciences, East China Normal University, Shanghai, China
- 8Department of Endocrinology and Metabolism, Zhongshan Hospital, Key Laboratory of Metabolism and Molecular Medicine, the Ministry of Education, Fudan University, Shanghai, China
- 9Department of Anesthesiology, Zhongshan Hospital, Fudan University, Shanghai, China
Introduction: Although multiple targeted treatments have appeared, hepatocellular carcinoma (HCC) is still one of the most common causes of cancer-related deaths. The immunosuppressive tumor microenvironment (TME) is a critical factor in the oncogenesis and progression of HCC. The emerging scRNA-seq makes it possible to explore the TME at a high resolution. This study was designed to reveal the immune-metabolic crosstalk between immune cells in HCC and provide novel strategies to regulate immunosuppressive TME.
Method: In this study, we performed scRNA-seq on paired tumor and peri-tumor tissues of HCC. The composition and differentiation trajectory of the immune populations in TME were portrayed. Cellphone DB was utilized to calculate interactions between the identified clusters. Besides, flow cytometry, RT-PCR and seahorse experiments were implemented to explore potential metabolic and epigenetic mechanisms of the inter-cellular interaction.
Result: A total of 19 immune cell clusters were identified and 7 were found closely related to HCC prognosis. Besides, differentiation trajectories of T cells were also presented. Moreover, a new population, CD3+C1q+ tumor-associated macrophages (TAM) were identified and found significantly interacted with CD8+ CCL4+T cells. Compared to the peri-tumor tissue, their interaction was attenuated in tumor. Additionally, the dynamic presence of this newly found cluster was also verified in the peripheral blood of patients with sepsis. Furthermore, we found that CD3+C1q+TAM affected T cell immunity through C1q signaling-induced metabolic and epigenetic reprogramming, thereby potentially affecting tumor prognosis.
Conclusion: Our study revealed the interaction between CD3+C1q+TAM and CD8+ CCL4+T cells and may provide implications for tackling the immunosuppressive TME in HCC.
1 Introduction
With a high mortality rate, liver cancer is the second leading cause of cancer-related deaths under the age of 80 worldwide (1, 2). Hepatocellular carcinoma (HCC) accounts for 70% to 85% of primary liver cancers (3). HCC is an inflammation-driven disease. Surgery, radiofrequency ablation (RFA), transcatheter arterial chemoembolization, and targeted therapies are the most common treatments (4). The tumor microenvironment (TME) of HCC is strongly immunosuppressive; thus, it is extremely important to illustrate the immune characteristics of TME and develop new immunotherapies for HCC.
Immune cell infiltration is a well-known significant regulator of HCC progression (5). The density of tumor-infiltrating CD8+T cells was proved to be an effective prognostic indicator in HCC or many solid tumors (6, 7). It can be affected by multiple regulatory procedures in HCC TME, including the secretion of transforming growth factor β (TGF-β) and interleukin 10 (IL-10), the recruitment of regulatory T cells (Tregs) and myeloid-derived suppressor cells (MDSCs), high levels of programmed cell death 1 (PD-1), programmed death ligand-1 (PD-L1) and so on (4, 5, 8, 9). All the factors cause a rather low response rate of HCC to immunotherapy (10). Therefore it is critical to illustrate the complex network of cell-to-cell interactions within the TME. However, the traditional immunological technology has certain limitations. For example, the chasm between in vitro and in vivo systems, low throuhput, and limited information on immunohistochemistry, et al.
ScRNA-seq is an emerging and powerful tool for investigating the cellular components even rare populations and their interactions in TME (11, 12). It also helps to illustrate connections between TME and clinical outcomes in cancers (13–15). Recently several studies have portrayed the landscape of HCC at single-cell level (16–20). For example, heterogeneity of exhausted T cells (Tex) has been reported (16, 17). The study aims to deepen our understanding of cell-cell interactions and molecular pathways in TME based on scRNA-seq data of paired HCC tissues, and discover new cell subsets which cannot be achieved by traditional methods (21, 22).
We revealed TILs differentiation trajectories and identified certain populations associated with HCC prognosis. The scRNA-seq found that macrophage populations were much more complex than the M1/M2 dichotomy. Although macrophages and T cells are generally considered to belong to different cell lineages, recently a novel macrophage sub-population expressing CD3 molecule was reported in infectious, inflammatory diseases (23). However, the specific role of CD3+ tumor-associated macrophages (TAM) is poorly understood. In the study, we found that a new population of CD3+C1q+TAM regulated the anti-tumor immunity of the tumor-infiltrated CD8+CCL4+T cells through the C1q signaling pathway and subsequent metabolic and epigenetic remodeling, thereby potentially affecting tumor prognosis. Our study supported that the versatile molecular C1q expressed in TAM has functions beyond the complement cascade. The data reveals the interaction between C1q and the metabolism of CD8+T cells and provides implications for regulating immunosuppressive TME.
2 Methods
2.1 Human specimens
Paired carcinomatous and para-carcinomatous tissues were from patients with HCC. No chemotherapy or radiation therapy was performed on patients before tumor resection. Peri-tumor tissue was 3cm away from the edge of the conjugated tumor tissue. Samples were then obtained for the subsequent CD3+cell sorting, single-cell RNA sequencing analysis, or in vitro testing implement.
Peripheral blood was collected from sepsis patients on the 1st, 3rd, 7th and14th day for flow cytometry. We obtained approval from the ethics committee of Zhongshan Hospital, Fudan University, and written informed consent from all HCC patients.
2.2 Sample preparation
We immersed specimens of fresh tumor and adjacent normal tissue in RPMI‐1640 (Gibco) containing 10% FBS, shredded, ground the blocks and then passed the suspension through 40 μm cell strainers. Next, the single cell suspension was centrifuged (1500 rpm, for 10 min) and the supernatant was removed. The bottom cell pellets were resuspended in erythrocyte lysis buffer (Solarbio), kept on ice for 5 min and then washed twice with 1 × PBS.
2.3 Cell isolation and scRNA-Seq
Single-cell suspension was stained with fluorescent-labeled anti- CD3 (0.5%, Biolegend, Cat No.300308, Clone HIT3a) anti-CD45RO (0.5%, Biolegend, Cat No.304210, Clone UCHL1) for 30 min at 4°C. Subsequently, cells were rinsed and resuspended for FACS sorting. CD3+CD45RO+ T cells were isolated through FACS sorting (BD FACS Aria II).
For scRNA-Seq, the sorted cells were counted with a hemocytometer and diluted to 700–1200 cells/μl with targeted cell viability (>70%). Single cells were separated on a Chromium controller (10×Genomics) as the manufacturer’s instructions. 20 cDNA libraries were prepared using Single-Cell 3′ Reagent Kits V2 (10xGenomics, Pleasanton, California) after single cell purification on a Chromium controller (10×Genomics) following the manufacturer’s instructions. Library sequencing was conducted via Illumina sequencer following stringent quality control by fragment analysis (AATI) and the output data was processed through the Cell Ranger pipeline (version 2.1.1, 10×Genomics) default. Cells expressing less than 200 genes or with an improperly high fraction (> 5%) of mitochondrial genes were removed. The raw data was normalized on a log scale, facilitating the following clustering and principal component analysis. CytoTRACE (https://cytotrace.stanford.edu) helped to predict cell differentiation fate based on single-cell RNA-sequencing data Cell development trajectory was performed through ScVelo, an extensible RNA velocity analysis toolkit.
2.4 CD8+ T cells culture
Spleens of wild-type mice were extracted, ground, and filtered to single-cell suspension. CD8+ T cells were purified through negatively magnetic sorting (Biolegend, Cat. 480008), activated with anti-CD3 (eBioscience, Cat.16-0031-85) and anti-CD28 (eBioscience, Cat.16-0281-85) for 3 days, and then cultured with 100 U/ml IL-2 (Peprotech, Cat.200-02-50) or 10 μg/ml IL-15 (R&D, Cat.247-ILB-005) for another 3 days while adding 10 ug/mL or 25ug/mL C1q (Sigma, Cat.C1740-1MG).
2.5 FACS analysis
For flowing staining, viable cells identified through Fixable Viability Stain 510(BD). And cells were stained with CD3 (BioLegend, clone: HIT3a), CD8 (BioLegend, clone: SK1), CD45 (BioLegend, clone: HI30), CD68 (BioLegend, clone: Y1/82A), CD11B (BioLegend, clone: ICRF44), CD80 (BioLegend, clone: 2D10), CD206 (BioLegend, clone: 15-2), CD163 (BioLegend, clone: RM3/1), C1QA (abcam, clone: EPR2980Y), C1QC (abcam, clone: EPR2984Y), PD-1 (BioLegend, clone: RMP1-30), PD-L1 (BioLegend, clone: 10F.9G2), BCL2 (BioLegend, clone: 100), Ki67 ((BioLegend, clone: 11F6), GZMA (Biolegend, clone: CB9),GZMB (Biolegend, clone: QA16A02), IFNγ (Biolegend, clone: XMG1.2), TNFα (Biolegend, clone: MP6-XT22), LDHA (biocompare), Acetyl-Histone H3 (Lys27) Antibody (Cell Signaling), Acetyl-Histone H3 (Lys9) (C5B11) Antibody (Cell Signaling). Intracellular markers were stained before cell re-stimulation by PMA/ionomycin (Biolegend, Cat. 423303). Flow cytometry was performed on BD FACS Aria III flow cytometer and the data was analyzed by FlowJo6 software.
2.6 Metabolic gene-expression analysis by RT-PCR
Relative expression levels of selected genes were quantified by qRT-PCR. cDNA was synthesized via PrimeScript RT Master Mix (TAKARA, Cat. RR036A). TB Green Premix Ex Taq II (TAKARA, Cat. RR820A) was used for qRT-PCR analysis on the Roche LightCycler 480 System with primer sets in Supplementary Table 1. Expression data of all candidate genes were normalized to the housekeeping gene 36b4.
2.7 Metabolic assay
XF-96 extracellular flux analyzer (Seahorse Bioscience) was used for the mitochondrial fitness test, Oxygen consumption rate (OCR) was measured at basal (1.5 × 10^5 cells/well), and after treatment with oligomycin (1 μM), protonophore carbonyl cyanide p-trifluoromethoxyphenylhydrazone (FCCP,1.5 μM), etomoxide (200 μM) in some cases, and rotenone (100 nM) plus antimycin A (1 μM) sequentially. The lactate assay kit (Life Technologies) was used to measure the lactate content in cells and cell cultures.
2.8 Statistical analysis
Post-acquisition analysis was carried out on FlowJo Software (v10.5.30) and representative plots were exhibited. Statistical analysis and corresponding quantitative plots were conducted using GraphPad Prism 8 software. Significance levels were assessed by Student t-test or ANOVA test, where appropriate.
2.9 Data availability
Relevant data in this study are available within the article and supplementary files. Additional information is available upon reasonable request from the corresponding author
3 Results
3.1 Cell populations identification in HCC samples
To identify the cellular diversity in HCC patients, we implemented ScRNA-Seq using the 10x Genomics single-cell 3′ V2 chemistry. We generated a single-cell suspension of 6 samples from the paired tumor (T) and peri-tumor tissue (P) and enriched CD3+T cells by magnetic bead sorting. Barcoded sequencing reads went to the corresponding cells and transcriptome, then individual mRNA molecules were counted through unique molecular identifiers (UMIs). After quality control, we acquired single-cell transcriptomes of 24267 cells. Multiple single-cell analysis was carried out via R package Seurat (version 3.0; https://satijalab.org/seurat/) and high-quality transcriptome and visualization cells were presented over uniform manifold approximation and projection (UMAP) (Figure 1A). A total of 19 immune cell clusters were identified, including CD8+T cells (5 clusters) and CD4+T cells (4 clusters), NKT cells (2 clusters), gamma delta T cells, CD3+monocytes and some clusters undergoing proliferating or differentiating process, et al. (Figure 1A; Table 1). In addition to the typical ones, including naïve T cells (TN), central memory T cells (TCM), effector memory T cells (TEM), recently activated effector memory or effector T cells (TEMRA/TEFF, designated TEMRA hereafter), gamma delta T cells, tumor-Treg cells, and exhausted T (Tex) cells, we also identified two virus-responsive cell clusters, with the CD3-IFI6 and CD3-IGHA1 clusters expressing some markers of type I interferon-stimulated genes.
Figure 1 Single CD3+ cell transcriptome analysis (A) UMAP visualization of T cells clusters based on 24267 single-cell transcriptomes of 6 samples (paired tumor and peri-tumor tissue of HCC), showing the formation of 19 cell clusters. (B) Histogram visualization shows the tissue distribution of 19 clusters. The percentage of cell number of each cluster was compared between tumor and peri-tumor. The means of the 2 groups were tested for significant differences based on the t-test. ns, no significance; *P<0.05 **P<0.01; ***P<0.001; ****P<0.0001. (C) Dot plot showed the expression of representative genes for each cluster. (D) Pie charts showed the proportion of different clusters. (E) Developmental trajectories of CD4+T clusters in UAMP space by scvelo (right) and CytoTRACE (left). Values per phenotype using boxplots in CD4+ clusters. Scores from 0 to 1 represent progressively lower differentiation potential. (F) Developmental trajectories of CD8+T clusters in UAMP space by scvelo (right) and CytoTRACE(left). Summarize the median and distribution of CytoTRACE values per phenotype using boxplots in CD8+T clusters. Scores from 0 to 1 represent progressively lower differentiation potential.
We found that 11 clusters were differentially distributed in P and T. Statistical analysis was performed on the number of cells in each cluster (Figure 1B). Except for the C8, other CD8+ T cell clusters (C0, C6, C11, C13) had a remarkable reduction in tumor tissue than adjacent tissue. Among the CD4+ T cell clusters, C5(Treg) was significantly accumulated in the tumors, centering in establishing and maintaining an immunosuppressive environment. We then listed representative genes for each cluster in Figure 1C. Each cluster exhibited a distinct distribution and characteristics. From C0 to C9, the top 10 clusters of most CD3+cells infiltrated in HCC samples are mainly composed of CD4 and CD8 T cells (Figure 1D). Other cell types, C3(NK) C7(NKT), and C12 (Gamma delta T) are notably reduced in tumors; on the contrary, the proportion of C14: CD3-STMN1 (proliferating cells), C17: CD3-IFI6 (type I interferon signaling pathway activated T cells) increased significantly in the tumor.
3.2 Time differentiation trajectory map of different cluster cells
Next, we explored a hypothetical differentiation trajectory with multiple intermediate states and gene expression gradients by using scVelo and CytoTRACE. The scVelo and CytoTRACE analysis respectively reflect the direction and potential of cell differentiation.
Firstly, we obtained the developmental trajectory of CD4+ and CD8+ T cells on UAMP space via scVelo. As listed in Table 1, we identified 4 clusters of CD4+T cells. It was found that C1 showed progress trends toward 3 other clusters (Figure 1E). Using CytoTRACE to predict the differentiation potential of 4 clusters, in order from high to low is C1> C4> C2> C5(from purple to red, Figure 1E; Supplementary Figure 1A), which is consistent with the results of scVelo. The results showed that C1 (GPR183hiIL7RhiCD4+T cells) was speculated to be a group of cells with “stem cell-like” multiple differentiation potential. The memory CD4 T C2 with up-regulated expression of SELL, CCR7, and LEF1, may differentiate into effector CD4 (C4) or regulatory T cells (C5). C4 has TEMRAs. C4 and C5 developed into two individual branches without further differentiation potentials. The top 10 genes related to stemness (red) or differentiation (blue) were listed in Supplementary Figure 1B by calculating the correlation with CytoTRACE. It is not surprising that multiple genes coding ribosome proteins were related to CD4+T cell differentiation (Supplementary Figure 1B, blue part). Ribosome biogenesis is critical for T-cell activation because of its rate-limiting role in cell growth and proliferation (24). More interestingly, several immune checkpoint genes (TIGIT, CTLA4) were linked to stemness character (Supplementary Figure 1B, red part), reflecting a paradox that dysfunctional TILs in situ may be capable of stem cell-like behavior (25).
By analyzing the differentiation trajectory of CD8+T cells using the same strategy, scVelo depicted that the differentiation direction started from C11 to other clusters (Figure 1F). CytoTRACE predicted that the differentiation potential of the clusters varied from high to low by C8> C0> C11> C13> C6 (Figure 1F, purple to red). Track analysis showed the same trends (Supplementary Figure 1C). Different from the results of CD4+T cells, multiple genes coding ribosome proteins were linked to cell differentiation and stemness of CD8+T cells (Supplementary Figure 1D, blue part, and red part respectively), pinpointing the importance of metabolism in CD8+T cells. Besides, other genes including long non-coding RNA were found associated with CD8+T cell differentiation (Supplementary Figure 1D, blue part), suggesting diverse factors are required for the process.
3.3 Survival analysis of each cluster for liver cancer patients
Based on the TCGA data of 368 liver cancer patients, the survival analysis results obtained 7229 genes related to liver cancer prognosis. In each cluster, we counted the number of survival-related genes. It was found that C14 had the most survival-related genes, followed by C16, and C18, while there were no survival-related genes in C9, C10, and C11. To further understand the role of different cell subsets in prognosis, we performed a K-means cluster analysis on TCGA data based on the survival-related characteristic genes of liver cancer for each cluster. We divided liver cancer patients into two groups and calculated the survival curve to judge the role of this cluster in predicting the prognosis of patients (Figure 2).
Figure 2 Survival analysis of clusters (A–G). Kaplan–Meier survival analysis curves compare the high and low expression of characteristic genes in each cluster based on the TCGA HCC cohort. .
According to survival analysis results, survival curves of 7 clusters (C1, C3, C12, C13, C14, C16, and C18) could be significantly separated. Patients highly expressing C3, C12, and C14-featured genes presented a lower survival probability than those with moderate levels, indicating that the prognosis of patients with liver cancer enriched with this type of cells in vivo is relatively poor (the upper panel, Figures 2A–C). On the contrary, those with high expression ofC1, C13, C16, and C18-featured genes showed improved survival probability, suggesting that these clusters may point to better prognosis in liver cancer patients (the lower panel, Figures 2D–G). Interestingly, cell cycle genes STMN1 and TUBA1B were found increasingly expressed on C14. As switch molecules of the cell cycle from G1 to S phase, they were known for their crucial rules in eukaryotic cell proliferation (26). As shown in Figure 1B, the proportion of C14 increased by 3 times in HCC tumors. Although the destiny of these proliferating cells remains unknown, survival analysis suggests that it is detrimental to the prognosis of cancer patients (Figure 2C).
3.4 The interaction between clusters
It is essential to comprehend the interactions between clusters. Inter-cellular interactions by various molecules, such as receptors, ligands, structural proteins, secreted proteins, hormones, cytokines, neurotransmitters, etc, mediate cell communication and cell activities. Cellphone DB provides comprehensive recourses of curated receptors, ligands, and interactions(https://pypi.org/project/CellPhoneDB/).
In the study, we utilized scRNA-Seq data and Cellphone DB (version 3.1.0) to calculate the interaction scores between clusters, visualized the receptor-ligands interaction quantity, and analyzed the regulatory relationship. The output results are counted and presented in Figures 3A, B. Figures 3A, B respectively showed the interaction of clusters within peri-tumor and tumor tissue. We found that in HCC, the network with the most notable changes involved clusters C5, C16, and C4. In para-cancerous tissue, C6 and C16 had stronger interactions (Figure 3A). However, C16 not only lost self-regulation but also had a weaker connection with C6 in the tumor. In contrast, tumor-infiltrated C6 and C14 gained enhanced interactions (Figure 3B). Overall, the changed connections between clusters of TME suggested the potential immune suppressive mechanisms which need further investigation.
Figure 3 The interaction of CD3+C1q+TAM and CD8+CCL4+T cells a-b. Interaction between clusters of HCC para-cancer (A) and tumor tissues (B). The arrow represents the interaction from this cluster to that cluster. The thickness of the line indicates the value of the interaction pair. The size of the point refers to the interaction of the cluster as an end. The color of the point refers to the interaction of the cluster as a start. (C) Flow cytometric gating of C1q-T (CD11B+ and CD68+) staining on CD45+CD3+ cells of HCC. Right panel shows the proportion of CD11B+ CD68+ in total CD45+CD3+ cells. (D, E). The proportion of M1 (CD86+CD80+) (D), M2 (CD163+CD206+) (E) in total C1Q-TAM (CD11B+ and CD68+) were compared at P and T. (F–K). The expression level of C1QA and C1QC of C1Q-TAM cells (F, G), M1 (H, I) and M2 (J, K) were detected by cell flow cytometry. ns, no significance; *P<0.05; **P<0.01; ***P<0.001.
A meaningful receptor-ligand interaction leads to a physiological response. We identified 141 pairs of ligand-receptor interactions enriched in tumor samples. And the top five were KLRB1_CLEC2D, CD74_MIF, MIF_TNFRSF14, CD8 receptor_LCK, and CD2_CD58. There were 156 pairs in the samples of adjacent tissues, and the top five with the highest frequency were CD74_MIF, CD55_ADGRE5, MIF_TNFRSF14, SELL_SELPLG, CCL4L2_VSIR. The role of these ligand-receptor interactions in shaping TME still needs further investigation.
Next, we defined C16 as CD3+C1q+TAM. The presentative marker genes for C16 included the genes encoding the C1q protein family (C1QB, C1QA, C1QC, and CD68) (Figure 1C). The C6 was CD8+CCL4+T cells identified as effector T cells. Therefore, we speculate that C16 possibly regulates C6 function through the complement signaling pathway. The cluster C16 (CD3+C1q+TAM) found in our study has not been reported before. To further evaluate the existence of C16, we stain CD3+ CD45+ CD68+ CD11B+ cells for single-cell suspension of paired tumor and peri-tumor tissue (Figure 3C). Using flow cytometry, we detected the presence of CD3+C1q+TAM in HCC and found the proportion of this population was significantly lower in tumors (Figure 3C). The observation is consistent with the shift revealed by scRNA-Seq analysis (Figures 3A, B). Furthermore, we labeled M1-LIKE(CD86+CD80+) and M2-LIKE(CD206+CD163+) to explore its entity (Supplementary Figure 2). The staining suggested that the population CD3+C1q+TAM was a mixture of M1 or M2-LIKE cells. There are more M1-LIKE cells in tumor-adjacent tissues (Figure 3D). In contrast, a remarkably increased number of M2-LIKE cells was observed in tumor tissues (Figure 3E). Since the C1q signaling pathway was significantly expressed in C16 (Figure 1C). C1QC expression was elevated in tumors, especially in M2 cells, and C1QA did not change significantly (Figure 3F–K).
To further verify the existence of this population, we collected peripheral blood from patients with sepsis to measure the proportion of CD3+C1q+macrophage in infection by flow cytometry. We found that the CD3+C1q+macrophage ratio increased in the early stage of sepsis (day 1) and decreased after day 3 (Figures 4A, B). Furthermore, we labeled M1-LIKE (CD80+) and M2-LIKE (CD163+) and the data indicated that M2-LIKE cells were significantly higher than M1-LIKE cells in the peripheral blood of sepsis patients (Figures 4C–F). At the same time, we also detected the C1q signaling pathway and found that the expression of C1QC was elevated in the peripheral blood of patients with sepsis, while C1QA had no significant change (Figures 4G–J). These results are consistent with the HCC results, so we believe that CD3+C1q+TAM plays a pro-inflammatory and anti-tumor role in the initial stage of infection and tumors.
Figure 4 The existence of CD3+C1q+macrophage in peripheral blood of sepsis patients (A) The histogram shows the proportion of CD11B+CD68+ in the total CD45+CD3+ cells in the peripheral blood of healthy controls and sepsis patients. (B) Flow cytometry gating of C1Q-M (CD11B+ and CD68+) staining on CD45+CD3+ cells in healthy controls and sepsis. c-f. The histogram shows the proportion of CD163-CD80+ (C), CD80+ (D), CD163+CD80- (E), CD163+ (F) cells in total C1Q-M (CD11B+ and CD68+) by flow staining. (G–J) The expression levels (G, I) and ratios (H, J) of C1QA and C1QC in C1Q-TAM cells (F–G), M1 (H–I) and M2 (J–K) were detected by flow cytometry. ns, no significance; *P<0.05; **P<0.01.
3.5 Effects of C1q on T cell metabolism
Since energy metabolism is essential for T cell differentiation and function, we next explored how C1q modulates T cell metabolism. We tested the expression panel of metabolic genes and found most were upregulated in C1q-treated Teff cells, including HK1, TPI1, GPD2, GPM2, ENO3, PDK1, and LDHA (Supplementary Figure 3). In Tmem cells, only the expression of ALDOA and LDHA were increased after C1q treatment (Supplementary Figure 3). For both Teff and Tmem, the MFI of Carnitine palmitoyl transferase 1A(CPT1α) and Lactate dehydrogenase A(LDHA) increased significantly after C1q treatment (Figures 5A, B). CPT1α and LDHA respectively center in the process of fatty acid oxidation and glycolysis. Besides, the expression of transcription factor HIF1α, an important oxygen sensor and metabolic modulator, also increased after C1q treatment (Figure 5C). The data suggested that C1q treatment can increase T cell metabolism, which was confirmed by the seahorse experiment (Figures 5D, E). Consistent with metabolic changes aforesaid, the lactate content was found increased in Teff and Tmem internally and culture medium of Teff after C1q treatment (Figures 5F, G). Increased glycolysis and higher lactate production further confirmed the effects of C1q treatment on T cell metabolism. Taken together, the results suggested that Teff was considerably responsive to C1q on metabolic regulation.
Figure 5 Effects of C1Q on T cell differentiation and function (A) The expression of LDHA in Teff and Tmem cells with and without C1q. (B) The expression of CPT1α in Teff and Tmem cells with and without C1q. (C) The mRNA expression of HIF1a in Teff and Tmem cells with and without C1q. (D, E). Oxygen consumption rate (OCR) spectra of Teff and Tmem cells with and without C1q using Seahorse XF96 analyzer. (F) The lactate levels in Teff and Tmem cells with and without C1q. (G) The lactate levels in Teff and Tmem cells culture medium. (H) The expression of H3K27ac in Teff and Tmem cells with and without C1q. (I) The expression of H3K9ac in Teff and Tmem cells with and without C1q. ns, no significance; *P<0.05; **P<0.01; ***P<0.001; ****P<0.0001.
Some metabolites have been proved to involve in epigenetic regulation as substrates, donors, cofactors, or antagonists of chromatin- and DNA-modifying enzymes, indicating cross-talk between cell metabolism and differentiation (27). Previous studies have found that LDHA could regulate cell differentiation through the epigenetic axis (28). We next focused on whether the alternation of epigenetic profiling occurred in CD8+ T cells after C1q treatment. Since lactate is generated in the cytoplasm which gives rise to acetyl-CoA and increases the level of histone acetylation, H3K27 and H3K9 acetylation levels were also elevated in both Teff and Tmem after C1q treatment (Figures 5H, I), which is consistent with our observation of increased lactate concentration in treated cells (Figures 5F, G).
3.6 C1q regulates T cell differentiation and function
Then we turned to the impact of C1q from C16 on T cell activation and proliferation. Naïve CD8+T cells were isolated from wild type mice, then stimulated by anti-CD3 and anti-CD28 antibodies, along with C1q to simulate the effect of CD3+C1q+TAM on T cells function. We found that C1q treatment promoted Teff and Tmem differentiation (Figures 6A, B) through T cell responses analysis by flow cytometry. Then we assessed the effects of C1q on the production of IFN-γ, tumor necrosis factor–α (TNF-α), Perforin (PRF), Granzyme A (GZMA) and Granzyme B (GZMB) in Teff and Tmem cells. For both Teff and Tmem, after C1q treatment, cells had increased expression of Ki67 and Bcl2, possibly highlighting a link between C1q and T cell proliferation and survival capacity although the effector function was only significantly enhanced in Teff (Figures 6C, D). The treatment did not affect PD1 or PDL1 expression of T cells, thereby excluding that C1q might induce T cell exhaustion (Supplementary Figure 4).
Figure 6 Effects of C1q on T cell differentiation and function (A) Flow cytometric gating of CD44 and CD62L in CD8+Teff (left) cells and quantification (right) on Teff cells with and without C1q. (B) Flow cytometric gating of CD44 and CD62L in CD8+Tmem (left) cells and quantification (right) on Tmem cells with and without C1q. (C) The expression of IFNγ, GZMA, PRF, GZMB, TNFα, BCL2, and KI67 in Teff cells with and without C1q. (D) The expression of IFNγ, GZMA, PRF, GZMB, TNFα, BCL2 and KI67 in Tmem cells with and without C1q. ns, no significance; *P<0.05; **P<0.01; ***P<0.001. IFNγ, Interferon Gamma, GZMA, Granzyme A, GZMB, Granzyme B, PRF, Perforin, TNFα, Tumor Necrosis Factor Alpha, BCL2, BCL2 Apoptosis Regulator, KI67, Marker Of Proliferation Ki-67.
4 Discussion
Understanding the tumor immune microenvironment of HCC is critical for developing individualized immunotherapy. The study performed cell clustering and annotated them based on scRNA-seq from paired HCC tumor or peri-tumor tissues. Combing the TCGA data, we established the correlation between each cluster and the prognosis of liver cancer patients (Figure 2). Some clusters (C1, C13, C16, and C18) predicted better survival; in contrast, C3, C12, and C14 were linked to poor prognosis.
In view of the prognostic value of C16, it is conceivable that these cells play a crucial role in anti-tumor immunity. C16 is a new population simultaneously expressing CD3, C1Q, CD68, APOEFTL, SELENOP, FCER1G, and HLA-DRA (Table 1). CD68 is a classic marker for macrophages. Complement component C1q is a marker of a particular subpopulation of tissue-resident macrophages and tumor-associated macrophages (TAM). It is hard to define this population as M1(pro-inflammation) or M2(anti-inflammation) based on their DEGs generated from scRNA-seq data. Then we verified its entity by staining well-known markers for M1 (CD80, CD86) and M2(CD206, CD163) signaling, and the results suggested it is a unique population of cells that possess both M1 and M2 features.
With high plasticity, the differentiation process of macrophages often produces heterogeneous results. Recently, a novel subset of macrophages called CD3+ macrophages was found accumulating in both infectious and non-infectious condition (23).They have been reported to originate from monocytic myeloid-derived suppressor cells (M-MDSCs) (29). Meanwhile, CD3+monocyte-derived macrophages can be further differentiated into TCRαβ+ and TCRαβ− sub-population. Adriana et al. found that Tumor necrosis factor (TNF) knockout caused lower frequency of CD3+ TCRαβ+ macrophage and higher mortality in Bacillus Calmette-Guérin (BCG) -infected mice (23). Consistently, the ϵ and ζ chains of CD3 were recently reported to express on RAW, a classical murine macrophage cell line, and present increased transcriptional levels after BCG infection (30). In this study, we identified CD3+C1q+macrophages, a new subset in both HCC tumors (Figure 1) and peripheral blood of sepsis patients (Figure 4). In sepsis samples, the population of CD3+C1q+macrophages had a sharp rise at Day1 after infection and subsequent decline when the infection subsided (Figure 4A). So, we believe that the tumor-infiltrated CD3+C1q+TAM mainly plays a pro-inflammatory and anti-tumor role in the initial stage of infection and tumors although the population is experiencing differential and heterogeneous characteristics. It also could explain why the increased CD3+C1q+TAM links to a better survival prognosis for liver cancer patients (Figure 2E).
We quantified the interaction information between clusters for visual analysis. The results showed that the cluster-cluster interactions differed in HCC tumors and para-cancerous tissue. To further understand the transcriptional heterogeneity of clusters, we discovered transcription factors related to the oncogenesis and development of HCC. For C6 and C16, the gene regulation network of specific high expression transcription factors was constructed. Some key transcription factors were found, such as KDM5B, UQCRB, and CREM in C6, at the same time DAB2 in C16(data not shown).
It has been reported that C1q+TAMs could modulate tumor-infiltrating CD8+T cells via multiple immunomodulatory ligands (31). However, in our study the new identified CD3+C1q+TAM hasn’t been reported before. Therefore, the explanation for the relationship between newly identified C16(CD3+C1q+TAM) and C6 (CD8+CCL4+T) cells also remains unknown. Based on multiple roles of C1q independent of complement activation (32), we explored how CD3+C1q+TAM affected the CD8+CCL4+T cells.
In vitro experiment suggested that this population can promote specific T cell functions by some mechanisms independent of complement activation. In the study (Figure 5), C1q treatment affected multiple metabolic genes expression and which was confirmed by a seahorse experiment. Besides, an abundant metabolite produced during glycolysis, lactate, was increased along with enhanced histone acetylation in both Teff and Tmem after C1q treatment (Figures 5H, I). In other studies, complement C1q binding protein (C1qbp) has been reported to promote differentiation of Teff via metabolic-epigenetic reprogramming (33). Besides, research showed that C1q could restrain responses to self-antigens in CD8+T cells by modulating mitochondrial metabolism (34). Altogether these data link C1q to CD8+T cell metabolism and activation and point the protective role of this newly identified population. Dynamic changes of histone acetylation lead to differential expression of effectors, facilitating rapid and robust responses in CD8+Tmem (35). The nature of incoming environmental factors and the diversion of inherent metabolic programs during T cell activation heavily impact the differentiation and function of effector T cell subpopulations. Thus, our findings depict how C1q influences the metabolism of CD8+T cells and highlights the importance of interaction between complement and immunometabolism in anti-tumor immunity.
With the development of gene detection technology, genomic network-based stratification (NBS)-guided treatments have made great progress clinically. Tumor mutations were proved to help identify subtypes in ovarian, uterine, endometrial and lung cancer and head and neck squamous cell carcinoma (HNSCC) (36–39). A metabolic network-driven approach facilitated to divide HCC into 3 subsets with distinct metabolic features and different prognosis (39). Besides, transcription regulators were also reported to label immunoreactive or silent subtypes in cancer (40). Genomic directed stratification made it possible to implement personalized medicine. Our data showed that C16 was significantly associated with patient outcomes. Together, our results identify C1q+CD3+TAM cells relate with cancer patients prognosis and highlight these cells, as well as C1q itself, as potential targets to address TME-associated immune dysfunction. One of the limitations of our study is the lackness of in vivo animal model to evaluate whether T cell function can regain after restoration of the attenuate crosstalk(C16-C6) within HCC TME. Besides, more researches await to be done to explore how to strengthen this connection pointedly in clinical practice.
Data availability statement
Data transparency. The datasets used or analyzed during the current study are available from the corresponding author on reasonable request. The data presented in the study are deposited in the China National Center for Bioinformation / Beijing Institute of Genomics, accession number HRA003563 (https://ngdc.cncb.ac.cn/gsa-human).
Ethics statement
The studies involving human participants were reviewed and approved by Zhongshan Hospital ethics committee. The patients/participants provided their written informed consent to participate in this study. The animal study was reviewed and approved by Zhongshan Hospital ethics committee.
Author contributions
YY designed and conducted the research, collected and interpreted the data, and wrote the paper. LS, ZC, FL, MM, YC, YL performed and assisted with the experiments, helped to analyze the data and edit the manuscript. QX and GC performed scRNA-seq data analysis. WL, HF and YS provided the patient samples and assisted in the analysis of experimental results. LS, FL, WL, QX, MM, YC, YL, GC and YS provided advice and edited the paper. DW oversaw the research and revised the paper. The authors declared that they had no conflicts of interest. All authors contributed to the article and approved the submitted version.
Funding
This study was supported by the National Natural Science Foundation of China (No.81771672), “Cross key project of mathematics and medical health” of the National Natural Science Foundation of China (No.12026608), Open fund project of Shenzhen BGI Institute of Life Science(No.BGIRSZ20200004) and Special Fund for Clinical Research of Zhongshan Hospital, Fudan University (No.2020ZSLC07), Key Subject Construction Program of Shanghai Health Administrative Authority (ZK2019B30), and Shanghai Engineering Research Center of Tumor Multi-Target Gene Diagnosis (20DZ2254300). The funding agencies had no role in the study’s design and conduct, collecting, analyzing, and interpreting the data or preparation.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2023.1033497/full#supplementary-material
Supplementary Figure 1 | Differentiation and stemness potential analysis of CD4 and CD8 cell subsets.
Supplementary Figure 2 | Flow cytometry gating of M1 (CD80+ CD86+), M2 (CD163+ CD206+) staining on C1q-TAM (CD11B+ CD68+) cells in HCC.
Supplementary Figure 3 | Effects of C1q on T cell metabolic genes expression.
Supplementary Figure 4 | Effects of C1q on T cell exhaustion induction.
Supplementary Table 1 | Primer sequences of genes related to glycolysis and effector functions are provided.
References
1. Llovet JM, Kelley RK, Villanueva A, Singal AG, Pikarsky E, Roayaie S, et al. Hepatocellular carcinoma. Nat Rev Dis Primers (2021) 7(1):6. doi: 10.1038/s41572-020-00240-3
2. Sung H, Ferlay J, Siegel RL, Laversanne M, Soerjomataram I, Jemal A, et al. Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin (2021) 71(3):209–49. doi: 10.3322/caac.21660
3. McGlynn KA, Petrick JL, El-Serag HB. Epidemiology of hepatocellular carcinoma. Hepatology (2021) 73 Suppl 1:4–13. doi: 10.1002/hep.31288
4. Leone P, Solimando AG, Fasano R, Argentiero A, Malerba E, Buonavoglia A, et al. The evolving role of immune checkpoint inhibitors in hepatocellular carcinoma treatment. Vaccines (2021) 9(5):532–50. doi: 10.3390/vaccines9050532
5. Hinshaw DC, Shevde LA. The tumor microenvironment innately modulates cancer progression. Cancer Res (2019) 79(18):4557–66. doi: 10.1158/0008-5472.CAN-18-3962
6. Fridman WH, Zitvogel L, Sautès-Fridman C, Kroemer G. The immune contexture in cancer prognosis and treatment. Nat Rev Clin Oncol (2017) 14(12):717–34. doi: 10.1038/nrclinonc.2017.101
7. Fridman WH, Pagès F, Sautès-Fridman C, Galon J. The immune contexture in human tumours: impact on clinical outcome. Nat Rev Cancer (2012) 12(4):298–306. doi: 10.1038/nrc3245
8. Kumar V, Patel S, Tcyganov E, Gabrilovich DI. The nature of myeloid-derived suppressor cells in the tumor microenvironment. Trends Immunol (2016) 37(3):208–20. doi: 10.1016/j.it.2016.01.004
9. Eggert T, Greten TF. Tumor regulation of the tissue environment in the liver. Pharmacol Ther (2017) 173:47–57. doi: 10.1016/j.pharmthera.2017.02.005
10. Xing R, Gao J, Cui Q, Wang Q. Strategies to improve the antitumor effect of immunotherapy for hepatocellular carcinoma. Front Immunol (2021) 12:783236. doi: 10.3389/fimmu.2021.783236
11. Wagner J, Rapsomaniki MA, Chevrier S, Anzeneder T, Langwieder C, Dykgers A, et al. A single-cell atlas of the tumor and immune ecosystem of human breast cancer. Cellα (2019) 177(5):1330–45.e18. doi: 10.1016/j.cell.2019.03.005
12. Qian J, Olbrecht S, Boeckx B, Vos H, Laoui D, Etlioglu E, et al. A pan-cancer blueprint of the heterogeneous tumor microenvironment revealed by single-cell profiling. Cell Res (2020) 30(9):745–62. doi: 10.1038/s41422-020-0355-0
13. Tirosh I, Izar B, Prakadan SM, Wadsworth MH 2nd, Treacy D, Trombetta JJ, et al. Dissecting the multicellular ecosystem of metastatic melanoma by single-cell RNA-seq. Sci (New York NY) (2016) 352(6282):189–96. doi: 10.1126/science.aad0501
14. Guo X, Zhang Y, Zheng L, Zheng C, Song J, Zhang Q, et al. Global characterization of T cells in non-small-cell lung cancer by single-cell sequencing. Nat Med (2018) 24(7):978–85. doi: 10.1038/s41591-018-0045-3
15. Savas P, Virassamy B, Ye C, Salim A, Mintoff CP, Caramia F, et al. Single-cell profiling of breast cancer T cells reveals a tissue-resident memory subset associated with improved prognosis. Nat Med (2018) 24(7):986–93. doi: 10.1038/s41591-018-0078-7
16. Zhang Q, He Y, Luo N, Patel SJ, Han Y, Gao R, et al. Landscape and dynamics of single immune cells in hepatocellular carcinoma. Cellα (2019) 179(4):829–45 e20. doi: 10.1016/j.cell.2019.10.003
17. Zheng C, Zheng L, Yoo JK, Guo H, Zhang Y, Guo X, et al. Landscape of infiltrating T cells in liver cancer revealed by single-cell sequencing. Cellα (2017) 169(7):1342–56 e16. doi: 10.1016/j.cell.2017.05.035
18. Lim CJ, Lee YH, Pan L, Lai L, Chua C, Wasser M, et al. Multidimensional analyses reveal distinct immune microenvironment in hepatitis b virus-related hepatocellular carcinoma. Gutα (2019) 68(5):916–27. doi: 10.1136/gutjnl-2018-316510
19. Younossi Z, Tacke F, Arrese M, Chander Sharma B, Mostafa I, Bugianesi E, et al. Global perspectives on nonalcoholic fatty liver disease and nonalcoholic steatohepatitis. Hepatologyα (2019) 69(6):2672–82. doi: 10.1002/hep.30251
20. Sun Y, Wu L, Zhong Y, Zhou K, Hou Y, Wang Z, et al. Single-cell landscape of the ecosystem in early-relapse hepatocellular carcinoma. Cellα (2021) 184(2):404–21 e16. doi: 10.1016/j.cell.2020.11.041
21. Zhang L, Li Z, Skrzypczynska KM, Fang Q, Zhang W, O'Brien SA, et al. Single-cell analyses inform mechanisms of myeloid-targeted therapies in colon cancer. Cellα (2020) 181(2):442–59.e29. doi: 10.1016/j.cell.2020.03.048
22. Ren X, Zhang L, Zhang Y, Li Z, Siemers N, Zhang Z. Insights gained from single-cell analysis of immune cells in the tumor microenvironment. Annu Rev Immunol (2021) 39:583–609. doi: 10.1146/annurev-immunol-110519-071134
23. Rodriguez-Cruz A, Vesin D, Ramon-Luing L, Zuñiga J, Quesniaux VFJ, Ryffel B, et al. CD3(+) macrophages deliver proinflammatory cytokines by a CD3- and transmembrane TNF-dependent pathway and are increased at the BCG-infection site. Front Immunol (2019) 10:2550. doi: 10.3389/fimmu.2019.02550
24. Tan TCJ, Knight J, Sbarrato T, Dudek K, Willis AE, Zamoyska R. Suboptimal T-cell receptor signaling compromises protein translation, ribosome biogenesis, and proliferation of mouse CD8 T cells. Proc Natl Acad Sci United States America (2017) 114(30):E6117–E26. doi: 10.1073/pnas.1700939114
25. Vodnala SK, Eil R, Kishton RJ, Sukumar M, Yamamoto TN, Ha NH, et al. T Cell stemness and dysfunction in tumors are triggered by a common mechanism. Science (2019) 363(6434):eaau0135. doi: 10.1126/science.aau0135
26. Lu C, Zhang J, He S, Wan C, Shan A, Wang Y, et al. Increased alpha-tubulin1b expression indicates poor prognosis and resistance to chemotherapy in hepatocellular carcinoma. Dig Dis Sci (2013) 58(9):2713–20. doi: 10.1007/s10620-013-2692-z
27. Fanucchi S, Dominguez-Andres J, Joosten LAB, Netea MG, Mhlanga MM. The intersection of epigenetics and metabolism in trained immunity. Immunityα (2021) 54(1):32–43. doi: 10.1016/j.immuni.2020.10.011
28. Peng M, Yin N, Chhangawala S, Xu K, Leslie CS, Li MO. Aerobic glycolysis promotes T helper 1 cell differentiation through an epigenetic mechanism. Sci (New York NY) (2016) 354(6311):481–4. doi: 10.1126/science.aaf6284
29. Zhang N, Gao X, Zhang W, Xiong J, Cao X, Fu ZF, et al. JEV infection induces m-MDSC differentiation into CD3+ macrophages in the brain. Front Immunol (2022). doi: 10.3389/fimmu.2022.838990
30. Ocaña-Guzmán R, Ramón-Luing LA, Rodríguez-Alvarado M, Voss TD, Fuchs T, Chavez-Galan L. Murine RAW macrophages are a suitable model to study the CD3 signaling in myeloid cells. Cellsα (2022) 11(10):1635–50. doi: 10.3390/cells11101635
31. Dong L, Chen C, Zhang Y, Guo P, Wang Z, Li J, et al. The loss of RNA N(6)-adenosine methyltransferase Mettl14 in tumor-associated macrophages promotes CD8(+) T cell dysfunction and tumor growth. Cancer Cell (2021) 39(7):945–57 e10. doi: 10.1016/j.ccell.2021.04.016
32. Thielens NM, Tedesco F, Bohlson SS, Gaboriaud C, Tenner AJ. C1q: A fresh look upon an old molecule. Mol Immunol (2017) 89:73–83. doi: 10.1016/j.molimm.2017.05.025
33. Zhai X, Liu K, Fang H, Zhang Q, Gao X, Liu F, et al. Mitochondrial C1qbp promotes differentiation of effector CD8(+) T cells via metabolic-epigenetic reprogramming. Sci Adv (2021) 7(49):eabk0490. doi: 10.1126/sciadv.abk0490
34. Ling GS, Crawford G, Buang N, Bartok I, Tian K, Thielens NM, et al. C1q restrains autoimmunity and viral infection by regulating CD8(+) T cell metabolism. Sci (New York NY) (2018) 360(6388):558–63. doi: 10.1126/science.aao4555
35. Araki Y, Fann M, Wersto R, Weng NP. Histone acetylation facilitates rapid and robust memory CD8 T cell response through differential expression of effector molecules (eomesodermin and its targets: perforin and granzyme b). J Immunol (Baltimore Md 1950) (2008) 180(12):8102–8. doi: 10.4049/jimmunol.180.12.8102
36. Hofree M, Shen JP, Carter H, Gross A, Ideker T. Network-based stratification of tumor mutations. Nat Methods (2013) 10(11):1108–15. doi: 10.1038/nmeth.2651
37. Sanati N, Iancu OD, Wu G, Jacobs JE, McWeeney SK. Network-based predictors of progression in head and neck squamous cell carcinoma. Front Genet (2018) 9:183. doi: 10.3389/fgene.2018.00183
38. Urick ME, Bell DW. Clinical actionability of molecular targets in endometrial cancer. Nat Rev Cancer (2019) 19(9):510–21. doi: 10.1038/s41568-019-0177-x
39. Bidkhori G, Benfeitas R, Klevstig M, Zhang C, Nielsen J, Uhlen M, et al. Metabolic network-based stratification of hepatocellular carcinoma reveals three distinct tumor subtypes. Proc Natl Acad Sci United States America (2018) 115(50):E11874–e83. doi: 10.1073/pnas.1807305115
Keywords: immunometabolism, C1q, tumor-associated macrophage, T cell, HCC
Citation: Yang Y, Sun L, Chen Z, Liu W, Xu Q, Liu F, Ma M, Chen Y, Lu Y, Fang H, Chen G, Shi Y and Wu D (2023) The immune-metabolic crosstalk between CD3+C1q+TAM and CD8+T cells associated with relapse-free survival in HCC. Front. Immunol. 14:1033497. doi: 10.3389/fimmu.2023.1033497
Received: 31 August 2022; Accepted: 04 January 2023;
Published: 09 February 2023.
Edited by:
Soumya R. Mohapatra, KIIT University, IndiaReviewed by:
Yinming Liang, Xinxiang Medical University, ChinaAntonio Giovanni Solimando, University of Bari Aldo Moro, Italy
Copyright © 2023 Yang, Sun, Chen, Liu, Xu, Liu, Ma, Chen, Lu, Fang, Chen, Shi and Wu. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Geng Chen, Y2hlbmdlbmc2NjY2NkAxNjMuY29t; Yinghong Shi, c2hpLnlpbmdob25nQHpzLWhvc3BpdGFsLnNoLmNu; Duojiao Wu, d3UuZHVvamlhb0B6cy1ob3NwaXRhbC5zaC5jbg==
†These authors have contributed equally to this work