Metabolomic Profiling of Open-Angle Glaucoma Etiologic Endotypes: Tohoku Multi-Omics Glaucoma Study

Abstract

Purpose: The purpose of this study was to investigate biologically meaningful endotypes of open-angle glaucoma (OAG) by applying unsupervised machine learning to plasma metabolites.

Methods: This retrospective longitudinal cohort study enrolled consecutive patients aged ≥20 years with OAG at Tohoku University Hospital from January 2017 to January 2020. OAG was confirmed based on comprehensive ophthalmic examinations. Among the 523 patients with OAG with available clinical metabolomic data, 173 patients were longitudinally followed up for ≥2 years, with available data from ≥5 reliable visual field (VF) tests without glaucoma surgery. We collected fasting blood samples and clinical data at enrollment and nuclear magnetic resonance spectroscopy to profile 45 plasma metabolites in a targeted approach. After computing a distance matrix of preprocessed metabolites with Pearson distance, gap statistics determined the optimal number of OAG endotypes. Its risk factors, clinical presentations, metabolomic profiles, and progression rate of sector-based VF loss were compared across endotypes.

Results: Five distinct OAG endotypes were identified. The highest-risk endotype (endotype B) showed a significant faster progression of central VF loss (P = 0.007). Compared with patients with other endotypes, those with endotype B were more likely to have a high prevalence of dyslipidemia, cold extremities, oxidative stress, and low OAG genetic risk scores. Pathway analysis of metabolomic profiles implicated altered fatty acid and ketone body metabolism in this endotype, with 34 differentially enriched pathways (false discovery rate [FDR] < 0.05).

Conclusions: Integrated metabolomic profiles identified five distinct etiologic endotypes of OAG, suggesting pathological mechanisms related with a high-risk group of central vision loss progression in the Japanese population.

Primary open-angle glaucoma (POAG) is the most prevalent form of glaucoma, resulting in irreversible blindness.¹ Its only confirmed risk and progression factor is increased intraocular pressure (IOP); however, growing evidence suggests that POAG is genetically and environmentally heterogeneous in nature.^2–4 Patients with glaucoma exhibit variations in genomics,⁵ transcriptomics,⁶ proteomics,^7,8 cytokine,⁹ and metabolomic profiles¹⁰; these variations might be linked to systemic profiles and could contribute to the heterogeneous pathophysiology in glaucoma. 

POAG can present vulnerability in differing retinal ganglion cell (RGC) subpopulations, leading to heterogeneity in clinical presentation and outcomes.^11,12 The POAG subtype with early paracentral visual field (VF) loss owing to damage to the papillomacular nerve fiber bundle directly leads to decreased visual acuity.¹³ Paracentral nerve fibers are anatomically vulnerable owing to shear forces and local blood flow from the optic disc and accompanying blood vessels.¹⁴ Additionally, RGC axons are unmyelinated and have the longest node of Ranvier in the body, leading to a high mitochondrial energy demand.¹⁵ These features indicate a particular role of altered metabolism in the etiology of POAG with central VF loss. 

Metabolomics is being increasingly used for glaucoma research.¹⁰ A systematic review of 18 metabolomic studies identified several biomarkers closely linked to distinct biological mechanisms.¹⁰ Furthermore, a recent large-scale nested case-control study suggested the role of dysregulated lipid metabolism in POAG.¹⁶ However, these studies investigated different metabolomes by comparing glaucoma and non-glaucoma populations; no study has yet characterized distinct glaucoma endotypes or determined their relationships with clinical outcomes. 

Here, we analyzed metabolomics data of patients with open-angle glaucoma (OAG), including those with POAG with IOP >21 millimeters of mercury (mm Hg; high-tension glaucoma [HTG]) and normal-tension glaucoma (NTG) with IOP ≤21 mm Hg, to identify biologically distinct OAG endotypes by applying unsupervised machine learning (ML) and to explore relationships of metabolomic-driven endotypes with clinical, genomic, and longitudinal VF loss progression. 

Materials and Methods

Study Design and Population

We retrospectively reviewed medical records and baseline clinical data of enrolled consecutive patients aged ≥20 years with glaucoma treated at Tohoku University Hospital (Miyagi, Japan) between January 2017 and January 2020. 

After the initial diagnosis by a glaucoma specialist (author T. Nakazawa), the patients underwent comprehensive ophthalmic evaluations, including slit-lamp, gonioscopic, and dilated funduscopic examinations; IOP assessments; and swept-source optical coherence tomography (SS-OCT) and full static threshold perimetry tests. Glaucoma was confirmed based on the presence of glaucomatous optic neuropathy accompanied by reliable VF defects. For OAG confirmation, the following criteria were required: (a) open angle with no signs of exfoliation material according to gonioscopy/slit-lamp evaluation and (b) optic disc abnormality (i.e. vertical cup-to-disc ratio > 0.7) or retinal nerve fiber layer defects and reproducible glaucomatous VF defects. For patients with bilateral glaucoma, the eye with the greater mean deviation (MD) plot in the perimetry was used in all analyses. 

The study protocol was approved by the institutional review boards of Tohoku University Graduate School of Medicine (approval numbers: 2021-1-184 and 2202-1-1184) and complied with the tenets of the Declaration of Helsinki. All participants gave written informed consent. 

Systemic and Ocular Examinations

Details are provided in the Supplementary Methods. Briefly, structured interviews and chart reviews were conducted at patient enrollment; and trained technicians conducted systemic and ocular examinations. 

Metabolomic Profiling

Details regarding blood sample collection were previously reported.¹⁷ Briefly, blood samples were collected immediately using vacutainer tubes containing EDTA-2Na (Venoject II; Terumo Corporation) and centrifuged at 2330 × g for 10 minutes at 4°C. The aliquoted plasma was stored at −80°C in a MATRIX 2D screw tube (Thermo Fisher Scientific). Metabolomic profiling of 200 µL plasma per sample was performed using a standard methanol extraction procedure¹⁷; processed samples were transferred to a 3-mm Bruker SampleJet nuclear magnetic resonance (NMR) tube. We performed NMR experiments at 298 K on a Bruker 600 MHz spectrometer. The large interfering NMR signal arising from water in all biofluids is eliminated by using the appropriate standard NMR solvent suppression methods.¹⁸ For internal standardization, we have used the use of a known concentration of DSS (2,2-dimethyl-2-silapentane-5-sulfonate sodium salt) in the aqueous media. Using the target profiling approach in the Chenomx Profiler module of Chenomx NMR Suite (Chenomx version 8.x), 45 metabolites were identified and quantified. Standard 1D NOESY and CPMG (Carr-Purcell-Meiboom-Gill) spectra were obtained for each plasma sample; spectra were acquired with 64 scans and 32 k of complex data points. Standard 1D NOESY spectra were utilized for the identification and quantification of metabolites. Additionally, 1D CPMG spectra were used to minimize the impact of residual proteins on the quantification process. To further verify metabolite identifications, 2D TOCSY and ¹H,¹³C-HSQC experiments were conducted on several samples. Spiking experiments were also carried out to confirm these identifications. 

Genotyping, Imputation, and Genetic Risk Score Calculation for OAG

Details regarding genotyping and imputation are provided in the Supplementary Methods. Briefly, genomic DNA extracted from whole-blood samples was genotyped using the Japonica Array NEO.¹⁹ 

The OAG Genetic Risk Score (GRS) was calculated as follows:  

\begin{eqnarray*} GR{S_i} = \mathop \sum \limits_j ({x_{i.j}} \times {w_j}) \end{eqnarray*}

where GRS_i represents the GRS of subject i; x_i.j represents the number of non-reference alleles carried by subject i of variant j; and w_j represents a weight parameter of variant j. These were derived from weight parameters of 98 of 127 OAG-associated variants reported by Gharahkhani et al.⁵ available in our dataset. Weight parameters based on those in the Biobank Japan Project OAG-genomewide association study (GWAS)²⁰ were used in meta-GWAS.⁵ 

Outcome Measures

The SITA standard 24-2 program of the Humphrey field analyzer (Carl Zeiss Meditec) was used to determine MD and total deviation (TD). As previously determined,¹⁴ the primary outcome was the rate of glaucomatous VF loss progression (averaged TD slope, dB/year) in each of the six classified sectors (Supplementary Fig. S1).²¹ To analyze sector-based VF progression, we restricted the cases with a reproducible glaucomatous VF detect, including ≥5 reliable VF measurements with ≥2 years of follow-up post-enrollment. We omitted cases with unreliable VFs/VF defects presumably due to other causes (i.e. posterior cortical cataract or retinal diseases; see the Supplementary Methods). In an exploratory analysis, we further analyzed the regional thickness changes (µm/year) of the circumpapillary retinal nerve fiber layer (cpRNFL; overall and quadrants) and macular ganglion cell layers (overall and upper/lower areas). Stratified analysis by glaucoma subtypes (HTG versus NTG) was further conducted. 

Statistical Analyses

Statistical methods were not used to predetermine the sample size, as this was an exploratory analysis. We performed log₂ transformation and autoscaling for each metabolite²² and computed a distance matrix using Pearson distance.²³ To derive biologically distinct OAG endotypes, we applied partitioning around medoids. To choose an optimal number of clusters, we computed various unsupervised ML methods including silhouette, elbow, and gap statistics (Supplementary Fig. S2).²⁴ The final number of molecular endotypes was determined based on gap statistics. To visualize the metabolome endotypes, the t-distributed stochastic neighbor embedding method was applied (Supplementary Fig. S3). 

Between-endotype differences in patient characteristics and clinical presentations were determined using analysis of variance, Kruskal–Wallis, chi-squared, and Fisher's exact tests, as appropriate. We also visualized relationships between major clinical characteristics, including a priori known glaucoma progression risk factors, and metabolomic-driven endotypes using a chord diagram, a Venn diagram, and an upset plot. We used the circlize, VennDiagram, and ComplexUpset packages for the chord diagram, Venn diagram, and upset plot, respectively. Each cutoff point represented the median value of participants. 

To determine the association of endotypes with VF loss progression risk, we compared the average TD slope in each sector using the Kruskal–Wallis test. 

To identify biologically meaningful pathways, we used MetaboAnalyst's biomarker, functional, and enrichment analysis modules to compare the highest risk of endotype and the other endotypes. A functional class scoring analysis was conducted to investigate if metabolites for specific biological pathways were enriched. Although the statistical power was limited, we examined the robustness of the enrichment analysis among those with longitudinal VF data (n = 173) as sensitivity analyses. 

To assess the robustness of the causal link between endotype B and central VF loss progression to unmeasured confounding, we additionally calculated E-values²⁵ (see the Supplementary Methods). We also used Bayesian methods^26,27 to determine the risk difference of central VF loss progression between endotype B and non-endotype B (see the Supplementary Methods). 

SAS (version 9.4; SAS Institute, Cary, NC, USA) and R version 4.0.3 (R Foundation for Statistical Computing, Vienna, Austria) were used for all statistical analyses. Two-tailed P values < 0.05 indicated statistical significance. The Benjamini–Hochberg false discovery rate (FDR) method was used for multiple testing corrections. 

Results

Among 544 participants with glaucoma with clinical and metabolomic data initially recruited (Fig. 1), we excluded those with exfoliation glaucoma (n = 2), primary angle-closure glaucoma (n = 3), and other types of glaucoma (n = 16). For endotyping analysis, we included 523 patients with OAG whose 45 metabolites were profiled. Across different numbers of endotypes, both the average elbow score and the network modularity were reasonable (k = 3–5; see Supplementary Fig. S2). Based on the gap statistics, we finally selected the five endotypes. 

Figure 1.

 

Analytic workflow of this study. (a) A total of 544 (280 men and 264 women) participants with glaucoma with clinical and metabolomic data were initially recruited. After excluding those with exfoliation glaucoma (n = 2), primary angle-closure glaucoma (n = 3), and other types of glaucoma (n = 16), including childhood and secondary glaucoma, 523 patients with open-angle glaucoma (OAG) were included in this study. (b) Using the 45 plasma metabolites from 523 patients with OAG, we computed a distance matrix using Pearson distance and derived biologically distinct OAG endotypes by applying the partitioning around medoids. We used gap statistics to select an optimal number of profiles and determined five OAG endotypes. We applied the t-distributed stochastic neighbor embedding (t-SNE) method to visualize the metabolome endotypes. (c) We examined the differences in major clinical and genomic variables according to the five endotypes. For visualization, we used chord diagram, Venn diagram, and upset plots. (d) For 173 patients (out of 523) with OAG who had been followed for ≥2 years with results of ≥5 reliable visual field tests, we assessed the progression rate of sector-based glaucomatous visual field loss according to the five endotypes. (e) To examine the relationship of the metabolite profiles of the highest-risk group (highest progression of central visual field loss; endotype B) compared with those of the other groups (non-endotype B), we performed a functional pathway analysis. OAG, open-angle glaucoma; t-SNE, t-distributed stochastic neighbor embedding.

After classifying five endotypes (see Supplementary Fig. S3), we included 173 participants with longitudinal reliable VF loss data (median follow-up, 2.98 ± 1.0 years; see Fig. 1). Patient characteristics were generally similar between those with and without longitudinal data, except for older age, lower rate of NTG, and higher systolic blood pressure and oxidative stress in those without longitudinal VF data (Supplementary Table S1). 

Table 1 summarizes baseline characteristics of the 173 OAG samples according to the five endotypes. Patients with endotype A were male dominant (84%) and were likely to be obese and have hypertension, sleep apnea, cardiovascular diseases, and diabetes but were less likely to have migraine and cold extremities. Most patients with endotype C were female patients and had a higher prevalence of migraine, cold extremities, and oral contraceptive use. There was a difference in ophthalmic parameters among the five endotypes. Axial length was the highest in endotype A and the lowest in endotype D. IOP was the highest in endotype C and the lowest in endotype D. The central corneal thickness was relatively higher in endotypes A and E and the lowest in endotype D. Additionally, oxidative stress levels were higher in endotypes B and C. There was a significant difference in OAG genetic risk across endotypes (P = 0.02). The highest GRS was 0.84 ± 0.41 for endotype D, and the lowest was 0.51 ± 0.39 for endotype E. Relationships between OAG progression risk factors and the five endotypes are shown in Figure 2. 

Table 1.

 

Baseline Characteristics and Clinical Features of Participants According to Open-Angle Glaucoma Endotype (n = 173)

 

Figure 2.

 

Relationship between previously known risk factors for open-angle glaucoma progression and endotypes characterized by metabolites (n = 173). (a) Chord diagram showing the previously known risk factors for open-angle glaucoma (OAG) progression by metabolic-driven endotypes. The ribbons connect individual endotypes to previously known clinical and endotypes. The widths of the ribbon represent the proportion of patients with OAG within the endotype with the corresponding clinical and metabolomic characteristics. Then, it was scaled to a total of 100%. (b) Venn diagram of previously known risk factors for OAG progression and their intersections. The Venn diagram illustrates the composition of five clinical variables and their intersections. The numbers correspond to the number of patients with OAG in each subset and intersection. The cutoff points for age (56 years), central corneal thickness (CCT; 511 µm), pattern standard deviation (PSD; 12%), and intraocular pressure (IOP; 14 mm Hg) were the median values. (c) Upset plot corresponding to the presented Venn diagram. The plot illustrates the composition of five previously known OAG risk factors and their intersections visualized based on the five endotypes. Vertical stacked bar charts reflect the number of patients with glaucoma within each subset and intersection colored according to the endotypes. Horizontal bars indicate the number of patients with glaucoma in each clinical variable set. Black dots indicate the sets of subsets and intersections, and connecting lines indicate relevant intersections related to each stacked bar chart. CCT, central corneal thickness; IOP, intraocular pressure; OAG, open-angle glaucoma; PSD, pattern standard deviation.

To determine the VF progression rate, the MD and TD slopes for each Garway-Heath sector (see Supplementary Fig. S1) were compared. There was a significant difference between the TD-inferior (sector 5) and TD-central (central sector) slopes, where the progression rate was the highest in endotype B (Table 2). Although the sample size was limited, this trend was generally similar, irrespective of NTG or HTG (Supplementary Table S2). Similarly, thinning of the superior cpRNFL and macular ganglion cell layer (GCL) in the upper region was the greatest in endotype B (Supplementary Table S3). 

Table 2.

 

Association of Open-Angle Glaucoma Endotypes With Progression Rate of Overall and Sector-Specific Visual Field Loss (n = 173)

 

To better capture differences between the highest-risk (endotype B) and other endotypes, metabolomic profiles were compared (Supplementary Table S4). The overall metabolomic profile of these 2 endotypes was significantly different, with 12 differentially expressed metabolomes (see Supplementary Table S4, Figs. 3a–c). Biologically meaningful pathways were identified, where endotype B had 34 differentially enriched pathways (FDR < 0.05), including upregulated fatty acid biosynthesis and ketone body metabolism. Among these pathways, the amino acid degradation pathway had the highest number of associated metabolites. We used the longitudinal cohort for sensitivity analysis (n = 173), and similar associations were found on comparing endotype B and non-endotype B (Supplementary Fig. S4a–S4c). 

Figure 3.

 

Differential metabolite expression and functional pathway analyses of endotype B v ersu s non-endotype B (n = 523). (a) Heatmap and volcano plot of differentially expressed metabolites. For example (left), we included 45 metabolites with the most significant P values (two-sided raw P values), and the color bar indicates the scaled value of variance-stabilizing transformation. For the volcano plot (right), the threshold of log2 fold change was |0.58| (i.e. ≥|1.5|-fold change) and that of false discovery rate (FDR) was <0.1. Twelve differentially expressed metabolites met these criteria. (b) Functional pathway analysis. There were 34 differentially enriched pathways with a threshold of P value for FDR <0.05. (c) Enrichment analysis of metabolome data. For the Wilcoxon pathway enrichment analysis, we selected the top 25 pathways with the most significant FDRs, and enrichment ratios were shown. FDR, false discovery rate.

The E-values for central VF progression of <−5.5 and <−6.0 dB/year for TD slope in endotype B versus non-endotype B were 7.99 ± 2.30 and 10.2 ± 1.00, respectively (Supplementary Table S5). Bayesian analyses estimated the risk difference for greater central VF progression (<−5.5 and <−6.0 dB/year for TD slope) in endotype B versus non-endotype B was ≥2.87% and ≥5.65%, respectively, with 99% probability (see Supplementary Table S5, Supplementary Fig. S5). 

Discussion

We identified five distinct OAG endotypes, with one demonstrating a significantly higher central VF loss progression rate. Pathway analysis of metabolomics profiling data revealed an upregulation of fatty acid synthesis and altered ketone body metabolism in this endotype. Although our study should be interpreted with caution, we first evaluated molecular endotypes in patients with OAG, suggesting that its heterogeneity partially results from different pathophysiological states. 

Understanding the complex interplay among host characteristics, metabolomics, and clinical morbidities is essential for uncovering disease pathophysiology and exploring new therapeutic targets.^10,16,28,29 With the advancement of ML techniques, integrated analyses of high-throughput omics and clinical outcome data have identified novel endotypes in several diseases, such as cardiovascular diseases,³⁰ diabetes,³¹ and bronchiolitis,³² indicating the ability to substratify high-risk populations. 

Using unsupervised ML, we found five distinct endotypes that could be clinically meaningful. Accumulating evidence suggests a role of vascular dysfunction and subsequent hypoxia affecting the optic nerve, which presumably underlies the relationship between sleep apnea syndrome and glaucoma.³³ Endotype A was associated with male sex, higher body mass index (BMI), and high prevalence of cardiovascular disease and sleep apnea. Aligning with the current findings, sleep apnea syndrome was more common in men, attributable to sex-related differences in upper airway anatomy, fat distribution, responses to neurochemical compounds, and sex hormones.³⁴ In a longitudinal study of 319 Japanese patients with OAG, male sex, high BMI, and sleep apnea syndrome were identified as independent risk factors for VF loss progression in the central-lower sectors.¹⁴ In contrast, endotype C was associated with female sex, lower BMI, cold extremities, and higher oxidative stress. Studies have suggested a link among impaired vascular autoregulation,^2,3,35,36 oxidative stress,⁴ and NTG progression,³⁷ typically observed in Flammer syndrome.³⁸ Vascular dysregulation in OAG was particularly susceptible to impaired nitric oxide signaling,³⁹ as supported by a genetic study showing that intergenic single nucleotide polymorphisms of CAV1/2 were associated with OAG,⁴⁰ higher IOP,⁴¹ and connective tissue diseases.⁴² 

However, we did not observe sex-related differences among the other three endotypes. An important factor associated with glaucoma progression was a lower central corneal thickness,⁴³ although the differences were not statistically significant. Based on the mean GRS, genetic susceptibility for OAG was higher in endotype D and lower in endotypes B and E, suggesting that acquired factors are more strongly affected by pathogenesis in endotypes B and E. Patients with endotype B were more likely to have a higher prevalence of dyslipidemia, cold extremities, and oxidative stress than those with endotype E, which may contribute to faster progression in endotype B than endotype E, although we were unable to determine any specific factors. 

In pathway analysis, fatty acid biosynthesis was significantly upregulated in the endotype with a higher progression rate of central VF defects (endotype B). Epidemiological studies have reported conflicting associations between blood lipid levels and glaucoma.^44–46 A systematic review reported significant positive associations of hyperlipidemia and hypertriglyceridemia with the risk of open-angle glaucoma,⁴⁷ whereas other studies found an inverse⁴⁴ or a null association.⁴⁸ 

Our findings suggest adverse effects of dysregulated lipid metabolism in endotype B, which were generally consistent with those of plasma metabolite studies of US populations and were also replicated in samples from the UK Biobank.¹⁶ Intriguingly, a stronger adverse association between triglycerides and OAG with paracentral VF loss further supports a stronger link between this endotype and systemic factors, particularly abnormal lipid metabolism. Hyperlipidemia contributes to increased blood viscosity,⁴⁸ thereby leading to elevated IOP.⁴⁷ Considering the significant role of blood flow regulation in the optic nerve in OAG pathogenesis, altered ophthalmic hemodynamics due to hypertriglyceridemia may be important. A recent GWAS identified OAG loci near genes related to cholesterol (ABCA1 and CAV1/2)^40,41 and lipid metabolism (ELOVL5).⁴⁹ Proper coordination of lipid metabolism is essential for RGC homeostasis, as supported by a murine study reporting that increased phospholipid synthesis substituting for triglyceride synthesis may promote axon regeneration.⁵⁰ Similarly, metabolomic studies have reported increased acetylcarnitine and phosphatidylcholine levels in optic nerve crush murine models⁵¹ and plasma samples of patients with OAG.⁵² Overall, altered lipid metabolism may play a role in the pathogenesis of glaucoma, particularly in subtypes with central VF defects. 

Owing to the high energy demand in papillomacular nerve fiber bundles, mitochondrial dysfunction caused by altered metabolism has been suggested as a contributing factor in the pathogenesis of OAG with central VF defects.^14,15 Our enrichment analysis indicated a more significant role of propanoate metabolism in endotype B compared with that in other endotypes. Abnormal carbohydrate metabolism has been reported in the trabecular meshwork of glaucomatous eyes.⁵³ Additionally, a recent mitochondrial gene-set analysis reported a stronger association of altered carbohydrate metabolism with NTG than with overall or HTG.⁵⁴ 

Another key process related to mitochondrial impairment in RGC is ketone body metabolism.⁵⁵ In this study, two metabolites (acetone and succinic acid) were differentially associated with endotype B. The neuroprotective effects of ketone bodies as a substitute energy source for glucose are increasingly being recognized in glaucoma because of the higher yield of adenosine triphosphate per unit of oxygen.⁵⁵ A long-term prospective study revealed a decreased incident OAG with paracentral VF loss associated with plant-based low carbohydrate diets.⁵⁶ A genetic pathway analysis study on glaucoma suggested the butanoate pathway related to producing precursors for ketone bodies.⁵⁷ Although the increase in plasma ketone bodies in the high-risk endotype B group remains controversial compared with findings of neuroprotective effects of ketone bodies in glaucoma^55,56 and a plasma metabolomic study of UK Biobank data with low ketone concentrations among patients with glaucoma,¹⁶ all these findings suggest a key role of dysregulated ketone body metabolism in the etiology of glaucoma with central VF progression. 

Altered amino acid metabolism was remarkable in endotype B and was the most prominent for branched-chain amino acid (BCAA) degradation. A recent systematic review of >140 different metabolites identified 7 amino acids, suggesting a major role of amino acid metabolism in the pathophysiology of glaucoma.¹⁰ The potential neuroprotective effects of BCAA include suppressing toxic glutamate levels, enhancing mitochondrial efficiency, and increasing energy production.⁵⁷ In murine glaucoma models, BCAA administration was associated with decreased RGC survival, thereby improving visual function.⁵⁸ A prospective cohort study reported an inverse association of OAG with paracentral VF loss in women.⁵⁹ 

This study had some limitations. First, the study sample did not include healthy controls as it aimed to describe the endotypes of OAG (i.e. descriptive analysis) but not to develop diagnostic biomarkers related to OAG (i.e. association analysis). Although the current descriptive study derives well-calibrated hypotheses, further experiments to examine the diagnostic accuracy of the metabolites observed in this study or association between metabolites and OAG are required. Second, plasma metabolites, which may or may not reflect ocular conditions, were assessed. However, metabolomic studies of glaucoma remain limited to date, and half of the currently available studies were conducted on blood samples.¹⁰ Third, objective measurements of lens opacity were not assessed; non-glaucomatous factors (e.g. cataracts) might affect central VF progression. Furthermore, as with any observational study, the observed association between endotype B and central VF loss progression might have been affected by unmeasured confounding factors. Nonetheless, all cases were examined by board-certified ophthalmologists in a stringent manner, and the relatively large E-values add some level of robustness to our findings. Fourth, our samples included only East Asians with mainly NTG, limiting data generalizability. Additional studies with other ethnic compositions and other types of glaucoma are warranted. Fifth, we could only measure a subset of the metabolites in a targeted approach, which may overlook contributions from unmeasured metabolites affecting other pathways. Nonetheless, the metabolites detected in this study were mostly replicated in published findings from another glaucoma cohort study.¹⁰ Furthermore, the NMR technique used excels in measuring low-molecular and/or high-concentration metabolites, which may minimize false positives. Nonetheless, future comprehensive studies including untargeted metabolomics combined with lipidomics in larger samples are warranted. Sixth, all patients had prevalent glaucoma, and we could not eliminate the possibility that certain metabolites were influenced by medications, particularly for systemic diseases (such as neuroleptics, diabetic medications, or statins). Notwithstanding, the number and type of glaucoma eye drops and oral glaucoma medications (e.g. oral carbonic anhydrase inhibitors) were generally comparable across the five endotypes. Although the inferential goal of this study was not exploring certain metabolomes as predictive biomarkers for glaucomatous progression, future studies are required to examine the causality under the free of systemic medications. Last, the current study lacks the external validation. Thus, we have used bootstrapping to investigate the stability of endotyping. Accordingly, the bootstrap-based stability (with 5000 resamplings) was 0.72 and the mean Rand index (with 5000 resamplings) was 0.86, supporting the robustness of our endotyping result (both values above 0.7 considered as stable). 

Nevertheless, our study has several strengths. First, a relatively large sample size of plasma metabolomic datasets analyzed using unsupervised ML techniques distinguished endotypes of OAG associated with clinical outcomes. Second, confirmation of all OAG cases was based on comprehensive ocular examinations, and longitudinal VF follow-up was performed stringently. 

In conclusion, we revealed distinct etiologic endotypes of OAG based on plasma metabolomic data. Pathway analysis of metabolomics implicated dysregulation in lipid and amino acid metabolism, providing clues for elucidating the pathophysiology of OAG and potential endotype-specific strategies for glaucoma prevention and treatment. 

Acknowledgments

The authors assume full responsibility for data analyses and interpretation. 

Data Sharing Statements: The source data for the findings of this study are available from the corresponding author upon reasonable request. Summary statistics supporting our findings are provided in the Supplementary Material. 

Supported by the Tohoku Medical Megabank Organization, Tohoku University, shared in the MEXT Project for promoting public utilization of advanced research infrastructure (Program for advanced research equipment platforms; grant number JPMXS0450100322), Japan Agency for Medical Research and Development (grant number 22ak0101110h0004), Japan Science and Technology Agency (grant number JPMJPF2201), and Katsuaki Kitazawa Memorial Glaucoma Research Fund. The funding source had no involvement in study design; in the collection, analysis, and interpretation of data; in the writing of the report; and in the decision to submit the paper for publication. 

Disclosure: A. Hanyuda, None; Y. Raita, None; T. Ninomiya, None; K. Hashimoto, None; N. Takada, None; K. Sato, None; J. Inoue, None; S. Koshiba, None; G. Tamiya, None; A. Narita, None; M. Akiyama, None; K. Omodaka, None; S. Tsuda, None; Y. Yokoyama, None; N. Himori, None; Y. Yamamoto, None; T. Taniguchi, None; K. Negishi, None; T. Nakazawa, None