Twenty Novel MicroRNAs in the Aqueous Humor of Pseudoexfoliation Glaucoma Patients

The microRNAs (miRNAs) are short non-coding RNAs (19–25 nt) that regulate the level of gene expression at the post-transcriptional stage. Altered miRNAs expression can lead to the development of various diseases, e.g., pseudoexfoliation glaucoma (PEXG). In this study, we assessed the levels of miRNA expression in the aqueous humor of PEXG patients using the expression microarray method. Twenty new miRNA molecules have been selected as having the potential to be associated with the development or progression of PEXG. Ten miRNAs were downregulated in PEXG (hsa-miR-95-5p, hsa-miR-515-3p, hsa-mir-802, hsa-miR-1205, hsa-miR-3660, hsa-mir-3683, hsa -mir-3936, hsa-miR-4774-5p, hsa-miR-6509-3p, hsa-miR-7843-3p) and ten miRNAs were upregulated in PEXG (hsa-miR-202 -3p, hsa-miR-3622a-3p, hsa-mir-4329, hsa-miR-4524a-3p, hsa-miR-4655-5p, hsa-mir-6071, hsa-mir-6723-5p, hsa-miR-6847-5p, hsa-miR-8074, and hsa-miR-8083). Functional analysis and enrichment analysis showed that the mechanisms that can be regulated by these miRNAs are: extracellular matrix (ECM) imbalance, cell apoptosis (possibly retinal ganglion cells (RGCs)), autophagy, and elevated calcium cation levels. Nevertheless, the exact molecular basis of PEXG is unknown and further research is required on this topic.


Introduction
Glaucoma is a group of multifactorial eye diseases that result in progressive irreversible damage to the optic nerve and cause vision loss, and is considered the second leading cause of blindness in the world [1]. The basic process in the pathophysiology of glaucoma is the loss of retinal ganglion cells (RGCs) and their axons, which leads to changes in the morphology of the optic disc and consequently to visual field defects. The main risk factor for glaucoma is increased intraocular pressure (IOP), caused by an imbalance between the production and outflow of aqueous humor (AH) from the anterior chamber [2]. The most important outflow pathway is trabecular meshwork (TM), which consists of trabecular cells surrounded by the extracellular matrix (ECM). Dysfunction in the normal homeostatic process leads to increased outflow resistance and elevated IOP [3]. Lowering IOP is the only treatment that prevents the progression of vision loss [4].
MiRNAs are short non-coding RNAs (19-25 nt) that regulate the level of protein expression at the post-transcriptional level. These molecules are responsible for the regulation of many physiological cellular processes, including division, cell differentiation, and regeneration, but also for many pathological processes such as neoplastic transformation. MiRNAs are found both in the cell and in extracellular fluids, i.e., plasma, urine, and aqueous humor [27][28][29][30]. MiRNAs are transcribed in a similar way to protein genes. The first step in the expression of a gene encoding miRNA is the synthesis of primary miRNA (several hundred nucleotides long), which are then modified in the nucleus at both ends of the RNA. The next stage of miRNA maturation takes place in the cytoplasm and, as a consequence, a pre-miRNA of about 70 nucleotides is formed, followed by the production of mature micro RNA. The functioning of micro RNA is closely related to the RNA-induced silencing complex (RISC). Both molecules linked together to interact with the target mRNA through partial complementarity of the sequence of one of the miRNA strands, leading to degradation of the target mRNA through the activity of the RISC complex. By analyzing the level of miRNA expression, it is possible to draw conclusions about the pathomechanism of the studied diseases and to use them as biomarkers of these diseases [27,31,32]. The most important features of biomarkers are: specificity for the disease entity, the possibility of quick, easy, and safe determination, and the possibility of obtaining accurate results [33][34][35][36].
The aim of the study was to assess the levels of miRNA expression in a group of patients with PEXG and to look for mechanisms caused by altered expression of miRNA that may be related to this disease.

Study Groups
Eighteen western descent patients underwent routine cataract surgery at the Department of Diagnostics and Microsurgery of Glaucoma, Medical University of Lublin (Poland) (nine PEXG and nine age-matched control patients). Prior to the surgery, all of the patients underwent the ophthalmic examination including best-corrected visual acuity (BCVA) assessed with Snellen charts, an IOP test measured by Goldmann applanation tonometry (GAT), gonioscopy, a slit-lamp examination, indirect ophthalmoscopy after pupil dilation with a stereoscopic optic nerve head assessment, optical coherence tomography (Carl Zeiss Cirrus HD-OCT 5000), and standard automated perimetry (SAP) (24-2 strategy using the Humphrey perimeter).
Inclusion criteria in the study group: (1) Incipient senile cataract: cloudy area in the lens that causes decreased vision. In both groups, other causes of RNFL thinning (i.e., myopia, optic disc anomalies, ischemic optic neuropathy, optic neuritis, hereditary optic neuropathy, traumatic optic neuropathy, multiple sclerosis, and degenerative diseases such as Alzheimer's and Parkinson's disease) were excluded. Moreover, all of the participants had no previous intraocular surgeries, no other systemic or ocular diseases known to affect the visual field (e.g., pituitary lesions, demyelinating diseases, diabetes mellitus, etc.), and no previous eye or head trauma in their medical history.

Aqueous Humor
Approximately 100 µL of AH were obtained from the eye's anterior chamber at the beginning of the cataract surgery with special care to avoid contamination with blood or tears.

RNA Isolation
The total RNA was isolated from the AH samples using an miRNeasy Serum/Plasma Kit (Qiagen, Valencia, CA, USA) according to the manufacturer's instructions. The isolated RNA was stored at −80 • C for further analysis. The RNA concentration was determined by using a NanoDrop 2000c spectrophotometer (Thermo Fisher Scientific, Waltham, MA, Cells 2023, 12, 737 4 of 21 USA). In addition, RNA analysis was performed using an Agilent Bioanalyzer 2100 (Agilent Technologies, Santa Clara, CA, USA) and a Pico RNA Kit according to the manufacturer's protocol.

The miRNA Profiling
A microarray system (GeneChip miRNA 4.0 Array chip, Affymetrix, Santa Carla, CA, USA) was used to determine the miRNA expression profiles. The RNA preparation and hybridization were performed according to the manufacturer's protocol, with one modification to extend the hybridization time to 42 h. The gene chips were scanned with an Affymetrix GeneChip Scanner 3000 (Affymetrix, Santa Carla, CA, USA).
The raw data, in a CEL format, were analyzed as log2-transformed intensities using the Affymetrix Transcriptome Analysis Console (TAC), following the software's guidelines to determine differentially expressed genes (DEGs) between the PEXG and control patients. The miRNAs that showed at least a one-fold difference between the control and PEXG groups, with a p-value of less than 0.01, were considered statistically significant and differentially expressed.
The analysis of the interactions of the selected miRNAs with genes was carried out in an R environment (version 4.1.2, https://www.r-project.org (accessed on 1 December 2021)) using the multiMiR 1.16.0 package (https://bioconductor.org/packages/release/ bioc/html/multiMiR.html, (accessed on 2 September 2022)) [39]. The enrichment analysis was performed using the mirNET database (https://www.mirnet.ca/miRNet/home.xhtml (accessed on 2 September 2022)), where a list of selected miRNAs and a list of genes with which they interact were introduced [40]. Terms of the Kyoto Encyclopedia of Genes and Genomes (KEGG), REACTOME, and Gene Ontology (GO) categories were searched for in PEX-related pathways.
The plot visualizing the enrichment analysis was generated using the ggplot2 3.3.0 package in the R environment.

Statistical Analysis
Results at p < 0.05 were considered statistically significant. The following software was used: GraphPad Prism 9 (Graph Pad Software, San Diego, CA, USA) and Statistica 13.3 (StatSoft, Krakow, Poland). The analyzes included the Shapiro-Wilk test (to assess the compliance of the examined variables with the normal distribution), the student's ttest, and ROC curves (area under the curve (AUC) were calculated) for variables on a quantitative scale (data presented as an average). Data on a qualitative scale are presented as numbers and percentages and were analyzed using the Chiˆ2 test.

Results
Analysis comparing the expression level of miRNAs using expression microarrays in the PEXG group compared to the cataract group showed 20 miRNAs as DEG's.
Using the DIANA miRPath v3.0 database, a functional analysis of each of the selected miRNAs was performed (the analysis was performed separately for miRNAs with reduced and increased expression in PEXG). Next, networks were constructed visualizing possible functional relationships between the selected miRNAs (Figures 3 and 4).    The functional network for downregulated miRNA functionally linked together nine miRNAs (leaving one type of miRNA without functional linkage-hsa-miR-3683). The terms with the most associations between the selected miRNAs were: "Ion binding", "organelle", "biosynthetic process", "cellular nitrogen compound metabolic process", "molecular function" and "molecular function". The functional network for downregulated miRNA functionally linked together nine miRNAs (leaving one type of miRNA without functional linkage-hsa-miR-3683). The terms with the most associations between the selected miRNAs were: "Ion binding", "organelle", "biosynthetic process", "cellular nitrogen compound metabolic process", "molecular function" and "molecular function".
Functional network for upregulated miRNA functionally linked eight miRNAs (leaving two types of miRNAs (has-miR-6723-5p, and hsa-miR-4655-5p) without functional links between other miRNAs and with each other). The terms that obtained the most connections between the selected miRNAs were: "Ion binding", "organelle", "biosynthetic process", "cellular nitrogen compound metabolic process", "steroid hormone biosynthesis" and "gene expression".
The next step was to perform an enrichment analysis to evaluate the biological processes regulated by the selected miRNAs and interacted genes using the mirNET database. Genomes), GO: Biological Processing (GO:BP), GO:Cellular Compartment (GO:CC), GO: Molecular Function (GO:MF), and Reactome. Genes indicated as targets of the expressed miRNAs were associated with intracellular and extracellular structures in addition to various processes. The top 10 terms for each category that may be potentially related to glaucoma and/or PEX are presented in Figures 5 and 6, separately for downregulated and upregulated miRNAs in PEXG.  for enrichment adjusted by the Benjamini correction for multiple-hypothesis testing. The number in brackets following the names of terms indicates the number of associated genes. The plot was generated using the ggplot2 3.3.0 package in the R environment. Within the group of downregulated miRNAs, the most frequent terms are related to the regulation of transcription and translation, internal cellular transport ("trans-Golgi network transport vesicle", "endocytosis"), and signal transduction ( Figure 5). In the upregulated miRNA group, the most frequently repeated terms are related to the regulation of transcription and translation, the regulation of the cell cycle, and the activity of cations (mainly calcium and zinc) ( Figure 6).
A PPI network of the genes regulated by selected miRNAs were constructed and clustering analysis was performed based on these PPI networks. In total, six clusters (sub-networks) were created for each PPI network of the genes regulated by downregulated miRNAs (Figure 7, Table 4). Clusters I and IV are functionally related to chromatin organization and epigenetic regulation of gene expression, moreover, the key genes in them are the EZH2 gene (in cluster I) and BAZ1B (in cluster IV). Cluster II is functionally related to nuclear chromosome segregation and cellular response to DNA , REACTOME, and KEGG (Kyoto Encyclopedia of Genes and Genomes) categories, revealed for glaucoma-and/or pseudoexfoliation syndrome (PEX)-associated genes targeted by upregulated miRNAs found in the current study in the aqueous humor (AH) of patients with PEXG and age-matched controls. p-value-EASE score for enrichment adjusted by the Benjamini correction for multiple-hypothesis testing. The number in brackets following the names of terms indicates the number of associated genes. The plot was generated using the ggplot2 3.3.0 package in an R environment.
Within the group of downregulated miRNAs, the most frequent terms are related to the regulation of transcription and translation, internal cellular transport ("trans-Golgi network transport vesicle", "endocytosis"), and signal transduction ( Figure 5). In the upregulated miRNA group, the most frequently repeated terms are related to the regulation of transcription and translation, the regulation of the cell cycle, and the activity of cations (mainly calcium and zinc) ( Figure 6).
A PPI network of the genes regulated by selected miRNAs were constructed and clustering analysis was performed based on these PPI networks. In total, six clusters (sub-networks) were created for each PPI network of the genes regulated by downregulated miRNAs (Figure 7, Table 4). Clusters I and IV are functionally related to chromatin organization and epigenetic regulation of gene expression, moreover, the key genes in them are the EZH2 gene (in cluster I) and BAZ1B (in cluster IV). Cluster II is functionally related to nuclear chromosome segregation and cellular response to DNA damage (the key gene: TIMELESS). Cluster III is functionally related to the voltage-gated calcium channel complex and MAPK signaling pathway (the key gene: CACNA2D1). Cluster V is functionally related to protein transport (the key gene: RAB7A). Cluster VI is functionally related to ubiquitinization (the key gene: UBE2D1).
Cells 2023, 11, x damage (the key gene: TIMELESS). Cluster III is functionally related to the vo calcium channel complex and MAPK signaling pathway (the key gene: CA Cluster V is functionally related to protein transport (the key gene: RAB7A). C functionally related to ubiquitinization (the key gene: UBE2D1).     In total, six clusters (sub-networks) were created for each PPI network of the genes regulated by upregulated miRNAs (Figure 8, Table 5). Cluster VII is functionally related to cytokine-cytokine receptor interaction and cell-cell adhesion (the key gene: IL10). Cluster VIII is functionally related to the activity of RNA polymerases II and III (the key gene: PLR2D). Cluster IX is functionally related to ubiquitinization (the key gene: PHC1). Cluster X is functionally related to L-type voltage-gated calcium channel complex (the key gene: CACNA1G). Cluster XI is functionally related to the regulation to the jnk cascade (the key gene: AXIN1). Cluster XII is functionally related to Rho GTPases activate ROCKs and phospholipase D signaling pathway (the key gene: RHOA). In total, six clusters (sub-networks) were created for each PPI network of the genes regulated by upregulated miRNAs (Figure 8, Table 5). Cluster VII is functionally related to cytokine-cytokine receptor interaction and cell-cell adhesion (the key gene: IL10). Cluster VIII is functionally related to the activity of RNA polymerases II and III (the key gene: PLR2D). Cluster IX is functionally related to ubiquitinization (the key gene: PHC1). Cluster X is functionally related to L-type voltage-gated calcium channel complex (the key gene: CACNA1G). Cluster XI is functionally related to the regulation to the jnk cascade (the key gene: AXIN1). Cluster XII is functionally related to Rho GTPases activate ROCKs and phospholipase D signaling pathway (the key gene: RHOA).   IFNG, TNFSF9, IL11, TNFSF4, TNFSF12, CXCL5,  IL16, CD40LG, TNFRSF1B, IL10 GO:BP Positive regulation of cell-cell adhesion 8.60 × IFNG, TNFSF9, TNFSF4, ICOS, CD40LG, IL10 Figure 8. Results of protein-protein interaction analysis. The six top-scored clusters separated from PPI network of the genes regulated by upregulated miRNAs found in the current study. The hub genes identified using the cytoHubba app included the EZH2, H2AFX in hub I, CACNA1D in hub IIa, and MYC in hub IIb (Figure 9). Enrichment analysis for these networks indicated that hub I is related to generic transcription pathway and cellular responses to stress, hub IIa is related to voltage-gated channel activity and MAPK signaling pathway, hub IIb is related to cellular senescence and Wnt, Hippo, JAK-STAT signaling pathways ( Table 6).
The hub genes identified using the cytoHubba app included the EZH2, H2AFX in hub I, CACNA1D in hub IIa, and MYC in hub IIb (Figure 9). Enrichment analysis for these networks indicated that hub I is related to generic transcription pathway and cellular responses to stress, hub IIa is related to voltage-gated channel activity and MAPK signaling pathway, hub IIb is related to cellular senescence and Wnt, Hippo, JAK-STAT signaling pathways (Table 6).  CACNG8, CACNA1D, CACNB2, CACNA1E, CACNA1B, CACNA2D2, CACNA2D3, CACNB4

Discussion
MiRNA molecules are found both inside the cell and in extracellular fluids. The miRNA in the AH may be derived from blood plasma and may be secreted by intraocular cells. These molecules have been shown to be of great importance in processes such as retinal homeostasis and its development [44,45].
Current knowledge about the pathological mechanisms of PEXG is still unknown. Many studies suggest the role of non-coding RNAs, which, through altered expression levels, may be associated with the development of PEXG. The role of miRNAs [26,[46][47][48] and snoRNA [49] is mainly suggested.
In our work, we examined the level of miRNA expression in a group of patients with PEXG and cataracts. After selecting the 20 DEGs of miRNAs (10 upregulated and 10 downregulated), we performed an enrichment and functional analysis.
One of the pathological mechanisms related to PEX contains too high of a concentration of calcium cations, which is directly related to the activity of the channels for these cations. Currently, three types of calcium channels are cited, i.e., CACNA1A [21], CACNB3 [46,49,50], and CACNB4 [50]. The calcium channel is composed of subunits, of which the α subunit forms a transmembrane channel, while the remaining subunits (β and α2/δ) perform regulatory functions. Additionally, in some channels, there may be a γ subunit. The calcium channels are divided into three subfamilies, i.e., CaV1, CaV2, and CaV3. The activity of the CaV1 subfamily channels are responsible for the initiation of contraction and regulates gene expression and transmission in synapses. The CaV2 subfamily channels are associated with synaptic transmission, while the CaV3 subfamily channels are associated with the production of action potentials in myocytes and neurons. The CaV1.4 transporter has been shown to be present in the retinal cell membranes of the eye and its activity is related to visual signaling [51,52]. The results of functional analysis in our study showed that hsa-miR-8074 influences the transport of calcium cations across biological membranes through the regulation of the ATP2B2 gene. Moreover, we found that the voltage-gated calcium channel activity is being regulated by the genes targeted by selected miRNAs. Additionally, the Wnt signaling pathway is regulated by the genes targeted by selected miRNAs, where one of the signaling pathways is the Wnt/Ca2 + pathway, the activity of which is related to G proteins and phospholipases.
An important regulatory structure of the IOP is HTM (human trabecular meshwork), which is made up of trabecular cells surrounded by the ECM. HTM is responsible for the outflow of AH-in a situation of difficult outflow, the IOP increases [53,54]. The literature provides two pathways (TGFβ1 and TGFβ2), in which the changes lead to an imbalance in the eye's ECM and, consequently, to an increase in IOP [53,54]. Among the results of the functional analysis in our work, we obtained the term "ECM receptor interactions" for four miRNAs that are regulated by genes: hsa-mir-802 (COL5A1, LAMC1, FN1, CD47), hsa-mir-3660 (LAMA3), hsa-mir -3662a-3p (ITGB8, SV2B, COL24A1, COL4A6, TNR), and hsa-miR-4774-5p (COL11A1, COL24A1). In addition, proteins are a very important element of the ECM. Disruption of protein metabolism (synthesis, post-translational modification, or transport) may lead to the development of PEX [55].
Functional analysis of selected miRNAs showed that a relatively large number of these miRNAs influence genes related to protein metabolism. The process of protein modification in the cell is regulated by genes whose expression is influenced by six out of twenty selected miRNAs (hsa-miR-95-5p, hsa-miR-1205, hsa-miR-3660, hsa-miR-515-3p, and hsa-miR-6071,22-3p). In addition, subsequent miRNAs affect the expression level of genes related to protein metabolism, i.e., "protein complex" (hsa-miR-4329, hsa-miR-6071, hsa-miR-1205, and hsa-miR-3660), "cytoskeletar protein binding" (hsa-miR-6847-5p, 1205), "post-translational protein modification" (hsa-miR-7843-3p), "endoplasmic reticulum-Golgi intermediate compartment" (has-miR-802), and "cytoplasmic sequestering of protein "(Hsa-miR-4655-5p). Moreover, hsa-miR6723-5p regulates the expression level of the UBE2H gene related to ubiquitin-mediated proteolis. Among the enriched functional terms of genes targeted by selected miRNAs in AH in PEXG patients, it is worth mentioning the "trans-Golgi network transport vesicle" as a process that influences protein transport. In addition, many times the results of this analysis were terms related to the process of protein biosynthesis ("positive regulation of translation", "regulation of translation") and transcription ("negative regulation of transcription from RNA polymerase II promoter", "negative regulation of transcription, DNA-dependent "). In a previous work, we suggested that some snoRNA molecules are related to PEXG [49]. One of the terms obtained in enrichment analysis is "spliceosomal complex" which may suggest that the selected miRNAs may also be responsible for this process.
The mechanisms discussed above largely explained where the elevated IOP comes from. Elevated IOP results in RGC apoptosis and, furthermore, optic nerve degeneration which is often observed in late PEXG [56][57][58]. Another mechanism that can lead to IOP increase and cell death in glaucoma is oxidative stress [59][60][61][62]. The hsa-miR-802 molecule by regulating the expression of MCL1, PPP2CB, and SOD2 genes influences the processes taking place in the mitochondria during apoptosis. Among the enriched functional terms of genes targeted by selected miRNAs in AH in PEXG patients, it is worth mentioning "G0 and early G1" and "G2/M transition of the mitotic cell cycle". In addition, "AKT phosphorylates targets in the cytosol" and "PI3K/AKT activation" are terms that suggest that selected miRNAs may regulate the survival of ocular neurons and thus may lead to optic nerve degeneration. The role of the PI3K/Akt pathway was previously mentioned [63][64][65][66][67].
Hsa-miR-802 shows an increased expression in primary human umbilical vein endothelial cells (HUVECs) in response to induced hypoxia [67]. As above, hsa-miR-202-3p is one of the miRNAs whose expression levels have been altered in lung adenocarcinoma A549 cells by hypoxia [68].
Additionally, a correlation between the risk of developing normal-tension glaucoma (NTG) and the level of retinol in the blood serum was demonstrated [69]. Moreover, another group of scientists reports that retinoic acid (a retinol derivative) lowers the expression level of the gene which codes the small heat shock protein B8 (HSPB8) [70], which is a chaperone protein functionally related to chaperone-assisted selective autophagy (CASA) [71]. One of the miRNA molecules selected as DEGs in our work, hsa-miR-6509-3p regulates the expression level of the HSPB8 gene. The above data may suggest a potential influence of the autophagy process in the development of PEXG.
The results of the PPI analysis suggest and confirm that a very important mechanism of PEXG may be disturbances in the transport of calcium cations in the eye. One of the hub genes that have been selected is the CACNA1D gene.
Our research had some limitations. First, the study groups were small, consisting of nine patients each. A study involving a larger group of patients is needed to confirm the results. Furthermore, PEXG patients had advanced neuropathy, which suggests that miRNA molecules may be the result of the pathological process and not a causative factor. However, we precluded that the same stage of the disease within the group would enable more unique results. Since the advanced stages of PEXG, comparison of miRNAs in the early stages of PEXG is needed. In addition, all enrolled patients were diagnosed with senile cataracts. Therefore, the impact of cataracts on the expression of miRNA in the PEXG group should be taken into consideration, but in the studied group all patients also had similar stages of cataracts. However, it is ethically impossible to collect aqueous humor in healthy eyes due to the invasive nature of the procedure.
Until now, the scientific literature has not provided such information, which means that this is the first time we are doing it. These mechanisms are possible, but their role in the pathogenesis of PEXG/PEX must be carefully studied.

Conclusions
Our results suggest 20 new miRNAs whose altered expression in the aqueous humor may have an influence on the pathomechanisms of PEXG. Functional and enrichment analysis confirms the mechanisms, such as high levels of calcium cations, ECM imbalance, apoptosis of RGC cells, and the optic nerve damage or autophagy as possible causes of PEXG. Moreover, the PPI analysis suggests that CACNA1D, EZH2, and MYC may be hub genes for PEXG. However, confirmation of these conclusions requires further study.