Comprehensive analysis of m6A methylation modification in chronic spinal cord injury in mice

Chronic spinal cord injury (CSCI) is a catastrophic disease of the central nervous system (CNS), resulting in partial or complete loss of neurological function. N6‐methyladenosine (m6A) is the most common form of reversible posttranslational modification at the RNA level. However, the role of m6A modification in CSCI remains unknown. In this study, we established a CSCI model using a water‐absorbable polyurethane polymer, with behavioral assessment, electrophysiological analysis, and histochemical staining for validation. Methylated RNA immunoprecipitation sequencing (meRIP‐seq) and messenger RNA sequencing (mRNA‐seq) were jointly explored to compare the differences between CSCI spinal tissue and normal spinal tissue. Furthermore, real‐time quantitative reverse transcription pcr (qRT–PCR), western blot analysis, and immunofluorescence staining were used to analyze m6A modification‐related proteins. We found that water‐absorbable polyurethane polymer simulated well chronic spinal cord compression. Basso mouse scale scores and electrophysiological analysis showed continuous neurological function decline after chronic compression of the spinal cord. meRIP‐seq identified 642 differentially modified m6A genes, among which 263 genes were downregulated and 379 genes were upregulated. mRNA‐seq showed that 1544 genes were upregulated and 290 genes were downregulated after CSCI. Gene Ontology terms and enriched Kyoto Encyclopedia of Genes and Genomes pathways were also identified. qRT–PCR, western blotting, and immunofluorescence staining showed that Mettl14, Ythdf1, and Ythdf3 were significantly upregulated after CSCI. Our study revealed a comprehensive profile of m6A modifications in CSCI which may act as a valuable key for future research on CSCI.


| INTRODUCTION
Chronic spinal cord injury (CSCI) is a catastrophic disease of the central nervous system (CNS) that usually causes partial or complete loss of neurological function; to date, there is no effective treatment. 1,2 The major cause for CSCI is chronic spinal cord compression due to intervertebral disc herniation, ligament sclerosis, and bone hyperplasia. 3,4 Although surgical treatment, such as decompression, can effectively relieve the symptoms of spinal cord compression, neuronal loss, apoptosis, and inflammatory factor infiltration persist due to the occurrence of secondary SCI. 5,6 Therefore, exploring the pathological changes and internal mechanisms of CSCI is critical for CSCI treatment.
Epigenetic modification is a widespread form of posttranscriptional modification in organisms and includes DNA modifications, histone modifications, chromatin remodeling, and RNA modifications. 7 Epigenetic regulation affects gene expression but does not alter the DNA sequence and plays a vital role in CNS diseases. 8 DNA methylation affects local chromatin structure to change gene expression and has been reported to regulate axon regeneration in CNS injuries. 9 Histone deacetylase 1 (HDAC1) is a crucial epigenetic factor that removes acetyl from histones.
Inhibiting HDAC1 was found to promote neuronal loss and DNA damage, which further deteriorated behavioral outcomes after stroke. 10 It has been reported that knockdown of the epigenetic factor Ubiquitously Transcribed tetratricopeptide repeat on chromosome X (UTX) can facilitate axon regeneration and angiogenesis after acute spinal cord injury by targeting microRNA-24. 11,12 Epigenetic modifications are also involved in angiogenesis after chronic compressive spinal cord injury through the Xist/miR-32-5p/Notch-1 axis. 13 Among epigenetic modifications, N6-methyladenosine (m6A) is the most common and abundant mRNA modification. 14,15 m6A modification is a reversible and dynamic biological process that is mediated by methyltransferases (writers, mainly METTL3, METTL14 and WTAP), demethylases (erasers, mainly FTO and ALKBH5), and RNA-binding proteins (readers, mainly YTHDF1 and YTHDF3). 16 m6A modifications have been reported to play a crucial role in diverse biological processes.
The development of various tumors is accompanied by m6A methylation modifications, including in colorectal cancer, 17 bladder cancer, 18 and prostate cancer. 19 Moreover, m6A modifications has been found in many CNS diseases. 15 m6A was found to promote the synthesis of a range of proteins in a YTHDF1-dependent manner. Knockout of YTHDF1 caused hippocampal deficits, resulting in damage to learning and memory. 20 Mettl14 deficiency in the substantia nigra suppressed tyrosine hydroxylase production in dopaminergic neurons, ultimately causing impairment of motor function. 21 Another study confirmed that microRNA-421-3p can directly target YTHDF1 to downregulate the expression of p65, a key component of the NF-κB signaling pathway, thus preventing inflammation in cerebral ischemia/reperfusion injury. 22 Recently, studies on m6A modifications have emerged focusing on acute SCI. m6A modifications have been reported in acute spinal cord injury, especially in regulating neuronal apoptosis. 23,24 Xing et al. 25 analyzed altered m6A modifications following acute SCI in zebrafish larvae. However, a comprehensive analysis of m6A modification changes in mice with chronic compressive spinal cord injury has never been reported.
This study aimed to explore the modification of m6A in chronically injured and normal mice tissues. We constructed a chronic compressive spinal cord injury model using a water-absorbable polyurethane polymer.
Methylated RNA immunoprecipitation sequencing (meRIP-seq) and mRNA-seq were used to explore the m6A modifications in injured and normal spinal tissue. We jointly analyzed the differentially expressed m6A-modified genes and identify the enriched GO terms and the KEGG pathways. Collectively, our research may provide new insights for the treatment of CSCI.

| Animals and ethics statement
Eight-week-old C57BL/6 female mice (weighing approximately 20-23 g) were purchased from Hunan SJA Laboratory Animal Co, Ltd and housed in the Department of Laboratory Animals, Central South University. All animals had free access to food and water and lived in a temperaturesuitable room under natural light. All procedures carried out were approved by the Animal Ethics Committee of Central South University.

| Establishment of chronic spinal cord injury
Mice were randomly divided into two groups: a control group and a CSCI group. Based on previous research, we applied a waterabsorbable polyurethane polymer (Guangzhou Fischer Chemical Co) to construct the CSCI model. 26 In CSCI group, the mice were anesthetized with 75 mg/kg sodium pentobarbital. The dorsal skin at the T10 level was incised, followed by fascia and muscle stripping to expose the T10 lamina. A portion of the lamina was removed, and a polyurethane polymer sheet (1.5 ×1.0 × 0.5 mm) was inserted into the epidural space.
While in the control group, after the exposure of the T10 lamina, a portion of the lamina was removed without placing a polyurethane polymer sheet in the epidural space. Then, the muscles, fascia and skin were sutured in layers. After the mice were awake, they were returned to their rearing cages. For the next 3 days, the mice were injected with penicillin daily to prevent infection.

| Histological and immunofluorescence staining
The mice were anesthetized and sequentially perfused from the left ventricle with heparinized saline and 4% paraformaldehyde. Subsequently, a 10 mm length of spinal cord in the injured area was harvested and dehydrated for sectioning. For histological staining, samples were dehydrated in gradient alcohol and embedded in paraffin. Spinal cord segments were transected into 8 µm thick slices and stained with HE staining kit (Solarbio) according to the manufacturer's instructions.
For immunofluorescence staining, samples were dehydrated with sucrose solution and embedded with optimal cutting temperature LI ET AL.   Table S1) were added to the cDNA mixture for qRT-PCR. CT values were detected using a quantitative PCR system (ABI). GAPDH was used as an internal reference, and the relative gene expression was calculated using the ∆∆ 2 C − t method.

| Western blotting
Total protein was extracted from the spinal cord using RIPA lysis buffer (Beyotime), and the protein concentration was detected with a BCA kit (Beyotime). The denatured proteins were separated via electrophoresis using 10% sodium dodecyl sulfate polyacrylamide gel electrophoresis gels and then electrotransferred to PVDF membranes
Neuroelectrophysiological analysis was performed at 28 days postinjury as described in our previous research. 28 Stimulating electrodes were placed on the surface of the skull (caudal/lateral to the bregma:

| Bioinformatics analysis
Raw reads acquired from the Illumina HiSeq sequencer were filtered using Trimmomatic (version 0.36) to remove low-quality and adaptercontaminated beads. Clean beads were harvested through unique identifier (UID) deduplication. Briefly, we performed cluster analysis on all sequences of the same UID to obtain different clusters. Then, sequence alignments were performed on the multiple sequences under each cluster, thereby obtaining a consensus sequence. The deduplicated clean data were mapped to the reference genome of mice using STAR software (version 2.5.3a). The specific region to which the protein binds (peak calling) was analyzed using exomePeak (version 3.8) software. After m6A peak annotation with bedtools (Version 2.25.0), the differentially methylated sites were identified with a Python script using Fisher's test. GO and KEGG pathway enrichment analyses were carried out using KOBAS software (version: 2.1.1). The threshold of significant enrichment was set as p < 0.05.

| Statistics
GraphPad Prism 8.0 was used to analyze the data. Data are presented as the mean ± standard deviation (SD). Unpaired Student's t test was applied to compare differences between two groups, such as in qRT-PCR, western blotting, and fluorescence quantitative analyses.
Repeated-measures two-way analysis of variance was used to analyze BMS scores. p < 0.05 was considered to indicate a significant difference (*p < 0.05, **p < 0.01).

| Validation of the CSCI model in mice
We implanted a water-absorbable polyurethane polymer to induce chronic compressive spinal cord injury according to previous studies. 13,26 Since our experimental subjects were mice, we adjusted the size of the polyurethane polymer to 1.5 × 1.0 × 0.5 mm before implantation. Four weeks after implantation, the size of the waterabsorbable material reached approximately 2.5 × 1.65 ×1.0 mm ( Figure 1A). The BMS scores showed that the locomotor function in the CSCI group was greatly declined compared to that in the control group. In the CSCI group, BMS scores continued to decline in the first 3 weeks, and became stabilized from the fourth week to the  Figure 1B). Electrophysiological experiments showed that mice in the CSCI group displayed lower motor evoked potentials, indicating severely impaired neurological conduction capacity ( Figure 1C,D). To further explore the effect of chronic compression on the spinal cord, we performed HE staining at 3, 7, 14, 21, 28 days, and 2 months post injury. Similar to the trend of BMS scores, in the CSCI group, the spinal morphology was destroyed, and neurons were lost in large numbers. Spinal cord tissue defects gradually increased during the first 3 weeks, and tended to stabilize from the fourth week to 2 months ( Figure 1E).

| Quality monitoring for meRIP-sequencing results
To explore the changes in m6A methylation modifications in spinal cord tissues after CSCI, we performed meRIP sequencing.
The raw data showed that the average data size of the "IP" samples was 12.96 Gb and that of the "Input" samples was 10.77 Gb (Table 1). After removing low-quality and adaptercontaining reads through UID deduplication, the clean reads of "IP" samples reached approximately 76.82% of raw reads, and the "Input" samples reached approximately 69.81% (Table 1). The Clean_Q30 (%) value of each sample, which represents the quality of sequencing, exceeded 98.50%, and the GC base content of clean reads exceeded 52% (Table 1). Clean data were further aligned to the reference genome of mice to acquire comprehensive transcript information. The proportion of sequences for which each sample had a unique alignment position on the reference sequence reached 94% ( Table 2). The numbers of reads mapped to the reference genome at both ends of the sequence were basically the same ( Table 2).

| General features of m6A methylation in normal and CSCI spinal tissues
We analyzed the read distribution in different regions of the reference genome and found that most of the reads were distributed in the coding sequences (CDSs), 3'UTR and 5'UTR regions in both groups ( Figure 2). By searching for protein binding regions, we found 22798 m6A peaks in the control group and 22866 m6A peaks in the CSCI group ( Figure 3A).
Compared to the control group, 263 m6A peaks were downregulated and 379 m6A peaks were upregulated in the CSCI group ( Figure 3B). The top 10 m6A hypermethylated or hypomethylated genes are listed in Table 3.  Figure 3G). Among them, 1353 genes were only detected in the control group, and 1183 genes were unique to the CSCI group ( Figure 3H).

| Functional analysis of differential m6A methylation and differentially expressed genes after CSCI
We performed functional analysis to explore the effect of m6A modification on CSCI via GO analysis and KEGG pathway analysis.
T A B L E 1 Detailed information of raw data obtained by sequencing in the CSCI group and control group

| Conjoint analysis of differentially expressed and methylated genes after CSCI
The mRNA-seq results showed that there were 290 downregulated genes and 1544 upregulated genes ( Figure 6B,C)  had upregulated mRNA expression, and 2 genes had downregulated mRNA expression ( Figure 6D).

| Relative expression of m6A-related proteins after CSCI
M6A is a reversible epigenetic modification that requires the participation of multiple enzymes and proteins containing writers (catalyzing adenylate methylation), erasers (demethylating adenylate) and resders (identifying mRNA domains). 29 Therefore, we performed qRT-PCR analysis to detect these genes including methylases (Mettl3, Mettl14, WTAP), demethylases (FTO and ALKBH5), and RNA-binding proteins (YTHDF1, YTHDF3). The results showed that mettl14, YTHDF1 and YTHDF3 were significantly upregulated after CSCI, while other proteins showed no significant differences ( Figure 7A). Western blot analysis results also verified that mettl14, YTHDF1 and YTHDF3 expression was significantly increased in the CSCI group ( Figure 7C). Immunofluorescence costaining indicated that mettl14 was mainly expressed in neurons and astrocytes, with upregulated expression after CSCI ( Figure 7B, Supporting Information: Figure S1A). After CSCI, expression of RNA-binding proteins YTHDF1 and YTHDF3 was increased, while methylases Mettl3 showed no obvious alteration (Supporting Information: Figure S2).

| DISCUSSION
CSCI is a catastrophic central nervous system disease which up to date there is no effective treatment. 1 The pathological process of CSCI is complex, often resulting in vasculature damage and neuronal loss. 30 Unlike ASCI, which causes a rapid decline or loss of neurological function and could be divided into acute phase, subacute phase, and chronic phase according to different pathological changes, CSCI tends to progress more slowly. Common modeling approaches for CSCI include using stainless steel screws, 31 expanding polymers, 30 silastic sheets 32 and water-absorbing materials 13 or using Tiptoewalking Yoshimura (twy/twy) mice, which spontaneously deposit calcium at the C1-C2 level. 33 In this study, we used a waterabsorbable polyurethane polymer to construct the CSCI model. We found that after spinal cord compression, the BMS scores decreased significantly in the first 3 weeks, with tissue defects and neuronal loss, and tended to be stabilized 4 weeks after injury. m6A methylation is the most common posttranscriptional modification at the RNA level. 34 An increasing number of reports have shown that m6A methylation modifications involved in multiple diseases of the central nervous system. 35 In a cerebral ischemia/ reperfusion (I/R) injury model, microRNA421-3p specifically targeted YTHDF1 to reduce its expression. YTHDF1 further recognized and bound to the p65 m6A site to promote p65 translation, thereby causing an inflammatory response. 22 Another study comprehensively analyzed the m6A methylation profile after acute SCI in zebrafish and found that a bulk of neural regeneration-related genes were hypomethylated but highly expressed. 25 However, the role of m6A modifications after CSCI has rarely been reported. In this study, we first performed meRIP-seq and mRNA-seq by isolating spinal cord tissues after chronic compression. We observed 22798 m6A peaks in the sham group and 22866 m6A peaks in the CSCI group. A total of 642 m6A differentially methylated peaks were identified, of which 263 peaks showed m6A downregulation and 379 showed m6A upregulation, indicating that m6A methylation changes did occur in spinal cord tissue after CSCI. We also found that the m6A peak distribution was similar in the two groups, with most peaks located in the CDS, 3'UTR and introns. Meanwhile, using mRNA-seq, we found that 290 genes were downregulated and 1544 genes were upregulated after CSCI. These results showed that m6A methylation modifications were likely to regulate the expression of a series of genes after CSCI to affect the pathological process of CSCI.
By association analysis of meRIP-seq and mRNA-seq, we found that 21 genes were hypomethylated, of which 19 genes were upregulated and 2 genes were downregulated; 34 genes were hypermethylated, of which 25 genes were upregulated and 9 genes were downregulated. Our data indicated that m6A regulation was complex. Both hypomethylated and hypermethylated m6A could induce or reduce mRNA expression, which was consistent with previous studies. 36,37 The differentially expressed genes with altered m6A methylation included small proline-rich repeat protein 1A (sprr1a), c-type lectin domain family 4 member D (clec4d) and apelin receptor (aplnr). Sprr1a is highly expressed in injured neurons and thought to be a regeneration-associated gene (REG) that promotes spontaneous axon growth. Macrophage-inducible clec4d was found to promote neuroinflammation and neuronal damage after traumatic brain injury (TBI). 40,41 Apelin and its receptor aplnr have been reported to alleviate neuropathic pain after compression spinal cord injury. 42 Our study confirmed that these differentially expressed genes after CSCI were regulated by m6A modification, which may shed new light on CSCI pathophysiology for further studies.
Through GO analysis and KEGG pathway analysis, we found that m6A differentially methylated genes were enriched in developmental process, ion binding, the wnt signaling pathway, neuroactive ligand-receptor interaction, the PI3K-AKT pathway and the NFkappa B pathway. Activation of the wnt signaling pathway is critical for neurological recovery after SCI. It has been reported that simulating the wnt signaling pathway can promote vascular regeneration 43 and axon regrowth 44 and reduce inflammatory responses. 45 The PI3K-AKT pathway participates in many fundamental cellular processes, including cell growth, cell proliferation, and cell motility.
Curcumin-and lipopolysaccharide-activated olfactory ensheathing cells were found to provide proangiogenic effects through the PI3K-AKT signaling pathway. 46 Another study verified that activating the PI3K-AKT signaling pathway facilitated nerve reconstruction and attenuated neuropathic pain after peripheral nerve injury. 47 The NFkappa B pathway is a classical pathway involved in inflammatory and immune responses. miR-182 reduced the expression of the upstream target of the NF-kappa B pathway IκB kinase β to inhibit apoptosis and the inflammatory response. 48 Although these classical signaling pathways have been widely reported in SCI, little is known about whether m6A methylation regulates these pathways. Our study indicated that m6A modifications contributed to activation or blockade of these signaling pathways after CSCI; however, the specific role and mechanisms still need further exploration in the future.   The process of m6A modification requires the interaction of methyltransferases (writers), demethylases (erasers), and RNAbinding proteins (readers). 34 Mettl3 and mettl14 are common methyltransferases in the central nervous system. Xing et al. found that mettl3 was increased after acute spinal injury in zebrafish and was mainly expressed in astrocytes and neural stem cells. 25 Other studies have shown that mettl3 is significantly increased in a mouse model of chronic inflammatory pain and is mainly located in neurons.
Downregulation of mettl3 was reported to alleviate inflammatory pain by regulating pri-miR-65-3p processing or 10-11 translocation methylcytosine dioxygenase 1 expression. 49,50 Mettl14 has been reported to increase neuronal apoptosis by inhibiting EEF1A2 expression 24 or promoting miR-375 maturation after acute spinal cord injury. 23 Another study demonstrated that mettl14 is essential for injury-or pten deletion-induced axon regeneration of DRG neurons. 51 However, m6A modifications in chronic compressive spinal cord injury have never been reported. Our study showed that mettl14, Ythdf1, and Ythdf3 were significantly overexpressed after CSCI and mainly expressed in neurons, which was further validated by qRT-PCR, western blot analysis, and immunofluorescence.
Our study is not without limitations. First, we obtained whole spinal tissue for meRIP-seq and mRNA-seq. Although the total m6A modifications were assessed, it was difficult to distinguish m6A methylation changes in specific cell types. In further in-depth