Expression and clinical significance of RAG1 in myelodysplastic syndromes

Objective To determine the expression level of RAG1 and its clinical significance in myelodysplastic syndromes (MDS). Methods To explore the candidate genes, the microarray datasets GSE19429, GSE58831, and GSE2779 were downloaded from the Gene Expression Omnibus database. The differentially expressed genes (DEGs) in MDS were screened using RStudio, and overlapped DEGs were obtained with Venn Diagrams. Gene Ontology and Kyoto Encyclopedia of Genes and Genomes pathway enrichment analyses, and protein–protein interaction network were performed. Quantitative real-time PCR (qRT-PCR) was employed to confirm the microarray results. Results This study identified 26 DEGs. Functional enrichment analyses indicated that these DEGs were significantly enriched in the immune response, and hematopoietic cell lineage. Eight core genes, for example, RAG1 and PAX5, were identified with a high degree of connectivity. The result of qRT-PCR showed that RAG1 was significantly down-regulated in MDS patients, which helped in distinguishing MDS patients from normal controls. The area under the curve of the receiver operator characteristic was 0.913 (P < 0.0001). MDS patients with low RAG1 expression level had a poor long-term survival (P = 0.031). What's more, the expression of RAG1 was significantly increased in the patients who received treatment. Conclusion The results showed that the expression of RAG1 was down-regulated in MDS patients. Lower RAG1 expression was associated with adverse clinical outcomes. RAG1 may be a potential prognostic biomarker for MDS.


Introduction
Myelodysplastic syndromes (MDS) are a kind of malignant heterogeneous clonal hematological diseases originating from bone marrow hematopoietic stem cells, which are characterized by morphological abnormalities of myeloid cells or multiple myeloid cells, ineffective hematopoiesis. There is a high risk of transformation to acute myeloid leukemia (AML) [1,2]. The pathogenesis of MDS is very complex, involving abnormal gene expression, gene mutations, epigenetic abnormalities, defects in immune mechanism, and changes in bone marrow microenvironment [1,3]. With the growth of geriatric population, the incidence of MDS has been increasing, but its clinical treatment is difficult and the cure rate is low. MDS has seriously threatened human health, and there is an urgent need to develop new therapeutic strategies to ameliorate disorders of MDS [4].
The plasticity of the acquired immune system in recognizing millions of possible antigens is largely due to the combinatorial joining of variable (V), diversity (D), and joining (J) gene segments that encode the antigen-binding regions of T cell receptors (TCRs) in T cells and B cell receptors (BCRs) in B cells, and the junctional diversity that can be introduced during the process of V(D)J recombination [5]. RAG1 and RAG2 proteins form a complex and initiate V(D)J recombination by introducing DNA double-strand breaks (DSBs) between the recombination signal sequences and the flanking V, D, or J gene segment. The human RAG1 protein consists of 1,043 amino acids, and the catalytic core (amino acids 387-1011) contains a nonamer-binding domain, a dimerization and DNA binding domain, a pre-RNase H and catalytic RNase H domain, 2 zinc-binding domains, and the carboxyterminal domain, which are all crucial for V(D)J recombination [6,7].
RAGs also can cleave cryptic RSSs (cRSSs) that are abundant in the genome. Various lines of evidence reveal that RAGs can cleave at non-B DNA structures resulting in cell death or chromosomal translocations and cancer, including leukemia and lymphoma [8]. Gazumyan, et al. reported that DNA breaks caused by RAG1 and activation-induced cytidine deaminase (AID) induce c-myc/immunoglobulin (Ig) heavy chain chromosomal translocations and thereby stimulated mice lymphomagenesis [9]. Mice, lacking the RAG1 protein, showed signs of earlier transformation to acute leukemia and reduced survival rates [10]. Hunter ZR, et al. performed next-generation transcriptional profiling in 57 Waldenström macroglobulinemia (WM) patients, compared with healthy donors, WM patient samples showed greatly enhanced expression of RAG1, its aberrant expression was associated with specific somatic mutation patterns such as deletions in BTG1 and ETV6 [11]. Ting Kang, et al. showed that RAG1 was down-regulated in gastric cancer, low RAG1 expression correlated with poor survival of gastric cancer patients [12]. Harriet E Gee, et al. proved that low RAG1 expression significantly correlated with local recurrence of Early Breast Cancer [13]. These studies suggest that RAG1 may play an important role in the occurrence and development of cancer. However, the effects of RAG1 on the treatment and prognosis of MDS patients have not been extensively studied.
In our study, three datasets containing transcriptome sequencing data from 356 MDS patients and 43 healthy controls from the Gene Expression Omnibus (GEO) database, 59 MDS bone marrow (BM) specimens, and 12 BM specimens of patients with MDS who progressed to secondary acute myeloid leukemia (sAML) from our department, were used to explore the diagnosis and prognostic value of RAG1.

GEO database
Gene chip data is downloaded from the GEO database based on GPL570 [HG-U133_Plus_2] gene expression profile chip data set of Affymetrix human genome U133A array platform: GSE19429 and GSE58831. The platform of GSE2779 is GPL96 [HG-U133A] Affymetrix human genome U133A array. These three datasets determined the microarray datasets of the CD34 + cells from MDS patients and normal samples, including 356 MDS patients and 43 healthy controls. All data are freely available online.

Screening and identification of differentially expressed genes (DEGs)
We employed the R package 'limma' to identify the significant DEGs between MDS patients and normal samples. We set fold change = 1 and pad j < 0.05 to find differentially expressed genes in MDS patients and normal control samples. The Venn diagram tool helped identify the intersecting part (http:// bioinformatics.psb.ugent.be/webtools/Venn/).

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses of DEGs
The Database for Annotation, Visualization, and Integrated Discovery (https://david.ncifcrf.gov/, DAVID), as an online biological information database, contributes to the enrichment analysis to investigate DEGs. The GO and KEGG pathway enrichment analyses of DEGs were performed by the database for annotation. The results were considered statistically significant, when the p-values were less than 0.05.

Construction of protein-protein interaction (PPI) network and module analyses
The STRING database was used to build the PPI network specific to DEGs. Cytoscape is a public software program that can graphically display, analyze and edit the network. Cytohubba and MCODE which are plugins in the Cytoscape software (Version 3.9.0), were applied to calculate the degree of every protein node, and modular analysis of PPI networks, respectively.  Table 4. All the participants had been given informed consent according to the Declaration of Helsinki. The study was approved by the Human Ethics Committees Review Board at Guangxi Medical University, Nanning, China (Approval Number: 2022-KY-E-062).

Statistical analysis
All statistical analyses were performed using SPSS (version 25.0; SPSS Inc, Chicago, IL). qRT-PCR validation was based on Mann-Whitney's U test. Mann-Whitney's U test and Pearson Chi-square test were used to evaluate the relationship between RAG1 and clinical parameters. Receiver Operator Characteristic (ROC) analysis was used to detect sensitivity, specificity, and the area under the curve (AUC) for RAG1. Survival curves were plotted using the Kaplan-Meier method and compared with log-rank tests. P < 0.05 is considered a significant difference.  (Figure 1). At last, 26 overlapping DEGs were obtained, including PAX5, and RAG1, all of which were significantly down-regulated in MDS.

GO and KEGG pathway enrichment analyses
The results of the GO and KEGG pathways were analyzed based on 26 overlapping genes. GO term enrichment analysis showed that the main biological processes involved by DEGs in MDS were: immunoglobulin production, immune response, cellular response to cytokine stimulus, chondroitin sulfate biosynthetic process, B cell proliferation, and external side of plasma membrane (Table 2). KEGG enrichment analysis demonstrated that these DEGs were significantly enriched in the Hematopoietic cell lineage and Glycosaminoglycan biosynthesis-chondroitin sulfate/dermatan sulfate (Table 3).

Quantitative real-time PCR confirmed the microarray results
Based on bioinformatics analysis, considering the value of RAG1, PAX5 in MDS remain undefined, we investigated the expression levels of these genes in BM-MNCs in MDS. RAG1 and PAX5 expressions were validated in BM of 20 MDS patients and 10 normal controls using qRT-PCR. As displayed in (Figure 3), the expression of RAG1 was significantly down-regulated in MDS compared with normal controls. There was no significant difference in the expression levels of PAX5 between patients and normal controls. Next, the authors expanded the sample size, and focused on the analysis of RAG1.

RAG1 expression of newly diagnosed MDS patients
The expression of RAG1 was validated by qRT-PCR in BM from 59 MDS patients and 29 normal controls. As displayed in (Figure 4), the expression of RAG1 was significantly down-regulated in each group of MDS patients compared with normal controls.

Association RAG1 expression with clinical characteristics of MDS
To analyze the relationship between the expression of RAG1 and MDS clinical features, the 59 MDS patients were divided into two groups according to the median RAG1 expression, in which the patients with RAG1 expression value below the median level (0.128) were divided into the low expression group (n = 30), while the RAG1 level was equal to or above the median level were assigned to the high expression group (n = 29). RAG1 expression and clinical features were summarized in Table 4.
To evaluate the value of RAG1 in diagnosis and prognosis of MDS. We used ROC curve analysis to analyze the diagnostic accuracy of RAG1. The AUC for RAG1 was as follows: 0.913 (95% CI 0.846-0.980, P < 0.0001) ( Figure 5A). This result indicated that RAG1 could discriminate MDS patients from normal subjects. To investigate the prognostic implications of RAG1 in MDS, we estimated the relationship between RAG1 and overall survival (OS). Using the Kaplan-Meier analysis, we discovered that the high expression group was significantly different from the low expression group in OS (P = 0.031, (Figure 5B)), which suggested that RAG1 expression was associated with the survival of MDS patients. Multivariate analysis also showed that low expression of RAG1 was an independent prognostic factor for OS. (P = 0.045, Table 5).

The expression of RAG1 in MDS that transformed into AML
To investigate whether the expression of RAG1 was changed during the disease progression of MDS patients, we examined the qRT-PCR result in 59 MDS patients and 12 sAML patients. Our data showed that RAG1 expression level was higher in sAML patients compared to MDS patients, but there was no statistically significant.

Changes of RAG1 expression in patients with MDS and sAML before and after treatment
In order to evaluate the dynamic change trend of RAG1 in the progression of MDS disease, we collected bone marrow samples from 19 patients with MDS and 3 patients with sAML, both groups received treatment and had the expression of RAG1 measured again. Compared with the expression of RAG1 before treatment, the results showed that the change of RAG1 expression in patients with MDS after treatment was generally upregulated (P < 0.0001, (Figure 6A)). And the change of RAG1 expression in the sAML patients after treatment was down-regulated (P > 0.05). Furthermore, among 14 patients who achieved complete response (CR) after treatment, the expression level of RAG1 has returned to the original level, although it was lower than that of the normal controls, there was no statistical difference (P > 0.05). However, the expression of RAG1 in these patients with non-remission (NR) after treatment were significantly lower than that in normal controls (P = 0.014, (Figure 6B)). The results revealed that RAG1 expression level in MDS patients

Discussion
V(D)J recombination is an important mechanism by which the immune system recognizes various antigens existing in nature and exerts normal immune function to protect the body. It plays an important role in the formation of the genetic diversity of the immune system. As one of the recombinases, the RAG1 protein participates in this process from the initial stage to ensure the smooth and normal progress of V(D)J recombination [14,15]. The disorder or lack of RAG1 function can cause severe immune function  deficiency [16,17] and grow into clinically apparent cancers. In recent years, the role of RAG1 in tumors has attracted attention. For instance, RAG1 plays an important role in hepatocellular carcinoma (HCC) development. Transarterial chemoembolization (TACE) treated HCC patients with low RAG1 expression showed shorter OS compared to the Deletion.Deletion genotype carriers and RAG1 may impact TACE efficiency via regulating antitumor immunity in HCC. RAG1 may be a target to improve therapeutic strategy for unresectable HCC patients [18]. Accumulated evidence suggested that the expression levels of RAG1 vary among different tumors. RAG1 was high expressed in lung cancer and may play an important role in cancer growth and metastasis [19]. RAG1 was up-regulated in adult B-cell acute lymphoblastic leukemia. High expression of RAG1 was observed to be correlated with the high proliferation marker [20]. RAG1 was abnormally high expressed in chronic myelogenous leukemia (CML) from a chronic phase (CP) to a blast crisis of B lymphoid (CML-LBC) and BCR-ABL1 acute lymphoblastic leukemia (BCR-ABL1 ALL). RAG1mediated recombination in the formation of previously identified and novel oncogenic rearrangement in CML, was associated with poor survival [21,22]. RAG1 inhibited apoptosis through DSB-induced up-regulation of antiapoptotic genes BCL-2, BCL-XL, and MCL-1 and down-regulation of proapoptotic gene BAD expression in BCR-ABL1 cells [22]. Li Juanjuan, et al. reported a rare case of acquired immunodeficiency associated with myelodysplastic syndromes. The results showed that a homozygous RAG1 mutation was detected in his bone marrow and blood, and no   other pathogenic driver mutations in MDS-related genes were found [23]. However, the expression of RAG1 and its prognostic value in MDS has been reported rarely. Myelodysplastic syndromes are hematopoietic stem cell (HSC) disorders defined by myeloid cell dysplasia and ineffective hematopoiesis [4,24,25]. The occurrence and progression of MDS is a multistep pathological process. The precise progressive mechanism of MDS has not been fully elucidated. Microarray technology allows us to detect genetic alterations in MDS and has proven to be an effective method to determine novel indicators of other diseases [26,27]. In the current study, bioinformatics analyses were performed to analyze the role of genetic alterations in the development of MDS [28]. In the present study, we identified, RAG1 as latent core genes for myelodysplastic syndromes by bioinformatics, and determined the expression of RAG1 in the bone marrow of MDS patients. We discovered that RAG1 was markedly down-regulated in MDS patients compared with normal control, after treatment, RAG1 was markedly up-regulated. ROC analysis showed that RAG1 lowexpression could distinguish MDS cases from normal controls. Thus, we believe that RAG1 may be a potential biomarker for the accurate diagnosis of MDS. Our preliminary work revealed that compared to patients with high RAG1 expression, patients with low RAG1 expression had a poorer OS.
In addition, we also observed an interesting phenomenon that compared with the newly diagnosed MDS, the expression level of RAG1 in sAML was relatively higher, and a similar phenomenon was found in the early studies of our group, which made us curious. The disadvantage is that our sample size is too small, and we did not measure the expression of RAG1 in primary AML (pAML) patients. In the future, we will expand the sample size and measure the expression of RAG1 in primary AML patients to further explore the expression and clinical significance of RAG1 in MDS, sAML, and pAML.
In conclusion, the results showed that the expression of RAG1 was down-regulated in MDS patients. Lower RAG1 expression in BM samples of MDS patients was associated with adverse OS. RAG1 may be a potential prognostic biomarker for MDS. RAG1 may be a target for improving treatment strategies for MDS patients.