Skip to main content

ORIGINAL RESEARCH article

Front. Immunol., 23 May 2024
Sec. Autoimmune and Autoinflammatory Disorders : Autoimmune Disorders

Identification and validation of biomarkers in membranous nephropathy and pan-cancer analysis

Yue YangYue YangGu-ming ZouGu-ming ZouXian-sen WeiXian-sen WeiZheng ZhangZheng ZhangLi ZhuoLi ZhuoQian-qian Xu*&#x;Qian-qian Xu*†Wen-ge Li*&#x;Wen-ge Li*†
  • Department of Nephrology, China-Japan Friendship Hospital, Beijing, China

Background: Membranous nephropathy (MN) is an autoimmune disease and represents the most prevalent type of renal pathology in adult patients afflicted with nephrotic syndrome. Despite substantial evidence suggesting a possible link between MN and cancer, the precise underlying mechanisms remain elusive.

Methods: In this study, we acquired and integrated two MN datasets (comprising a single-cell dataset and a bulk RNA-seq dataset) from the Gene Expression Omnibus database for differential expression gene (DEG) analysis, hub genes were obtained by LASSO and random forest algorithms, the diagnostic ability of hub genes was assessed using ROC curves, and the degree of immune cell infiltration was evaluated using the ssGSEA function. Concurrently, we gathered pan-cancer-related genes from the TCGA and GTEx databases, to analyze the expression, mutation status, drug sensitivity and prognosis of hub genes in pan-cancer.

Results: We conducted intersections between the set of 318 senescence-related genes and the 366 DEGs, resulting in the identification of 13 senescence-related DEGs. Afterwards, we meticulously analyzed these genes using the LASSO and random forest algorithms, which ultimately led to the discovery of six hub genes through intersection (PIK3R1, CCND1, TERF2IP, SLC25A4, CAPN2, and TXN). ROC curves suggest that these hub genes have good recognition of MN. After performing correlation analysis, examining immune infiltration, and conducting a comprehensive pan-cancer investigation, we validated these six hub genes through immunohistochemical analysis using human renal biopsy tissues. The pan-cancer analysis notably accentuates the robust association between these hub genes and the prognoses of individuals afflicted by diverse cancer types, further underscoring the importance of mutations within these hub genes across various cancers.

Conclusion: This evidence indicates that these genes could potentially play a pivotal role as a critical link connecting MN and cancer. As a result, they may hold promise as valuable targets for intervention in cases of both MN and cancer.

1 Introduction

Membranous nephropathy (MN) encompasses a spectrum of disorders distinguished by the thickening of the glomerular basement membrane (GBM) and the accumulation of immune complexes beneath the epithelial cells on the outer surface of the GBM (1). MN can be categorized into primary MN and secondary MN based on its underlying causes. The causative factors for MN are multifaceted, encompassing infections, autoimmune diseases, malignancies, pharmaceutical agents, heavy metals, and more.

For decades, a connection between MN and cancers has been established (2), tracing back to 1966 when Lee et al. (3) disclosed that 11% of nephrotic syndrome patients had carcinoma, the solid tumors most commonly associated with MN are lung and gastric cancers, followed by renal cell carcinoma, prostate cancer, and thymoma (4), these glomerular lesions are thought to be paraneoplastic. Nonetheless, the precise prevalence of cancer among MN patients remains elusive, with estimates ranging from 5% to 22%. A meta-analysis indicated that the cancer prevalence in MN patients was approximately 10% (5). This association is frequently observed in patients aged over 60, and most cancer cases are identified either before or concurrent with the diagnosis of MN (6), unfortunately, in most cases, the exact pathogenesis is unclear.

Microarray-based gene expression profiling is a widely utilized, high-throughput technique for investigating complex disease mechanisms (7). This approach has facilitated the identification of diagnostic and prognostic biomarkers, disease classification, monitoring of treatment responses, and understanding of disease pathogenesis (8, 9). Recently, numerous bioinformatics studies have aimed to elucidate the pathogenesis of membranous nephropathy. These studies have included searching for biomarkers (10, 11), examining their association with immune infiltration (12), identifying hub genes involved in disease mechanisms (13), and investigating the miRNA-mRNA regulatory networks related to podocyte autophagy, lipid metabolism, and renal fibrosis (14). Such bioinformatics analyses enhance our understanding of membranous nephropathy, allow for personalized molecular assessments of patients, and identify potential therapeutic targets.

Among various types of cancers, solid cancers originating from the lung, prostate, gastrointestinal tract, and breast, as well as certain hematological cancers, exhibit a closer relationship to MN (5). Some researchers have suggested that T-cell responses, especially those triggered by tumor antigens, may play a significant role in the interaction between cancers and MN (15), but the precise mechanism underlying the occurrence of cancer-associated MN remains unidentified (16). It is worth highlighting that the latest Kidney Disease: Improving Global Outcomes (KDIGO) clinical practice guidelines (17) recommend the use of the monoclonal antibody rituximab, originally employed in lymphoma treatment, as a 1B recommendation for intermediate-to-high-risk MN treatment. Furthermore, newer monoclonal antibodies, such as ocrelizumab, obinutuzumab, and ofatumumab have also been explored for MN treatment (18). As antineoplastic agents assume an increasingly pivotal role in treating immune-mediated non-neoplastic conditions (19, 20), the relationship between immune-related disorders and cancer has garnered more attention. We contend that numerous potential connections between MN and cancer remain undisclosed. Our study aims to employ bioinformatics methods to analyze the potential link between MN and pan-cancer.

2 Methods

2.1 Data source

We conducted a search in the Gene Expression Omnibus (GEO) for gene expression datasets. Firstly, we searched for “membranous nephropathy” as the keyword, selecting “series” for Entry type and “Homo sapiens” for Organisms, and then filtered the single-cell RNA sequencing dataset and bulk RNA sequencing dataset from these datasets. In cases where multiple datasets met the aforementioned criteria, we opted for the dataset with a larger sample size and a greater count of differentially expressed genes (DEGs). In addition, RNA sequencing and clinical data for 33 distinct cancer types were obtained from The Cancer Genome Atlas (TCGA).

2.2 Data processing of the single-cell RNA-seq dataset

The single-cell RNA-seq data underwent filtering and analysis using the Seurat R package. The filtering criteria were established with nFeature_RNA falling within the range of 300 to 7500. To mitigate batch effects among samples, the Harmony R package was employed. The ScaleData function was utilized for data normalization, ensuring zero-centered data for principal component analysis (PCA). For data dimensionality reduction, the RunUMAP function was applied, and the FindAllMarkers function was employed to identify DEGs across distinct clusters. Clustering was carried out at a resolution of 0.8.

2.3 DEG screening

Quantile normalization was conducted using the preprocessCore R package, and DEGs were identified through the utilization of the limma R package. Genes with an adjusted p-value < 0.05 were exclusively considered. The upregulated and downregulated DEGs from the two datasets were subjected to an intersection process to derive the shared upregulated and downregulated DEGs.

2.4 Enrichment analyses and senescence-associated DEGs

We conducted enrichment analyses utilizing the clusterProfiler R package for both Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses. A total of 318 genes associated with senescence were identified by utilizing the cellular senescence pathway (Homo sapiens) within the KEGG database (https://www.kegg.jp/pathway/hsa04218) and the human gene set REACTOME_CELLULAR_SENESCENCE from the REACTOME database (https://gsea-msigdb.org/gsea/msigdb/cards/REACTOME_CELLULAR_SENESCENCE). The intersection between DEGs and senescence-related genes was considered as the set of senescence-related DEGs.

2.5 Protein-protein interaction networks

Constructing protein-protein interaction (PPI) networks for the 13 senescence-related DEGs were constructed using the STRING Database (https://string-db.org/).

2.6 Machine learning and identification of hub genes

we utilized least absolute shrinkage and selection operator (LASSO) regression and random forest to further refine the selection of DEGs. The kernlab R package was employed for LASSO regression, while the randomForest R package was utilized for the random forest algorithm. Subsequently, by intersecting the gene sets derived from LASSO regression and random forest, we identified the central hub genes.

2.7 Characterization and functional analyses of hub genes

We generated receiver operating characteristic (ROC) curves for the hub genes using the pROC R package. Correlation analyses among the hub genes were conducted using the circlize R package. To assess immune cell infiltration, we employed the ssGSEA function from the GSVA R package. Visualization of immune cell correlations, expression differences in immune cells across different groups, and correlations between hub genes and immune cells was achieved using the ggplot2 R package. Furthermore, leveraging the top 50 genes that exhibited positive correlations with the hub genes, we conducted gene set enrichment analysis (GSEA) utilizing the Reactome database. The outcomes were presented through a ridge plot showcasing the 20 most significant pathways.

2.8 Prediction of upstream transcription factors and miRNAs

We utilized the Regnetwork database (https://regnetworkweb.org) to predict the upstream transcription factors (TFs) and microRNAs (miRNAs) associated with the hub genes. The resulting network was visualized using Cytoscape software.

2.9 mRNA differential expression analysis of hub genes in pan-cancer

mRNA sequencing data and clinical information were sourced from the TCGA database. For analysis, we included only 14 cancer types that featured more than ten pairs of cancer and normal samples. These cancer types encompassed BLCA, BRCA, COAD, ESCA, HNSC, KICH, KIRC, KIRP, LIHC, LUAD, LUSC, PRAD, STAD, and THCA. The mRNA expression values from TCGA are presented as normalized RSEM values. The fold change was computed as the ratio of the mean expression in cancer samples to the mean expression in normal samples. To determine statistical significance, a t-test was employed, and the resulting p-values were adjusted using the false discovery rate (FDR).

2.10 Survival analysis of hub genes in pan-caner

The mRNA expression data of hub genes and the corresponding clinical survival data across 33 cancer types were integrated for survival analysis. The cancer samples were categorized into high and low expression groups based on the median gene RSEM value. The survival R package was employed to model the survival time and status for these two groups, encompassing disease-free interval (DFI), disease-specific survival (DSS), overall survival (OS), and disease-free survival (DFS).

2.11 Pan-cancer mutations of hub genes

Copy number variation (CNV) data were collected from the TCGA database across 33 different cancer types. CNV is classified into homozygous and heterozygous types, which include amplifications and deletions. These variations indicate the presence of CNVs on either one or both chromosomes. The frequency (percentage) of single nucleotide variation (SNV) mutations within each gene’s coding region was calculated using the following formula: Number of Mutated Samples/Number of Cancer Samples. An SNV oncoplot was generated using the maftools package. Additionally, CNV profiles were examined, and the correlation between CNV and mRNA expression was assessed. Statistical significance was determined using p-values adjusted through the FDR correction. The relationship between paired mRNA expression and methylation levels was evaluated using Pearson’s product-moment correlation coefficient, followed by a t-distribution test. The resulting p-values were then adjusted using FDR correction. Genes with an adjusted FDR of 0.05 or lower were retained for further analysis.

2.12 Drug sensitivity analysis

The small molecules were sourced from the Genomics of Drug Sensitivity in Cancer (GDSC) and Cancer Therapeutics Response Portal (CTRP) databases. The Spearman correlation coefficient was employed to indicate the potential correlation between gene expression and drug sensitivity. A positive correlation suggests that genes with high expression levels confer resistance to a drug, while a negative correlation suggests that genes with high expression levels render sensitivity to a drug.

2.13 Validation by immunohistochemistry

The experimental procedures conducted in this study were approved by the Ethics Committee of China-Japan Friendship Hospital (Ethics Approval Number: 2021–75-K43). For the control group, normal kidney tissue adjacent to the cancerous tissue removed during kidney cancer surgery was selected. All patients in the MN group were confirmed through pathological examination.

Immunohistochemistry was employed to validate the protein expression levels of the hub genes. Paraffin-embedded tissues were re-sectioned, underwent antigen retrieval, and were then treated with 5% serum to prevent nonspecific binding. Endogenous peroxidase activity was also blocked. Primary antibodies (PI3KR1, CAPN2, TERF2IP from Santa Cruz Biotechnology; CCND1, TXN from Cell Signaling Technology; SLC25A4 from Affinity Biosciences) were applied and incubated at 4°C to evaluate their expression in kidney tissue. Subsequently, secondary antibodies were used, followed by standard incubation, staining, and observation procedures.

3 Results

3.1 Selecting of MN dataset

According to the search conditions, we finally screened and identified a single-cell RNA-seq dataset (GSE171458) and a bulk RNA-seq dataset (GSE108109) for analysis, and a search flow (Supplementary Figure 1) has been drawn up.

3.2 Processing of single-cell RNA-seq data

The single-cell RNA-seq dataset GSE171458 comprises 6 MN patients and 2 healthy subjects (Supplementary Figure 2). Following filtration, a total of 25,223 genes and 14,357 cells were retained. In the MN group, there were 543 genes with up-regulated expression and 193 genes with down-regulated expression (Supplementary Figure 3).

3.3 DEG identification

GSE108109 contains 44 MN patients and 6 healthy subjects, all patients with MN had a clinical presentation of nephrotic syndrome, after quantile normalization of the gene expression matrix (Figures 1A, B), DEGs were identified (Figures 1C, D). The numbers of upregulated genes in GSE171458 and GSE108109 is 2632 and 2389, respectively, while the numbers of downregulated genes were 2632 and 2171, respectively. Intersections of the up- and downregulated DEGs of the two datasets were performed separately to obtain 222 upregulated and 144 downregulated genes (Figures 1E, F).

Figure 1
www.frontiersin.org

Figure 1 Identification of differentially expressed genes. (A, B) Box plot of data before and after quantile normalization from GSE108109; (C) Volcano plot showing the identification of upregulated and downregulated genes. (D) Heatmap displaying the top 20 upregulated and downregulated genes, where upregulated genes are highlighted in light red, and downregulated genes in light blue. (E, F) Venn diagrams illustrating the overlapping upregulated and downregulated DEGs from GSE171458 and GSE108109, revealing 222 and 144 common genes, respectively.

3.4 Enrichment analyses of senescence-associated genes

To comprehensively explore the biological functions and pathways associated with the 366 common DEGs, GO and KEGG pathway enrichment analyses were conducted (Supplementary Figure 4). In the GO enrichment analysis, the top three biological processes exhibiting significant enrichment were cellular cation homeostasis cellular cation homeostasis (p = 5.01E-8), cellular metal ion homeostasis (p = 2.77E-7), and response to metal ion (p = 9.44E-7). In the KEGG enrichment analysis, the top three pathways were protein processing in the endoplasmic reticulum (p = 1.03E-8), amyotrophic lateral sclerosis (p = 0.000453), and Parkinson’s disease (p = 0.000102).

We performed an intersection analysis between the 318 senescence-related genes and the 366 common DEGs, resulting in the identification of 13 senescence-related DEGs (Figure 2A). In the GO enrichment analysis, in terms of biological processes, genes exhibited significant enrichment in functions such as negative regulation of organelle organization (p = 8.61E-5), regulation of DNA binding (p = 6.98E-5), and negative regulation of mitochondrial outer membrane permeabilization involved in the apoptotic signaling pathway (p = 3.39E-5) (Figure 2B). The top three pathways in the KEGG enrichment analysis were cellular senescence (p = 9.64E-17), human T-cell leukemia virus 1 infection (p = 1.41E-7), and necroptosis (p = 3.85E-5) (Figure 2C).

Figure 2
www.frontiersin.org

Figure 2 Enrichment analyses with senescence-related differentially expressed genes. (A) The Venn diagrams show an overlap of 13 genes; (B) GO analysis; (C) Top 5 pathways of KEGG analysis and the genes they contained.

3.5 Protein-protein interaction networks

PPI networks predictions for 13 senescence-related DEGs (Supplementary Figure 5).

3.6 Expression of senescence-related DEGs

The expression patterns of the 13 senescence-related DEGs in GSE108109 are depicted in the volcano plot and heatmap (Supplementary Figures 6A-B). When contrasted with healthy controls, all of the senescence-related DEGs exhibited significant expression differences, except for TERF2IP (Supplementary Figure 6C).

3.7 Identification of hub genes

To further narrow down the genes, we employed machine learning techniques, namely LASSO regression and random forest. Initially, LASSO regression identified seven genes (Figure 3A) and subsequently, the random forest method was used to obtain the top ten genes ranked by their importance (Figure 3B). The intersection of these two sets yielded the final selection of six hub genes (Figure 3C): PIK3R1, CCND1, TERF2IP, SLC25A4, CAPN2, and TXN.

Figure 3
www.frontiersin.org

Figure 3 Identification of hub genes. (A) LASSO regression; (B) Top ten genes ranked by importance in random forest; (C) Venn diagram of LASSO regression and random forest; (D) ROC curve and AUC of hub genes; (E) Correlation analysis between hub genes, red symbolizes positive and green symbolizes negative correlation, and darker color means stronger correlation.

The ROC analysis revealed that the AUC values for PIK3R1 (0.996), SLC25A4 (0.922), CCND1 (1.000), and CAPN2 (0.936) exceeded 0.9, indicating a robust predictive classification capability of these hub genes for distinguishing between MN and healthy controls (Figure 3D). Moreover, visual representations were used to display the correlations among the hub genes (Figure 3E).

3.8 Immune infiltration of hub genes

The correlations among immune cells were intricate and widespread (Supplementary Figure 7A). In comparison to healthy controls, MN patients exhibited significantly elevated levels of gamma delta T cells, macrophages, mast cells, myeloid-derived suppressor cells, monocytes, natural killer cells, natural killer T cells, cytoid dendritic cells, regulatory T cells, follicular helper cells, type 1 T helper cells, and type 2 T helper cells (Supplementary Figure 7B). Additionally, all hub genes showed close associations with immune cell infiltration (Supplementary Figure 7C).

3.9 GSEA of hub genes

In order to gain a deeper insight into the significance of hub genes in MN and to anticipate their potential functions, we initially identified the top 50 genes in the dataset that were positively correlated with the hub genes (Supplementary Figure 8). Subsequently, GSEA was conducted for each individual hub gene. Focal adhesion ranked first in CAPN2 and CCND1, while intrinsic component of organelle membrane, Golgi apparatus, actin cytoskeleton organization, and regulation of anatomical structure morphogenesis ranked first in PIK3R1, SCL25A4, TERF2IP, and TXN, respectively (Supplementary Figure 9).

3.10 Prediction of the upstream regulation network

Using the Regnetwork database, we identified TFs and miRNAs upstream of the six hub genes (Figure 4).

Figure 4
www.frontiersin.org

Figure 4 Network between transcription factors, miRNAs and hub genes. Red indicates hub genes, and blue indicates transcription factors and miRNAs.

3.11 Hub genes expression and survival analysis of hub genes in pan-cancer

Hub genes were significantly differentially expressed in multiple types of cancer. Among them, SLC25A4 was significantly differentially expressed in 11 cancers (highly expressed in KICH, and lowly expressed in HNSC, ESCA, BLCA, STAD, LUSC, KIRP, COAD, PRAD, LUAD, and KIRC), followed by PIK3R1, CAPN2, TXN, TERF2IP and CCND1 in 9, 8, 7, 7, and 6 types of cancer, respectively (Figure 5A). The expression of hub genes correlates with the prognostic indicators (including progression-free internal survival, disease-specific survival, overall survival and disease-free survival) of many types of cancer (Figure 5B). In survival analysis, a Hazard Ratio (HR) less than 1 and a p-value less than 0.05 signify that heightened gene expression diminishes the risk of patient mortality and enhances survival—an advantageous outcome. Conversely, an HR exceeding 1 with a p-value under 0.05 indicates that increased gene expression elevates the hazard of patient demise and diminishes survival—a detrimental scenario. When HR equals 1 and p-value is less than 0.05, it suggests that heightened gene expression has no discernible impact on patient survival.

Figure 5
www.frontiersin.org

Figure 5 Hub gene expression and survival analysis in pan-cancer. (A) mRNA expression profiles of hub genes across diverse cancer types; (B) Assessment of the prognostic significance of hub genes in pan-cancer; DFI, progression free internal; DSS, disease specific survival; OS, overall survival; DFS, disease free survival. In this representation, red indicates elevated gene expression in cancer, blue signifies reduced expression, and a solid circle denotes statistical significance.

3.12 Mutations of hub genes in pan-cancer

First, we analyzed the distribution of CNV types, and the CNV pie chart revealed that the primary CNV types were heterozygous amplification and deletion (Figure 6A). Subsequently, we assessed the mutation frequencies of hub genes in pan-cancer, with PIK3R1 exhibiting the highest frequency of SNV (Figure 6B). The mutation landscape indicated that missense mutations were the predominant type, with the hub genes displaying mutation frequencies in the following descending order: PIK3R1, CAPN2, CCND1, TERD2IP, SLC25A4, and TXN, with mutation percentages of 69%, 17%, 13%, 8%, 6%, and 2%, respectively (Figure 6C).

Figure 6
www.frontiersin.org

Figure 6 Pan-cancer mutations of hub genes. (A) Distribution of copy number variation (CNV) in 33 cancers. This pie chart illustrates the proportion of various CNV types for a single gene in one cancer, with different colors representing distinct CNV types. Hete Amp, heterozygous amplification; Hete Del, heterozygous deletion; Homo Amp, homozygous amplification; Homo Del, homozygous deletion; None, no CNV. (B) Mutation frequency of hub genes. The numbers indicate the count of samples with the corresponding mutated gene in a given cancer. ‘0’ signifies no mutation in the gene coding region, and the absence of a number indicates no mutation in any region of the gene. (C) Single nucleotide variation (SNV) oncoplot. This chart depicts the distribution of mutations in hub genes and categorizes SNV types.

CNV percentage analysis revealed that hub genes displayed heterozygous amplification or deletion in nearly all cancer types, and in some cancers, they also exhibited homozygous amplification or deletion (Figure 7A). Correlation analysis unveiled a close relationship between the mRNA expression of hub genes and pan-cancer CNV. These genes displayed significant positive or negative correlations across multiple cancer types (Figure 7B). These findings suggested that CNV in hub genes mediated their aberrant expression, which could potentially play a crucial role in cancer progression. Methylation and mRNA expression correlation analysis indicated that, for the most part, the expression levels of hub genes were negatively correlated with their methylation levels, especially PIK3R1, CCND1, and CAPN2. Conversely, only TXN in BLCA, KICH, UCS, and UVM, as well as SLC25A4 in CESC, MESO, and TGCT, exhibited a positive correlation between methylation and gene expression (P < 0.05, Figure 7C).

Figure 7
www.frontiersin.org

Figure 7 Correlation analyses of mRNA expression, CNV and drug sensitivity in pan-cancer. (A) Heterozygous and homozygous CNVs in pan-cancer. Genes with > 5% CNV in a given cancer are represented as data points on the figure. (B) Correlation between CNV and mRNA expression. (C) Correlation between methylation and mRNA gene expression. (D, E) Gene set drug resistance analysis from Genomics of Drug Sensitivity in Cancer (GDSC) and Cancer Therapeutics Response Portal (CTRP) (Top 30). The size of the data points indicates the statistical significance, with larger dots indicating higher statistical significance. The false discovery rate (FDR) was used for correction.

Genomic aberrations have a significant impact on the clinical response to both chemotherapy and targeted therapy treatments. To investigate the role of hub genes in chemotherapy and targeted therapy, we analyzed drug sensitivity and gene expression profiling data from cancer cell lines in GDSC and CTRP. The correlation analysis specifically refers to the relationship between gene expression and the half-maximal inhibitory concentration (IC50) of a drug, which is usually used to assess antitumor activity. A lower IC50 value indicates greater drug potency. Thus, a positive correlation implies that higher gene expression weakens the drug’s inhibitory effect, while a negative correlation suggests that higher gene expression strengthens it. The mRNA expression level of TXN, CCND1, and CAPN2 showed a positive correlation with the sensitivity to most drugs, except for 17-AAG and docetaxel. On the other hand, the mRNA expression level of TERF2IP and PIK3R1 exhibited a negative correlation with the sensitivity to most drugs in GDSC and CTRP, again with the exceptions of 17-AAG and docetaxel (Figures 7D, E).

3.13 Immunohistochemical verification of hub gene expression in the kidneys

We performed immunohistochemical staining on 22 MN kidney specimens and 3 control kidney specimens (Figure 8). In the control group, PIK3R1 and SLC25A4 were expressed in both glomeruli and tubules, while CCND1 and TXN were primarily expressed in tubules, with minimal expression in glomeruli. CAPN2 was expressed in glomeruli and in some renal tubular epithelial cells. In contrast, in the MN group, PIK3R1 exhibited reduced expression in glomeruli, while TXN showed significantly decreased expression in tubules. CCND1 and CAPN2 displayed markedly increased expression in glomerular podocytes, and SLC25A4 exhibited increased expression along the glomerular basement membrane. TERF2IP, on the other hand, showed minimal expression in renal tissues in both the control and MN groups.

Figure 8
www.frontiersin.org

Figure 8 Immunohistochemical Validation of Hub Gene Expression in Human Kidney Specimens (400×). (A, B) PIK3R1, (C, D) CCND1, (E, F) TERF2IP, (G, H) SLC25A4, (I, J) CAPN2, (K, L) TXN. (A, C, E, G, I, K) are normal kidney tissues, while (B, D, F, H, J, L) are membranous nephropathy kidney tissues.

4 Discussion

In recent years, an increasing body of research has highlighted a potential correlation between immune-mediated diseases and cancer (21, 22). Possible mechanisms involved include: 1. Cancer-associated antigens induce the host to produce antibodies, and these antigens may become lodged beneath glomerular epithelial cells, forming in situ immune complexes that ultimately mediate renal injury (23). 2. Certain oncogenic viruses, such as hepatitis B virus, cytomegalovirus, and Epstein-Barr virus, infect the host, resulting in the production of antibodies against viral antigens. These viral antigen-antibody complexes may then deposit in the glomerulus, activate the complement system, and cause renal injury (24). 3. Necrotic tumors release substantial amounts of DNA, which triggers the body to produce anti-DNA antibodies. These immune complexes may accumulate in the kidneys, subsequently causing kidney damage (25).

In this study, we aimed to explore the potential connection between MN and pan-cancer. Our analysis led to the identification of PIK3R1, CCND1, TERF2IP, SLC25A4, CAPN2, and TXN as hub genes, which could potentially serve as targets for interventions in both MN and cancer.

The PI3K enzymes constitute a conserved family of lipid kinases, consisting of a catalytic subunit and a regulatory subunit. Phosphatidylinositol-3-kinase regulatory subunit 1 (PIK3R1) exhibits low expression levels in the majority of cancers and is believed to function as a cancer suppressor. Downregulation of PIK3R1 is associated with poor survival outcomes for most cancer patients (26). Conversely, PIK3R1 has been observed to be significantly upregulated in rats with adriamycin-induced chronic glomerulonephritis (27). Simultaneously, phospholipase A2 receptor (PLA2R), the primary target antigen in MN, has been demonstrated to activate the upregulated PI3K/AKT/mTOR pathway (28).

CCND1 plays a pivotal role as a critical regulator of the cell cycle and holds a central position in the development of cancer by driving uncontrolled cellular proliferation. Its activity is significantly heightened in various cancer contexts (29), and the expression of CCND1 is indispensable for the survival and proliferation of cancer cells (30). In glomerular intrinsic cells, CCND1 was found to be expressed in podocytes in both the Heymann nephritis model of rats (31) and cases of FSGS in humans (32). It is primarily enriched in actively proliferating podocytes, as opposed to those in a quiescent state, and its expression increases following injury in passive Heymann nephritis rats (31).

TERF2IP, the most highly conserved component of the shelterin complex, plays a multifaceted role in the regulation of various cellular processes, encompassing cell metabolism, DNA damage response, and NF-κB signaling (33). It has been shown to play a part in oncogenesis, cancer progression, and the development of resistance to chemotherapy in human cancers, with multiple mutations and diverse expression patterns of TERF2IP reported in cancer contexts (34, 35). On the other hand, TERF2IP serves as a central signaling hub within podocytes (36). In murine disease models and kidney biopsies from glomerulosclerosis patients, injured podocytes displayed reduced activation of TERF2IP within the glomeruli. Notably, severe glomerulosclerosis manifests in mice with diminished podocyte expression of TERF2IP, leading to early mortality from renal failure by 8 weeks of age. Furthermore, podocyte-specific TERF2IP haploinsufficiency also resulted in significant podocyte damage, including signs of podocyte detachment (37).

The solute carrier protein 25 (SLC25) family, which is the largest gene transporter family, consists of membrane proteins that regulate the transport of various solutes in and out of cells. They play crucial roles in essential physiological processes such as cellular material transport, energy transmission, signal transduction, and nutrient metabolism. Pan-cancer analysis of the SLC25 family suggests that SLC25A4 is linked to multiple oncogenic pathways, including the PI3K-AKT-MTOR pathway, MYC-TARGETS-V1 pathway, MYC-TARGETS-V2 pathway, and MTORC1 pathway (38). Regrettably, there have been no reported associations between SLC25A4 and membranous nephropathy.

CAPN2 (calpain-2) is a prototypical classical isoform of the calpain family of calcium-activated cysteine proteases. Its substrate proteins are involved in a wide range of cellular processes, including transcription, survival, proliferation, apoptosis, migration, and invasion. Dysregulated calpain activity has been linked to tumorigenesis, suggesting that calpains may hold promise as therapeutic targets (39). Interestingly, researchers have observed that inhibiting both calpain 1 and 2 in cell cycle protein G-related kinase knockout mice mitigated podocyte injury. This finding establishes a direct correlation between calpain-1/-2 activity and podocyte injury, proteinuria, and glomerulosclerosis (40).

TXN (thioredoxin-1) is a multifunctional protein with a molecular weight of 12 kDa, primarily localized within the cytosol. TXN plays a pivotal role in a diverse range of cellular functions, including cell proliferation, the maintenance of redox homeostasis, DNA synthesis, gene expression regulation, and the regulation of apoptosis-mediated cell death. TXN is indispensable for the normal functioning of both organs and tumors. It is strongly associated with various diseases, notably cancer, and ample evidence has been presented to underscore its significance in influencing the phenotype and prognosis of lung, gastrointestinal, and urological cancers (41). On the other hand, urinary TXN is regarded as a biomarker for diagnosing tubular redox dysregulation (42). Furthermore, recombinant long-acting TXN has been shown to ameliorate the transition from acute kidney injury to chronic kidney disease by modulating renal oxidative stress and inflammation (43).

Although a causal association between lung cancer and MN was not found in the Mendelian randomization-based study by Yang et al. (44), we don’t think they are the opposite of our conclusions. We speculate that the negative result may be related to the fact that the investigators used the primary MN dataset for their analyses. In recent years, it has become increasingly evident that the development of MN is associated with a variety of target antigens present on podocytes (45). Among the six central genes identified through bioinformatics analysis, subsequent immunostaining of renal tissues revealed a notable increase in the expression of CCND1 and CAPN2 within glomerular podocytes in the MN group compared to the control group. Moreover, these two genes have been substantiated to have significant implications in podocyte injury (32, 40). We believe that further exploration of these two genes in our research endeavors may unveil novel and intriguing findings in the future.

Certainly, this study has limitations. First, MN is the most prevalent pathologic type of cancer-associated nephropathy, though other types such as minimal change disease, focal segmental glomerulosclerosis, IgA nephropathy, and membranoproliferative glomerulonephritis have also been documented (5). Our study focused exclusively on MN, and thus it remains unclear whether the six hub genes identified are specific to MN. The intricate relationship between kidney disease and cancers warrants further investigation. Second, while elevated expression of CCND1 and CAPD2 has been noted in MN glomerular podocytes, their roles in other podocytopathies like minimal change disease are yet to be defined. Whether these genes could serve as specific markers of podocyte damage in MN or their association with MN target antigens (e.g., PLA2R, NELL1, THSD7A) require more extensive research. Third, SLC25A4, and TXN have not been previously reported in the context of MN or podocyte injury. Investigating the mechanisms underlying their roles in MN requires more research.

5 Conclusions

In conclusion, we identified PIK3R1, CCND1, TERF2IP, SLC25A4, CAPN2, and TXN as potential markers associated with both cancer and MN. This discovery enhances our understanding of the potential connection between cancers and MN. Furthermore, these genes represent potential therapeutic targets for MN as well as various types of cancers.

Data availability statement

Publicly available datasets were analyzed in this study. This data can be found here: Gene Expression Omnibus at https://www.ncbi.nlm.nih.gov/geo/, reference number: GSE171458 and GSE108109.

Ethics statement

The studies involving humans were approved by the Ethics Committee of China-Japan Friendship Hospital (Ethics Approval Number: 2021-75-K43). The studies were conducted in accordance with the local legislation and institutional requirements. The participants provided their written informed consent to participate in this study.

Author contributions

YY: Data curation, Methodology, Validation, Writing – original draft. GZ: Methodology, Validation, Writing – original draft. XW: Funding acquisition, Methodology, Software, Writing – original draft. ZZ: Data curation, Methodology, Writing – original draft. LZ: Writing – original draft. QX: Writing – review & editing. WL: Funding acquisition, Project administration, Writing – review & editing.

Funding

The author(s) declare financial support was received for the research, authorship, and/or publication of this article. This research was funded by National High Level Hospital Clinical Research Funding (2023-NHLHCRF-YS-01) and College-level Projects in China-Japan Friendship Hospital (2017-2-QN-19).

Acknowledgments

We gratefully acknowledge Dr. Ji Wang for her invaluable assistance with the immunohistochemistry experiments.

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.2024.1302909/full#supplementary-material

References

1. Ronco P, Beck L, Debiec H, Fervenza FC, Hou FF, Jha V, et al. Membranous nephropathy. Nat Rev Dis Primers. (2021) 7:69. doi: 10.1038/s41572-021-00303-z

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Beck LH Jr., Salant DJ. Membranous nephropathy: from models to man. J Clin Invest. (2014) 124:2307–14. doi: 10.1172/JCI72270

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Lee JC, Yamauchi H, Hopper J Jr. The association of cancer and the nephrotic syndrome. Ann Intern Med. (1966) 64:41–51. doi: 10.7326/0003-4819-64-1-41

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Rosner MH, Jhaveri KD, McMahon BA, Perazella MA. Onconephrology: The intersections between the kidney and cancer. CA Cancer J Clin. (2021) 71:47–77. doi: 10.3322/caac.21636

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Leeaphorn N, Kue A, Pai P, Thamcharoen N, Ungprasert P, Stokes MB, et al. Prevalence of cancer in membranous nephropathy: a systematic review and meta-analysis of observational studies. Am J Nephrol. (2014) 40:29–35. doi: 10.1159/000364782

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Plaisier E, Ronco P. Screening for cancer in patients with glomerular diseases. Clin J Am Soc Nephrol: CJASN. (2020) 15:886–8. doi: 10.2215/CJN.09000819

CrossRef Full Text | Google Scholar

7. Hephzibah Cathryn R, Udhaya Kumar S, Younes S, Zayed H, George Priya Doss C. A review of bioinformatics tools and web servers in different microarray platforms used in cancer research. Adv Protein Chem Struct Biol. (2022) 131:85–164. doi: 10.1016/bs.apcsb.2022.05.002

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Balasundaram A, Udhaya Kumar S, George Priya Doss C. A computational model revealing the immune-related hub genes and key pathways involved in rheumatoid arthritis (RA). Adv Protein Chem Struct Biol. (2022) 129:247–73. doi: 10.1016/bs.apcsb.2021.11.006

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Kumar SU, Kumar DT, Siva R, Doss CGP, Zayed H. Integrative bioinformatics approaches to map potential novel genes and pathways involved in ovarian cancer. Front Bioengineering Biotechnol. (2019) 7:391. doi: 10.3389/fbioe.2019.00391

CrossRef Full Text | Google Scholar

10. Zhang P, Geng Y, Tang J, Cao Z, Xiang X, Yang K, et al. Identification of biomarkers related to immune and inflammation in membranous nephropathy: comprehensive bioinformatic analysis and validation. Front Immunol. (2023) 14:1252347. doi: 10.3389/fimmu.2023.1252347

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Han M, Wang Y, Huang X, Li P, Liang X, Wang R, et al. Identification of hub genes and their correlation with immune infiltrating cells in membranous nephropathy: an integrated bioinformatics analysis. Eur J Med Res. (2023) 28:525. doi: 10.1186/s40001-023-01311-3

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Han M, Wang Y, Huang X, Li P, Shan W, Gu H, et al. Prediction of biomarkers associated with membranous nephropathy: Bioinformatic analysis and experimental validation. Int Immunopharmacol. (2024) 126:111266. doi: 10.1016/j.intimp.2023.111266

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Cai XY, Wang ZF, Ge SW, Xu G. Identification of hub genes and immune-related pathways for membranous nephropathy by bioinformatics analysis. Front Physiol. (2022) 13:914382. doi: 10.3389/fphys.2022.914382

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Hou Y, Li Y, Wang Y, Li W, Xiao Z. Screening and analysis of key genes in miRNA-mRNA regulatory network of membranous nephropathy. J Healthcare Engineering. (2021) 2021:5331948. doi: 10.1155/2021/5331948

CrossRef Full Text | Google Scholar

15. Holdsworth SR, Kitching AR, Tipping PG. Th1 and Th2 T helper cell subsets affect patterns of injury and outcomes in glomerulonephritis. Kidney Int. (1999) 55:1198–216. doi: 10.1046/j.1523-1755.1999.00369.x

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Hoxha E, Wiech T, Stahl PR, Zahner G, Tomas NM, Meyer-Schwesinger C, et al. A mechanism for cancer-associated membranous nephropathy. N Engl J Med. (2016) 374:1995–6. doi: 10.1056/NEJMc1511702

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Kidney Disease: Improving Global Outcomes Glomerular Diseases Work Group. KDIGO 2021 clinical practice guideline for the management of glomerular diseases. Kidney Int. (2021) 100:S1–S276. doi: 10.1016/j.kint.2021.05.021

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Deng L, Xu G. Update on the application of monoclonal antibody therapy in primary membranous nephropathy. Drugs. (2023) 83:507–30. doi: 10.1007/s40265-023-01855-y

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Kaegi C, Wuest B, Crowley C, Boyman O. Systematic review of safety and efficacy of second- and third-generation CD20-targeting biologics in treating immune-mediated disorders. Front Immunol. (2021) 12:788830. doi: 10.3389/fimmu.2021.788830

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Wechsler ME, Ruddy MK, Pavord ID, Israel E, Rabe KF, Ford LB, et al. Efficacy and safety of itepekimab in patients with moderate-to-severe asthma. N Engl J Med. (2021) 385:1656–68. doi: 10.1056/NEJMoa2024257

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Zhao Z, He S, Yu X, Lai X, Tang S, Mariya M, et al. Analysis and experimental validation of rheumatoid arthritis innate immunity gene CYFIP2 and pan-cancer. Front Immunol. (2022) 13:954848. doi: 10.3389/fimmu.2022.954848

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Zhao S, Wu Y, Wei Y, Xu X, Zheng J. Identification of biomarkers associated with CD8+ T cells in coronary artery disease and their pan-cancer analysis. Front Immunol. (2022) 13:876616. doi: 10.3389/fimmu.2022.876616

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Pascal RR, Slovin SF. Tumor directed antibody and carcinoembryonic antigen in the glomeruli of a patient with gastric carcinoma. Hum Pathol. (1980) 11:679–82. doi: 10.1016/S0046-8177(80)80082-7

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Birkeland SA, Storm HH. Glomerulonephritis and Malignancy: a population-based analysis. Kidney Int. (2003) 63:716–21. doi: 10.1046/j.1523-1755.2003.00771.x

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Higgins MR, Randall RE, Still WJ. Nephrotic syndrome with oat-cell carcinoma. Br Med J. (1974) 3:450–1. doi: 10.1136/bmj.3.5928.450

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Liu Y, Wang D, Li Z, Li X, Jin M, Jia N, et al. Pan-cancer analysis on the role of PIK3R1 and PIK3R2 in human tumors. Sci Rep. (2022) 12:5924. doi: 10.1038/s41598-022-09889-0

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Gao JR, Qin XJ, Jiang H, Wang T, Song JM, Xu SZ. Screening and functional analysis of differentially expressed genes in chronic glomerulonephritis by whole genome microarray. Gene. (2016) 589:72–80. doi: 10.1016/j.gene.2016.05.031

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Chiou TT, Chau YY, Chen JB, Hsu HH, Hung SP, Lee WC. Rapamycin attenuates PLA2R activation-mediated podocyte apoptosis via the PI3K/AKT/mTOR pathway. BioMed Pharmacother. (2021) 144:112349. doi: 10.1016/j.biopha.2021.112349

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Montalto FI, De Amicis F. Cyclin D1 in cancer: A molecular connection for cell cycle control, adhesion and invasion in tumor and stroma. Cells. (2020) 9(12):2648. doi: 10.3390/cells9122648

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Laphanuwat P, Likasitwatanakul P, Sittithumcharee G, Thaphaengphan A, Chomanee N, Suppramote O, et al. Cyclin D1 depletion interferes with oxidative balance and promotes cancer cell senescence. J Cell Sci. (2018) 131(12):jcs214726. doi: 10.1242/jcs.214726

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Petermann A, Hiromura K, Pippin J, Blonski M, Couser WG, Kopp J, et al. Differential expression of d-type cyclins in podocytes in vitro and in vivo. Am J Pathol. (2004) 164:1417–24. doi: 10.1016/S0002-9440(10)63228-2

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Wang S, Kim JH, Moon KC, Hong HK, Lee HS. Cell-cycle mechanisms involved in podocyte proliferation in cellular lesion of focal segmental glomerulosclerosis. Am J Kidney Dis. (2004) 43:19–27. doi: 10.1053/j.ajkd.2003.09.010

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Deregowska A, Wnuk M. RAP1/TERF2IP-A multifunctional player in cancer development. Cancers (Basel). (2021) 13(23):5970. doi: 10.3390/cancers13235970

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Kotla S, Vu HT, Ko KA, Wang Y, Imanishi M, Heo KS, et al. Endothelial senescence is induced by phosphorylation and nuclear export of telomeric repeat binding factor 2-interacting protein. JCI Insight. (2019) 4(9):e124867. doi: 10.1172/jci.insight.124867

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Anuja K, Kar M, Chowdhury AR, Shankar G, Padhi S, Roy S, et al. Role of telomeric RAP1 in radiation sensitivity modulation and its interaction with CSC marker KLF4 in colorectal cancer. Int J Radiat Biol. (2020) 96:790–802. doi: 10.1080/09553002.2020.1721609

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Zhu B, Cao A, Li J, Young J, Wong J, Ashraf S, et al. Disruption of MAGI2-RapGEF2-Rap1 signaling contributes to podocyte dysfunction in congenital nephrotic syndrome caused by mutations in MAGI2. Kidney Int. (2019) 96:642–55. doi: 10.1016/j.kint.2019.03.016

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Potla U, Ni J, Vadaparampil J, Yang G, Leventhal JS, Campbell KN, et al. Podocyte-specific RAP1GAP expression contributes to focal segmental glomerulosclerosis-associated glomerular injury. J Clin Invest. (2014) 124:1757–69. doi: 10.1172/JCI67846

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Liu AR, Liu YN, Shen SX, Yan LR, Lv Z, Ding HX, et al. Comprehensive analysis and validation of solute carrier family 25 (SLC25) and its correlation with immune infiltration in pan-cancer. BioMed Res Int. (2022) 2022:4009354. doi: 10.1155/2022/4009354

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Shapovalov I, Harper D, Greer PA. Calpain as a therapeutic target in cancer. Expert Opin Ther Targets. (2022) 26:217–31. doi: 10.1080/14728222.2022.2047178

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Tian X, Inoue K, Zhang Y, Wang Y, Sperati CJ, Pedigo CE, et al. Inhibiting calpain 1 and 2 in cyclin G associated kinase-knockout mice mitigates podocyte injury. JCI Insight. (2020) 5(22):e142740. doi: 10.1172/jci.insight.142740

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Oberacker T, Kraft L, Schanz M, Latus J, Schricker S. The importance of thioredoxin-1 in health and disease. Antioxid (Basel). (2023) 12(5):1078. doi: 10.3390/antiox12051078

CrossRef Full Text | Google Scholar

42. Kasuno K, Yodoi J, Iwano M. Urinary thioredoxin as a biomarker of renal redox dysregulation and a companion diagnostic to identify responders to redox-modulating therapeutics. Antioxid Redox Signal. (2022) 36:1051–65. doi: 10.1089/ars.2021.0194

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Nishida K, Watanabe H, Murata R, Tokumaru K, Fujimura R, Oshiro S, et al. Recombinant long-acting thioredoxin ameliorates AKI to CKD transition via modulating renal oxidative stress and inflammation. Int J Mol Sci. (2021) 22(11):5600. doi: 10.3390/ijms22115600

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Yang K, Ding X, Liu J, Liu S, Liu Q, Li J, et al. Two-sample mendelian randomization reveals a causal association between membranous nephropathy and lung cancer. Commun Biol. (2023) 6:887. doi: 10.1038/s42003-023-05111-7

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Sethi S, Madden B. Mapping antigens of membranous nephropathy: almost there. Kidney Int. (2023) 103:469–72. doi: 10.1016/j.kint.2023.01.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: membranous nephropathy, pan-cancer, differential expression gene, potential link, biomarkers

Citation: Yang Y, Zou G-m, Wei X-s, Zhang Z, Zhuo L, Xu Q-q and Li W-g (2024) Identification and validation of biomarkers in membranous nephropathy and pan-cancer analysis. Front. Immunol. 15:1302909. doi: 10.3389/fimmu.2024.1302909

Received: 27 September 2023; Accepted: 10 May 2024;
Published: 23 May 2024.

Edited by:

Mario Ollero, INSERM U955 Institut Mondor de Recherche Biomédicale (IMRB), France

Reviewed by:

Shuwang Ge, Huazhong University of Science and Technology, China
Udhaya Kumar, Baylor College of Medicine, United States

Copyright © 2024 Yang, Zou, Wei, Zhang, Zhuo, Xu and Li. 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: Qian-qian Xu, xqq816926@163.com; Wen-ge Li, wenge_lee2002@126.com

These authors have contributed equally to this work and share last authorship

Disclaimer: 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.