Role of N6-methyladenosine-related lncRnas in pseudoexfoliation glaucoma

ABSTRACT

To explore the role of lncRNA m⁶A methylation modification in aqueous humour (AH) of patients with pseudoexfoliation glaucoma (PXG). Patients with open-angle PXG under surgery from June 2021 to December 2021 were selected. Age- and gender-matched patients with age-related cataract (ARC) were chosen as control. Patients underwent detailed ophthalmic examinations. 0.05–0.1 ml AH were extracted during surgery for MeRIP-Seq and RNA-Seq. Joint analysis was used to screen lncRNAs with differential m⁶A methylation modification and expression. Online software tools were used to draw lncRNA-miRNA-mRNA network (ceRNA). Expression of lncRNAs and mRNAs was confirmed using quantitative real-time PCR. A total of 4151 lncRNAs and 4386 associated m⁶A methylation modified peaks were identified in the PXG group. Similarly, 2490 lncRNAs and 2595 associated m⁶A methylation modified peaks were detected in the control. Compared to the ARC group, the PXG group had 234 hypermethylated and 402 hypomethylated m⁶A peaks, with statistically significant differences (| Fold Change (FC) |≥2, p < 0.05). Bioinformatic analysis revealed that these differentially methylated lncRNA enriched in extracellular matrix formation, tight adhesion, TGF- β signalling pathway, AMPK signalling pathway, and MAPK signalling pathway. Joint analysis identified 10 lncRNAs with differential m⁶A methylation and expression simultaneously. Among them, the expression of ENST000000485383 and ROCK1 were confirmed downregulated in the PXG group by RT-qPCR. m⁶A methylation modification may affect the expression of lncRNA and participate in the pathogenesis of PXG through the ceRNA network. ENST000000485383-hsa miR592-ROCK1 May be a potential target pathway for further investigation in PXG m⁶A methylation.

KEYWORDS:

• Pseudoexfoliation glaucoma
• methylation
• aqueous humour
• lncRNA

Introduction

Pseudoexfoliation glaucoma (PXG) is the most common type of secondary glaucoma, with abnormal extracellular fibrillar material accumulated in the anterior chamber as its typical clinical presentation [1]. PXG has strong genetic susceptibility and shows specific regional and ethnic aggregation.The variation of LOXL1 gene is found to be mainly responsible for the fibrillar material and disfunctioning elastin to increase the outflow resistance of the anterior trabecular meshwork [2]. Other genes, such as the CACNA1A, TBC1D21 and ATXN2 are also found in patients with PXG according to the previous reports [3].

In recent years, epigenetics has drawn increasing attention. Accumulating evidence suggests that epigenetic factors may be involved in the pathogenesis of PXG. Two sounding studies at home and abroad confirmed that LOXL1 promoter methylation was significantly increased in patients with PXG compared to control patients, resulting in a decrease in the expression level of LOXL1 mRNA and protein product [4,5]. Noncoding RNAs are involved in gene expression regulation as another epigenetics. A meta-analysis showed that CDKN2B-AS1 lncRNA contributed to the development of PXG via regulating the expression of CDKN2B, a member of the TGF-β signalling pathway [6]. LOXL1-AS1 lncRNA localized to the nucleus and was altered in response to oxidative stress and cyclic mechanical stress. Dysregulation of its expression was critical for the global gene expression in ocular cells of PXG [7]. Five miRNAs (miR-122-5p, miR-3144-3p, miR-320a, miR-320e, and miR-630) from patients with PXG were identified and involved in potential glaucoma pathways, including focal adhesion, tight junctions, and TGF-β signalling. Among them, miR-122-5p may target three glaucoma-associated genes: OPTN, TMCO1, and TGF-β [8]. The joint participation of multiple risk genes is an important reason for the occurrence and development of PXG, but a functional role for these variants has not been established.

RNA methylation, a kind of post-transcriptional mRNA modifications, is further classified as N6-methyladenosine (m⁶A), 5-methylcytosine, N1-Methyladenosine, with the former as the most prevalent [9]. m⁶A plays an important role in regulating RNA processing, localization, transportation, splicing, and translation. Compared with oncology, m⁶ modification research in Ophthalmology is still in its early stage. Recently, emerging evidence has indicated that the abbrerant modification of m⁶A is associatited with progression of pterygium [10], Graves’ ophthalmopathy [11], high myopia [12], and experimental fungal keratitis [13]. Decreased m⁶A levels of HINT2 mRNA, a tumour suppressor, indicated a poor prognosis and progression in ocular melanoma [14]. Although the role of m⁶A modification has been investigated in many diseases, data are lacking in PXG.

In our previous work, we showed significantly higher m⁶A levels in the PXG aqueous humour and abnormal methylation in target mRNAs [15]. Here, we hope to conduct in-depth research on the impact of m⁶A methylation on ceRNA networks and explore the mechanism of m [6]A methylation modification in PXG.

Methods

Participants

The Ethics Committee of the First People’s Hospital of Kashi Prefecture approved the study (2021ksyd-22). Informed consents were acquired by the participants. Besides, this study adhered to the tenets of the Declaration of Helsinki. PXG was diagnosed by elevated intraocular pressure (IOP >21 mmHg), visible extracellular fibrillar material in the anterior segment structures, glaucomatous optic nerve head changes, and visual field defects.

In total, 30 patients with PXG were enrolled and the exclusion criteria were those with eye surgeries, or systemic diseases, such as hypertension, diabetes, coronary heart disease, and autoimmune diseases. Age- and sex-matched ARC patients were recruited as control. All participants were Uyghur. They belong to the yellow and white mixed race and constitute a minority population in China.

AH samples were acquired during surgery for both groups. They were immediately collected in a cryopreservation tube, snap-frozen in liquid nitrogen, and stored at − 80°C until RNA isolation. Table 1 summarized the demographic data of the patients.

Table 1. Patient characteristics of PXG and ARC.

	Discovery phase			Validation phase
	PXG	ARC		PXG	ARC
Characteristics	(n = 9)	(n = 9)	P-value	(n = 21)	(n = 21)	P-value
Sex (male/female)	7/2	7/2	1.000	15/6	9/12	0.118
Age, years, mean±SD	76.44 ± 6.33	73.56 ± 5.20	0.306	72.67 ± 7.82	71.81 ± 8.95	0.743
Eye (right/left)	5/4	3/6	0.637	11/10	14/7	0.530
BCVA, LogMAR, mean±SD	0.82 ± 0.55	0.64 ± 0.23	0.397	1.13 ± 0.60	0.72 ± 0.46	0.016*
IOP, mmHg, mean±SD	30.11 ± 6.79	13.56 ± 1.74	0.000*	30.00 ± 6.24	13.06 ± 3.19	0.000*

Data are shown as number of cases or mean±SD, as appropriate. P-values were calculated using Student’s t-tests (age, BCVA, and IOP) and χ2-tests (sex and eye). Abbreviations: PXG, pseudoexfoliation glaucoma; ARC, age-related cataract; n, number of patients; SD, standard deviation; BCVA, best-corrected visual acuity; IOP, intraocular pressure. * Statistically significant difference (p < 0.05).

Workflow and RNA isolation

The discovery phase consisted 9 pairs of samples for MeRIP and RNA-seq analysis. The rest 21 pairs of samples were independently validated using real-time PCR (Figure 1). We extracted RNA using TRI-Reagent (Cat. no.15596026, Thermo Fisher Scientific, Waltham, MA, USA) and examined the quality using a NanoDrop ND-1000 spectrophotometer (Thermo Fisher Scientific). Only the samples with OD260/OD280 ratios within 1.8–2.1 were selected for analysis.

Figure 1. Workflow of the study protocol.

MeRIP sequencing

The m⁶A-IP-Seq was analysed by CloudSeq Biotech Inc. (Shanghai, China). Samples in the discovery phase were immunoprecipitated by GenSeq® m⁶A-IP Kit (GenSeq Inc.). Specifically, RNA samples (n = 3 pools per group) were subjected to immunoprecipitation using the GenSeq® m6A-IP Kit (GenSeq Inc., Shanghai, China) according to the manufacturer’s instructions. Briefly, the RNA was randomly fragmented to ~ 200 nucleotides using RNA fragmentation reagents. Protein A/G beads were coupled to the m6A antibody via rotation at 25°C for 1 h. Next, the RNA fragments were incubated with the bead-linked antibodies at 4°C for 4 h with rotation. This was followed by the washing of the RNA/antibody complexes several times, after which the captured RNA was eluted and purified. Using the GenSeq® Low Input Whole RNA Library Prep Kit (GenSeq Inc.), RNA libraries were constructed through immunoprecipitated and input samples. Thereafter, the libraries were qualified using an Agilent 2100 Bioanalyzer and sequenced on a NovaSeq platform (Illumina, San Diego, CA, USA). Paired-end reads were then obtained using an Illumina NovaSeq 6000 sequencer and checked for errors in base calling using a quality score of Q30. Cutadapt software (v1.9.3) was used for 3´-adaptor trimming and low-quality read removal. Thereafter, clean reads of the input libraries were aligned to the reference genome (hg19) using STAR software. Next, clean reads from all the libraries were aligned to the reference genome using HISAT2 software (v2.0.4). The methylated sites on the RNAs (peaks) were identified using MACS software; differentially methylated sites were identified using diffReps. The peaks identified by both programmes that overlapped with mRNA exons were chosen by homemade scripts.

High-throughput sequencing

rRNA was removed using the GenSeq® rRNA Removal kit (GenSeq Inc.) according to the manufacturer’s instruction. The rRNA-depleted samples (n = 3 pools per group) were then subjected to library construction using the GenSeq® Low Input RNA Library Prep kit (GenSeq Inc.) according to the manufacturer’s instruction. The libraries thus obtained were subjected to quality checks and quantification using the Bioanalyzer 2100 system (Agilent Technologies, Inc., Santa Clara, CA, USA). Library sequencing was performed on an Illumina NovaSeq 6000 sequencer using 150-bp paired-end reads; errors in base calling were checked using a quality score of Q30. After 3´-adaptor trimming and low-quality read removal using Cutadapt (v1.9.3), the high-quality trimmed reads were aligned to the reference genome (Hg19) using HISAT2 software (v2.0.4). Next, HTSeq software (v0.9.1) was used to obtain mRNA raw counts for expression profiling, while edgeR (v3.16.5) was used for normalization, and lncRNAs with p < 0.05 and |fold change| ≥ 2 were considered differentially expressed.

Conjoint analysis and lncRNA-mediated ceRNA interactions identification

Conjoint analysis was used to explore the relationship between the RNA-Seq and MeRIP-Seq data. Ten lncRNAs with differentially m⁶A-modified expression were included for prediction. The mircode online tool (http://www.mircode.org/) was used to predict the lncRNA and miRNA interaction. The top five miRNAs associated with each lncRNA were obtained. The putative miRNA-mRNA interactions were collected from miRDB (http://mirdb.org/), TargetScan (http://www.targetscan.org/vert_72/), and miRTarBase (http://mirtarbase.mbc.nctu.edu.tw/php/index.php) databases. The top 100 mRNAs of each miRNA, predicted by all three databases, were included for the intersection with the differentially expressed mRNAs through high-throughput sequencing.

The lncRNA-mediated ceRNA network was constructed on the basis that ceRNA can bind to miRNA through MREs. This ceRNA network was visualized using Cytoscape software (Version 3.6.0). In this network, nodes and edges represented biological data in a direct manner in which each node represented a biological molecule, and the edges represented interactions between nodes. LncRNAs, mRNAs, and miRNAs in the ceRNA network were presented as blue diamonds, green ellipses, and red triangles, respectively.

Bioinformatics analyses

We conducted the gene ontology (GO) and pathway enrichment analysis according to differentially methylated or expressed lncRNAs. GO and pathway enrichment analysis were performed using the Database for Annotation, Visualization, and Integrated Discovery. The ontology covers three parts: cellular component (CC), molecular function (MF), and biological process (BP). The p-value denotes the significance of GO term enrichment of the genes. Pathway enrichment analysis is a functional analysis that maps genes to the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways. The Fisher p-value denotes the significance of the pathway correlated to the conditions.

Validation of differentially expressed lncRnas and mRnas

The expression levels of 6 selected lncRNAs and 13 mRNAs were validated by RT-qPCR in validation phase. The total RNA was extracted and processed to performed a RT-qPCR process as described previously [15]. All experiments were conducted in triplicate and were repeated thrice. The primer sequences are listed in Table 2.

Table 2. Primer sequences for RT-qPCR.

Gene	Forward Primer (5´→3´)	Reverse Primer (5´→3´)
ENST00000429216	AGAGGTCGAAGTGGTTCTGG	TCCCCACTGCCATCATATCC
ENST00000569046	TCGGCCTTCTTCCTCTTCTG	TCTCTCCCAGTTTCTTCGCA
ENST00000478357	TGAGACCTGGTTCTGCTGTT	TCCTCAAGTCCCTCAAGCAG
ENST00000462917	GGGGTCCGCTTTTAACTTGG	CTTGGCTGAAACATGACGGT
ENST00000422040	TAGCTGCCCAATCTCTTCCC	CCTCTGCTGGTTGGTTTCAG
ENST00000485383	CTTCCTCTTCCTCCTCCACC	CACTGAAACCAACCATGCCA
MMP14	CAAGATTGATGCTGCTCTCTTC	ACTTTGATGTTCTTGGGGTACT
TRPM7	TGCTGTGGAAGAATGTATGACT	CAAATTTTCACCAACACCTTGC
COL5A2	ACGATCAAGCTAAGAACCTCAA	CACATTTCCATTCCGCTTAGAG
P2RX7	GAGCATGAATTATGGCACCATT	GTCTAAAAAGCTGCAAAGGGAA
HDAC6	GTGTCACTTCGAAGCGAAATAT	CCACGATTAGGTCTTCTTCCAT
TGFBR1	CTCTTCAAAAACTGGGTCTGTG	CATCAACATGAGTGAGATGCAG
VEGFB	TAGATGTGTATACTCGCGCTAC	TTTTAGGTCTGCATTCACACTG
UCP2	GTGGTCAAGACGAGATACATGA	CTGCTCATAGGTGACGAACATC
GJA1	CTTTCGTTGTAACACTCAGCAA	CGCATCACATAGAACACATGAG
ROCK1	GCTGAGCAATATTTCTCGACAC	TCCAACTGAGTAGCAAGAGTTT
SPHK2	GGGGAGGTGTGTCTTTCAGT	GTTGGGTGATGATAACGGGC
LOXL2	ACAGAATGTGAAGGAGACATCC	TGATGTTGTTGGAGTAATCGGA
SIRT1	TATACCCAGAACATAGACACGC	CTCTGGTTTCATGATAGCAAGC
ACTB	GTGGATCAGCAAGCAGGAGT	AAAGCCATGCCAATCTCATC

Statistical analysis

Statistical Package for Social Sciences version 23.0 software (IBM Corp., Armonk, NY, USA) was used in this study. Data are presented as means ± standard deviation (SD). Differences between two groups were analysed by two-tailed Student’s t-test. Statistical significance was set at p < 0.05.

Results

Alterations of m⁶A-related lncRnas in the AH of PXG patients

A total of 4386 m⁶A modified peaks were identified in the PXG group, which were related to 4151 transcripts of lncRNAs. Similarly, 2595 m⁶A modified peaks among 2490 lncRNAs transcripts were detected in the ARC group. Compared to the ARC group, 234 hypermethylated and 402 hypomethylated m⁶A peaks were identified in the aqueous humour of PXG patients, which presented a significant alteration, respectively (|FC|≥2 and p < 0.05). Moreover, hierarchical cluster analysis revealed interrelationships between the aqueous humour samples, which were also classified based on the similarities of the m ⁶A methylation levels (Figure 2).

Figure 2. The alterations of m⁶A-related lncRnas in the aqueous humour of PXG and ARC. (a) Venn diagrams show 289 overlaps of m ⁶A peaks and 414 overlaps of m⁶A- modified transcripts. (b) Hierarchical clustering presents the different lncRNA m⁶A methylation levels of PXG and ARC. Abbreviations: PXG, pseudoexfoliation glaucoma; ARC, age-related cataract.

Important biological functions and pathways of altered m⁶A-modified lncRnas

To shed light on the pathological significance of m⁶A modification in PXG, GO and KEGG pathway analyses were performed via systematically screening these different and common peaks and the related genes.

With GO analysis, the hypermethylated lncRNAs were mainly involved in a variety of biological processes, including extracellular matrix organization (ontology: biological process), extracellular matrix (ontology: cellular component), and extracellular matrix structural constituent (ontology: molecular function). The hypomethylated lncRNAs were significantly associated with positive regulation of cytoplasmic translation (ontology: biological process), cytoplasm (ontology: cellular component), and RNA binding (ontology: molecular function). GO enrichment analysis revealed that the extracellular matrix was associated with the pathological processes of PXG. Outstandingly, KEGG pathway analysis showed the upregulated peaks in PXG were significantly associated with the TGF-β signalling pathway, and the downregulated peaks were significantly associated with tight adhesion, AMPK signalling pathway, and MAPK signalling pathway (Figure 3).

Figure 3. Biological function and pathway predictions of statistically significantly altered m⁶A-modified lncRnas in PXG. (a) GO terms were significantly enriched for the hypermethylated lncRnas. (b) GO terms were significantly enriched for the hypomethylated lncRnas. (c) The top ten significantly enriched KEGG pathways for the hypermethylated lncRnas. (d) The top ten significantly enriched KEGG pathways for the hypomethylated lncRnas. Abbreviations: PXG, pseudoexfoliation glaucoma; GO: gene ontology; KEGG: kyoto encyclopedia of genes and genomes.

lncRNA expression profiles associated with m⁶A methylation

Compared with ARC, 502 lncRNAs were differentially expressed in PXG (|FC|≥2 and p < 0.05, Figure 4a), including 422 upregulated genes and 80 downregulated genes. Next, according to the conjoint analysis of RNA-seq and MeRIP-seq data (|FC|≥2 and p < 0.05), 10 altered lncRNAs with three modes of interaction were found: 17 hypermethylated and downregulated patterns; 1 hypermethylated and upregulated pattern; 16 hypomethylated and upregulated patterns. However, m⁶A hypomethylated and downregulated patterns were not identified (Figure 4c). All differently expressed m⁶A methylation-modified lncRNAs are listed in Table 3.

Figure 4. Conjoint analysis of lncRNA m⁶A modification and lncRNA expression in PXG (|FC|≥2 and p < 0.05). (a) Scatter plots showing the differentially expressed lncRNA between PXG and ARC. (b) Scatter plots showing the differentially m⁶A-methylated lncRNA between PXG and ARC. (c) Four-quadrant graph exhibiting the relationship of the m⁶A and expression level of lncRnas between the two groups. Abbreviations: PXG, pseudoexfoliation glaucoma; ARC, age-related cataract.

Table 3. 10 differently expressed m⁶A-modified lncRnas.

Transcript ID	Gene Name	Chromosome	Fold change	P-value	Methylation Regulation	Expression Regulation
ENST00000429216	HNRNPA1P35	Chr2	508.7	3.27681E–09	up	down
ENST00000569046	RP11-673C5.1	Chr15	495.9	3.03534E–09	up	down
ENST00000422040	AC008781.7	Chr5	18.9	0.015988172	up	down
ENST00000469685	RN7SL381P	Chr16	173.9	8.75635E–14	up	up
ENST00000478357	ZMIZ1	Chr10	332	4.07902E–09	down	up
ENST00000462917	EIF1	Chr17	36.26153846	1.21724E–08	down	up
ENST00000485383	HNRNPKP4	Chr3	107	1.88217E–07	down	up
ENST00000461347	NCL	Chr2	504.6	5.95846E–09	down	up
ENST00000442006	RP11-334L9.1	Chr1	113.4	3.14364E–09	down	up
ENST00000397627	PTMAP2	Chr5	8.399399399	1.93934E–09	down	up

ceRNA network construction and analysis

Ten PXG-specific lncRNAs that putatively targeted 50 PXG-specific miRNAs and the comparisons of 50 PXG-specific miRNAs and 5000 PXG-specific mRNAs were involved in the ceRNA network. As displayed in Figure 5, 10 lncRNA nodes, 50 miRNA nodes, 5000 mRNA nodes, and their edges comprised the lncRNA-miRNA-mRNA network. Among these mRNA, ROCK1, GJA1, P2RX75, and HDAC6 have been confirmed to play important roles in the pathogenesis of glaucoma by affecting extracellular matrix composition, cell ageing and death, retinal stress, and microglial cell proliferation.

Figure 5. The lncRNA-miRNA-mRNA competing endogenous RNA network in PXG. Diamonds, rectangles, and spheres represent lncRnas, miRnas, and mRnas, respectively. The black line represents the interaction between RNAs. The red nodes represent upregulated RNAs, while the blue nodes represent downregulated RNAs. Abbreviations: PXG, pseudoexfoliation glaucoma.

Validation of differentially expressed lncRNA and mRNA by RT-qPCR

To verify the transcriptome sequencing results, we used RT-qPCR to examine the expression of ENST00000485383 from the conjoint analysis results and ROCK1 from the intersectional mRNAs. The PCR results were consistent with the RNA-Seq results (Figure 6), showing significant differences (p < 0.05).

Figure 6. ENST00000485383 and ROCK1 were examined in patients with PXG (n = 21) and controls (n = 21) by quantitative PCR and normalized to ACTB. The student’s t-test was used to assess the differences in each gene between patients in PXG and controls. Abbreviations: PXG, pseudoexfoliation glaucoma; ARC, age-related cataract.

Discussion

Our study is the first to report preliminary findings on m⁶A-related pathway in pseudoexfoliation glaucoma. M⁶A modification disorder was found in the pathogenesis of PXG, and affecting lncRNA expression through the ceRNA network. LncRNA ENST00000048538 and ROCK1 May be the potential key factor behind the m⁶A modification disorder. This implies that the m [6]A-related pathway and ceRNA regulation network may be a potential therapeutic target for PXG.

First, our study found that compared with the control group, there were 502 differentially expressed lncRNAs in the aqueous humour of PXG patients, of which 422 were upregulated and 80 were downregulated. LncRNA is a DNA transcription fragment with a length exceeding 200 nucleotides. It lacks an open reading frame and only has weak protein coding [16] function [17,18]. However, it can interact with replication and transcription related elements, mRNA, ribonuclease Dicer, and other biological molecules to regulate the transcription and expression of mRNA, proteins, and other molecules, playing a key role in physiological and pathological changes in organisms [19]. LncRNA plays an important role in maintaining normal vision and multiple stages of visual impairment. For example, Six3os1, MIAT, lncRNA BB283400, lncRNA TUG1Vax2os, ENSMUST0000013486 and other lncRNAs are involved in various physiological and pathological processes of the retina [20]. Previous epigenetic studies on PXG have shown that CDKN2B-AS1 lncRNA affects TGF-β. The expression of CDKN2B, a member of the signalling pathway, promotes the development of PXG [6]. LOXL1-AS1 lncRNA is located in the nucleus and changes with oxidative stress and cyclic mechanical stress. The imbalance of its expression is crucial for the overall gene expression in PXG eye cells [21]. These differentially expressed genes may play an important role in the pathogenesis of PXG, but their specific mechanisms of action require further experimental confirmation.

Second, we detected a total of 4386 m⁶A methylation modification sites associated with lncRNA in the aqueous humour of PXG patients. Compared with the control group, there were 234 sites with differential hypermethylation and 402 sites with differential hypomethylation in the aqueous humour of PXG patients. GO analysis shows that these differentially m⁶A modified lncRNAs are mainly enriched in extracellular matrix formation, tight adhesion, AMPK signalling pathway, MAPK signalling pathway, and TGF- β signalling pathway. RNA methylation is an important post transcriptional modification event that includes modifications such as N6 methyladenosine (m⁶A), 5-methylcytosine, and N1 methyladenosine. As the most common and abundant gene modification method in eukaryotes, m⁶A methylation modification is involved in regulating RNA stability, translation, nuclear output, and decay, and participates in many biological processes without changing the base sequence [9,22]. M⁶A methylation modification not only occurs in mRNA, but also in many non coding RNA, such as lncRNA, circRNA, etc [23]. There have been a series of reports on the relationship between the methylation of lncRNA m⁶A and tumours, such as glioblastoma [24], pancreatic cancer [25], gastric cancer and liver cancer [26]. In the retinal microvascular endothelial cells of patients with diabetes retinopathy, the researchers found that the level of m⁶A methylation and expression of lncRNA SNHG7 were significantly down regulated, and the endothelial mesenchymal transformation in diabetes retinopathy was regulated through the KHSRP/MKL1 axis, providing a new target for the treatment of diabetes retinopathy [27]. In our study, the relationship between m⁶A methylation modification of lncRNA and PXG was first reported. This suggests that lncRNAs affected by m⁶A methylation modification may affect the occurrence and development of PXG in the aforementioned pathways.

However, the RNA stability mediated by m⁶A methylation modification is a complex process that may involve other regulatory factors, and changes in m⁶A methylation modification in genes may not necessarily affect gene expression [28]. In order to further understand the impact of m [6]A methylation modification on lncRNA expression in the aqueous humour of PXG patients, we conducted a joint analysis of these differentially expressed lncRNAs and differentially methylated lncRNAs of m⁶A. A total of 10 lncRNAs were identified, including ENST000000429216, ENST000000569046, ENST000000422040, ENST00000046985, ENST000000478357, ENST000000462917, ENST000000485383, ENST000000461347, ENST00000442006, and ENST000000397627. Among them, ENST000000429216, ENST000000569046, and ENST000000422040 showed upregulation of m⁶A methylation modification and downregulation of expression; ENST000000469685 undergoes upregulation of m⁶A methylation modification and downregulation of expression; ENST000000478357, ENST000000462917, ENST000000485383, ENST00000461347, ENST000000442006, and ENST000000397627 showed downregulation of m⁶A methylation modification and upregulation of expression. In this study, we did not identify lncRNAs that underwent both downregulation of m⁶A methylation modification and downregulation of expression. The RT-qPCR results further confirmed the differential expression of four lncRNAs: ENST000000429216, ENST000000569046, ENST000000422040, and ENST00000048538. The expression changes of the above 10 lncRNAs may be caused by m⁶A methylation modification, and they are the most likely candidate target genes for future research.

In addition, lncRNA can participate in the transcription process of RNA as a competitive endogenous RNA (ceRNA). The vast majority of ceRNAs contain potential miRNA response elements (MREs) that can competitively bind to miRNAs, thereby acting as miRNA sponges and inhibiting the regulatory effect of miRNAs on downstream RNA transcription. They regulate gene expression at the epigenetic, transcriptional, and post transcriptional levels, forming a regulatory network for gene transcription [29]. This network of non-coding RNA and coding RNA interacting together has been reported in various complex diseases, especially neurodegenerative diseases [19,30], immune diseases [31], and cancer [32]. In 2020, Zhou et al. [33] constructed a ceRNA network map related to primary open angle glaucoma using genetic data from patients with primary open angle glaucoma in a publicly available database, and identified hub lncRNAs such as OIP5-AS1, DNAJC27-AS1, AF121898, and SNX29P2, providing new insights for further exploration of the potential mechanisms of primary open angle glaucoma. This study found that there were 10 differentially methylated and expressed lncRNAs of m⁶A in the aqueous humour of PXG patients. Further prediction revealed that four lncRNAs, including ROCK1, GJA1, P2RX7, and HDAC6, influenced the expression of four downstream mRNA through the ceRNA network. RT-qPCR has confirmed their expression levels. Among them, ROCK1 participates in the pathogenesis of glaucoma by affecting the inflammatory response of retinal pigment epithelial cells [34] and the fibrous proliferation of trabecular meshwork cells [35]. GJA1 is a gap junction protein present mostly in glial cells of retina and optic nerve. Its channel blockers is a novel therapy for glaucomatous neuropathies [36]. It also plays a role in the influence of trabecular meshwork cells on mechanical stretching [37]. P2RX7 affects the proliferation and migration of retinal microglia [38] and participates in the death of retinal ganglion cells [39]. HDAC6 is involved in the process of excitatory toxic degeneration of retinal ganglion cells [40]. We speculate that the direction of changes in lncRNA and mRNA expression is consistent, based on the ceRNA network pattern, as lncRNA promotes mRNA expression by binding to miRNA. Among the above four mRNA, only the ENST000000485383-hsa miR592-ROCK1 pathway conforms to this pattern. From this, we speculate that in the aqueous humour of PXG patients, low methylation modification of m⁶A transcripts leads to overexpression of lncRNA ENST000000485383. Through the sponge mechanism, a large number of miRNAs are bound in the ceRNA network, relieving the inhibitory effect of miRNAs on downstream mRNA expression, thereby promoting the expression of downstream ROCK1, and participating in the pathogenesis of PXG by affecting the inflammatory response of retinal pigment epithelial cells and the fibrous proliferation of trabecular meshwork cells. In clinical practice, the efficacy of ROCK inhibitors in reducing intraocular pressure in glaucoma has been well established, pointing towards the involvement of ROCK in glaucoma pathogenesis. Our study proposes that aside from targeting ROCK, the intervention of lncRNAs (ENST000000485383) and miRNAs (miR592) may also have a potential impact on lowering intraocular pressure.

This study has some limitations. First, considering the low incidence of PXG in the Chinese population, the sample size was small and all the enrolled patients with PXG were Uyghur elderly. Thus, there may be a selection bias, implying that further studies are required to determine whether these results are generalizable. Secondly, this study utilized MeRIP-Seq and RNA-Seq to analyse lncRNAs and m6A methylation modification. While these are important tools as their high-throughput sequencing, quantitative analysis, detection of heterogeneity, no prior knowledge requirement, and data reproducibility, it’s essential to acknowledge that these methods have limitations, such as potential biases in RNA-seq and the need to validate MeRIP-Seq results. We suggest that the authors consider these limitations in future research and propose possible improvements. Thirdly, no relevant cell or animal experiments have been conducted in this study, and the evidence of ceRNA network pathogenesis is not sufficient. Further experimental studies are needed to confirm the specific mechanism of methylated lncRNA in the pathogenesis of PXG.

In summary, this study proposes for the first time that m [6]A methylation modification may be involved in the pathogenesis of PXG by affecting lncRNA expression and thus through the ceRNA network. Formation of extracellular matrix, tight adhesion, AMPK signalling pathway, MAPK signalling pathway, and TGF- β signalling pathway may be a key target for lncRNA methylation modification. ENST00000048538 and ROCK1 May be potential target genes for future research.

Disclosure statement

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.

Data availability statement

Data are available at reasonable request to the corresponding authors.

Additional information

Funding

This work was supported by the National Key R&D Project of China (2020YFA0112701), the National Natural Science Foundation of China (82171057), the Science and Technology Program of Guangzhou, China (202206080005) Guangdong Provincial Medical Science and Technology Research Fund (A2022301), and Kashi Regional Science and Technology Project (KS2021048). The funding bodies had no role in the design of the study; nor in the collection, analysis, or interpretation of the data; nor in writing the manuscript.