Firstly, the chemical structure of ALA was obtained from the PubChem database, and its potential protein targets were predicted using six databases, including SwissTargetPrediction, PharmMapper, Comparative Toxicogenomics Database (CTD), TargetNet, HIT, and Similarity Ensemble Approach (SEA). Secondly, the genes related to AMD were retrieved from the GeneCards (https://www.genecards.org/), DisGeNET (https://disgenet.com/), and OMIM (https://www.omim.org/) databases using the keyword “age-related macular degeneration.” Ferroptosis-related genes were obtained from the FerrDb V2 database (Zhou et al., 2023). All gene sets were merged and deduplicated to generate a comprehensive list of ALA-ferroptosis-AMD interactions. Additionally, all predicted targets were standardized using UniProt and mapped to official gene symbols via the UniProt ID mapping tool to ensure data consistency (Kim et al., 2021). The intersection of ALA-associated targets with AMD- and ferroptosis-related genes was visualized using a Venn diagram, and overlapping genes were retained as candidate targets for further analysis.

To explore the interactions between ALA and AMD- or ferroptosis-related genes, we constructed an ALA-ferroptosis-AMD network using Cytoscape (version 3.9.1) (Shannon et al., 2003), with nodes representing ALA, the ferroptosis-related genes, and AMD-related genes, while edges indicated predicted associations.

The protein-protein interaction (PPI) network was established using the STRING database with a minimum interaction score threshold of 0.4 (Szklarczyk et al., 2019). The key genes were identified using the CytoHubba plugin based on topological parameters such as degree centrality (Chin et al., 2014). Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were performed using the clusterProfiler package in R (v4.2.2) (Wu M. et al., 2021). Pathways with a false discovery rate (FDR) < 0.05 were considered statistically significant. Enrichment results were visualized as bubble plots using the ggplot2 package in R.

Transcriptomic data related to human AMD was retrieved from the GEO database (GSE135092), comprising 104 AMD samples and 433 controls. Probe annotation and gene symbol conversion were performed using custom Perl scripts. The raw expression data were then subjected to log2 transformation and quantile normalization. To remove potential batch effects, the ComBat algorithm from the sva package was applied, yielding a standardized expression matrix for downstream analysis. The differentially expressed genes (DEGs) between AMD and control samples were identified using the limma package (Ritchie et al., 2015). Genes with a fold change greater than 2 (|log₂FC| > 1) and an adjusted p-value <0.05 were considered statistically significant. Overlapping genes between the DEGs and previously identified ALA-ferroptosis-AMD targets were defined as significant differential core genes (SDECGs) and visualized using boxplots and heatmaps. The chromosomal localizations of these SDECGs were extracted using Perl scripts, and circular diagrams were generated using the RCircos package. Pearson correlation analysis (r > 0.6) was performed to evaluate the co-expression relationships among SDECGs, and gene-gene interaction networks were constructed accordingly.

CIBERSORT and single-sample gene set enrichment analysis (ssGSEA) were performed on the transcriptomic dataset to assess the immune landscape associated with AMD and its potential modulation by ALA-related targets (Chen et al., 2018; Subramanian et al., 2005). Specifically, CIBERSORT was run 1,000 times on the samples to estimate the relative abundance of 22 immune cell types in each sample. The ssGSEA scores were calculated with predefined immune-related gene sets, representing the enrichment level of immune cell signatures or pathways across individual samples. Subsequently, differential analysis of immune cell infiltration between the AMD and control groups was performed using the Wilcoxon rank-sum test. Correlation analysis was performed between the previously identified SDECGs and ssGSEA scores to identify immune-gene associations. The correlation results were visualized using the ComplexHeatmap package.

Two machine learning algorithms, including Least Absolute Shrinkage and Selection Operator (LASSO) regression and Random Forest (RF), were used to identify key ferroptosis-related genes. LASSO regression was applied with 10-fold cross-validation to identify the optimal penalty parameter (λ = 0.01), resulting in a subset of genes with non-zero coefficients. In parallel, the RF algorithm constructed 500 decision trees and ranked gene features based on their Gini importance, selecting those with scores exceeding 9.0. The intersection of both methods was considered high-confidence candidate genes for further validation. A nomogram model was subsequently constructed to visualize and quantify their diagnostic contribution.

The three-dimensional structures of ALA and the key target proteins were retrieved from the PubChem and Protein Data Bank (PDB) databases, respectively. Molecular docking was performed using AutoDock Vina (v1.1.2) to predict binding conformations and calculate binding affinities. Heatmaps of docking scores were generated, and the docking poses were visualized in three dimensions using PyMOL. Subsequently, molecular dynamics (MD) simulations were performed using GROMACS 2022 with the Amber99SB-ILDN force field to assess the stability of the ALA-protein complexes. Each complex was solvated in a cubic TIP3P water box (10 × 10 × 10 nm³) with a minimum margin of 1.2 nm from the protein surface. Counterions were added to neutralize the system. Long-range electrostatic interactions were treated using the particle-mesh Ewald (PME) method. The cutoff distances for Coulombic and van der Waals interactions were set to 1.0 nm. Energy minimization was performed using the steepest descent algorithm for 50,000 steps. Equilibration was conducted under the NVT ensemble at 300 K and the NPT ensemble at 1 bar, using the Langevin thermostat and Berendsen barostat, respectively. A 100-nanosecond production run was then performed under periodic boundary conditions.

MD analyses included root mean square deviation (RMSD) to assess global conformational stability, radius of gyration (Rg) to reflect protein compactness, root mean square fluctuation (RMSF) to evaluate residue flexibility, and solvent-accessible surface area (SASA) to measure solvent exposure. The binding free energy was estimated using the gmx MMPBSA tool, with more negative values indicating stronger affinity. Free energy landscape (FEL) analysis was performed using RMSD and Rg as collective variables, with blue regions indicating the most stable conformational states of the system.

Male C57BL/6J mice (6–8 weeks old, 20–25 g) were obtained from Beijing Vital River Laboratory Animal Technology Co., Ltd. (Beijing, China) and housed in a pathogen-free facility under standard conditions (12-h light/dark cycle, 22°C ± 2°C, 50%–60% humidity) with free access to food and water. All animal experiments were approved by the Institutional Animal Care and Use Committee of Wuhan University and Hubei Medical University, and performed in accordance with the ARVO Statement for the Use of Animals in Ophthalmic and Vision Research.

Mice were randomly assigned to three groups: control, sodium iodate (SI), and SI + ALA groups (n = 6 per group). Mice in the SI + ALA group received intraperitoneal injections of ALA (10 mg/kg/day; T5625, >99% purity, Sigma-Aldrich, St. Louis, MO, United States) dissolved in sterile PBS for 7 consecutive days. On day 8, SI (25 mg/kg; S4007, Sigma Industries) dissolved in sterile PBS was administered via tail vein injection. Mice in the SI group received PBS for 7 days followed by SI injection, while the control mice received PBS for the first 7 days and on day 8. After SI administration, mice in the SI + ALA group continued to receive ALA (10 mg/kg/day) via intraperitoneal injection for an additional 15 days, whereas mice in the SI and control groups received equivalent volumes of saline.

At the end of the treatment, animals were euthanized, and the eyes were enucleated for downstream analyses. For histological analysis, eyes were fixed in 4% paraformaldehyde (PFA) for 24 h, dehydrated through a graded ethanol series, cleared in xylene, and embedded in paraffin. For cryo-sectioning, retinas were fixed in 4% PFA for 30 min, dehydrated in a sucrose gradient (10%–30%), embedded in optimal cutting temperature (OCT) compound, and cryo-sectioned at 10 μm using a Leica cryostat.

Following overnight dark adaptation (12 h), electroretinography (ERG) was performed using the Reti-MINER IV system on day 5, day 10, and day 15 after SI injection (Wu T. et al., 2021). Mice were anesthetized and maintained at 37°C. A reference electrode was inserted subcutaneously into the forehead, a ground electrode was placed at the base of the tail, and a recording electrode was positioned on the cornea using a contact lens. Scotopic ERG responses were recorded at increasing flash intensities ranging from 0.01 to 10 cd·s/m². Visual acuity was assessed using the optomotor response (OMR) test on day 15 after SI injection (Wu, et al., 2021b). Mice were placed on a platform surrounded by a virtual cylinder of rotating vertical gratings (100% contrast, 12°/s). Head-tracking responses to the rotating gratings were recorded across spatial frequencies from 0.042 to 0.336 cycles/degree.

Paraffin-embedded retinal sections (5–8 μm) were stained with hematoxylin and eosin (H&E) for morphological evaluation. Iron accumulation was assessed using Prussian blue staining with Perl’s reagent (Solution A and B, mixed 1:1) (Wuhan Service Bio-Technology Co. Ltd.) for 1 h at room temperature, followed by nuclear counterstaining. After dehydration, the sections were observed under a light microscope (BX53, Olympus, Japan). Iron-positive areas were quantified using ImageJ software with red particle enhancement as previously described (Zhang et al., 2018). Reactive oxygen species (ROS) levels were evaluated in frozen retinal sections using dihydroethidium (DHE) staining. Sections were incubated with 1 μmol/L DHE at 37°C for 30 min in the dark. Fluorescence intensity was measured under a fluorescence microscope (BX53, Olympus, Japan) and used as an indicator of ROS accumulation.

Retina and serum were obtained from each group of mice, weighed, and centrifuged. Malondialdehyde (MDA) and GSH levels were detected using commercial kits (Beyotime Institue of Biotechnology, Jiangsu, China) following the manufacturer’s instructions. The absorbance at 532 nm and 412 nm (for MDA and GSH, respectively) were measured by a microplate reader (Multiskan GO; Thermo Scientific, MA, United States).

Total RNA was extracted from retinal tissues using TRIzol reagent (Invitrogen, United States) according to the manufacturer’s instructions. Complementary DNA (cDNA) was synthesized using PrimeScript™ RT Master Mix (Takara, Japan) following the recommended protocol. Quantitative real-time PCR (qRT-PCR) was performed in triplicate using SYBR Green Master Mix on a QuantStudio™ 5 Real-Time PCR System (Applied Biosystems, United States). The primer sequences are listed in Supplementary Table S1. GAPDH was used as the internal control, and relative gene expression levels were calculated using the 2^−ΔΔCT method. All qPCR experiments were performed in three biologically independent replicates, and four retinal sections were analyzed per replicate.

All statistical analyses were performed using GraphPad Prism 9.0 (GraphPad Software, United States) and R software (v4.2.2). Data were expressed as mean ± standard deviation (SD) unless otherwise specified. Comparisons between two groups were conducted using unpaired Student’s t-test. One-way analysis of variance (ANOVA) followed by Dunnet’s post hoc test was used for comparisons among multiple groups. The Benjamini-Hochberg false discovery rate (FDR) correction method was applied to adjust p-values for multiple comparisons. A p-value <0.05 was considered statistically significant. All experiments were performed with at least three biological replicates unless otherwise stated.

Pharmacophore-based target prediction identified 667 potential targets of ALA across six databases. A total of 8,548 AMD-associated genes and 484 ferroptosis-related genes were collected from GeneCards, CTD, OMIM, and FerrDb V2, respectively. A three-way intersection analysis yielded 71 overlapping targets shared among ALA, AMD, and ferroptosis (Figure 1A). With these genes, an ALA-ferroptosis-AMD interaction network was then constructed in Cytoscape, comprising 74 nodes and 213 edges, with an average node degree of 5.757, network heterogeneity of 2.330, density of 0.079, and centrality of 0.919 (Figure 1B). Furthermore, The STRING-derived PPI network (815 edges) revealed a key iron homeostasis module (HMOX1-FTL-SLC40A1) and a redox signaling hub (MAPK1-NOS2-JUN), as well as the hub genes mediating the crosstalk between GSH metabolism and mitochondrial function (Figure 1C; Table 1).

GO and KEGG functional enrichment analyses revealed a total of 2,667 enriched terms. Notably, top terms in the enriched biological processes included cellular response to metal ion, cellular response to reactive oxygen species and oxidative stress (Figure 1D). Molecular functions featured superoxide-generating NAD(P)H oxidase activity and MAP kinase activity (Figure 1E). Cellular component terms were enriched for NADPH oxidase complex and chromatin silencing complex (Figure 1F). KEGG pathway analysis highlighted the AGE-RAGE signaling pathway, the reactive oxygen species signaling pathway, the lipid and atherosclerosis pathway (Figure 1G).

Using the AMD transcriptomic dataset GSE135092 from the GEO database, which comprises 433 control and 104 AMD samples (Supplementary Table S2), the expression profiles of the 71 previously identified intersection genes (Supplementary Table S3) were examined. Among them, 18 core differentially expressed genes (DEGs) were identified (Figure 2A): FABP4, BCAT2, AKR1C3, AHCY, DHODH, NQO1, HMOX1, G6PD, MAPK3, KEAP1, SIRT2 were significantly downregulated in AMD compared to the control, whereas MAPK8, GSK3B, MAPK1, NOS2, SIRT1, MFN2, SREBF1 were upregulated (p < 0.05), suggesting that inflammation and lipid metabolism disorders may play critical roles in AMD development. Chromosomal mapping revealed that two key ferroptosis regulators, KEAP1 and HMOX1, are located at Chr19p13 and Chr22q12, respectively (Figure 2B). Hierarchical clustering and co-expression network analysis showed high inter-gene correlation, indicating their potential cooperative roles in AMD progression (Figures 2C–E). Immunoinformatic analysis using CIBERSORT and ssGSEA revealed significant reductions in CD8⁺ T cells, regulatory T cells (Tregs), activated mast cells, and dendritic cells in the Control and AMD samples (p < 0.05; Figures 2F,G, Supplementary Data S1). Notably, HMOX1 expression was positively correlated with M1 macrophage polarization (r = 0.596), and MAPK8 was strongly associated with the recruitment of follicular helper T cells (r = 0.77, p < 0.001; Figure 2H). These results collectively suggest that ferroptosis-related pathways, as well as immune dysregulation, may be involved in AMD pathogenesis.

To further identify ferroptosis-related genes with high diagnostic and therapeutic potential, we employed two mainstream machine learning algorithms, LASSO regression and RF, to the 18 core differentially expressed genes (DEGs) identified in prior multi-omics analysis. Specifically, LASSO regression was used for feature compression and variable selection, yielding 16 genes with non-zero coefficients, indicating their substantial contribution to AMD risk prediction (Figure 3A; Table 2). Meanwhile, the RF classification model identified six key genes with importance values >9.0, including AHCY, DHODH, MAPK1, MAPK8, NOS2, and HMOX1, significantly exceeding those of other candidate genes (Figure 3B; Table 3). The intersection between the LASSO and RF-identified targets confirmed these six genes as overlapping top-ranked features (Figure 3C). Subsequently, a multivariate nomogram model was constructed by assigning weights to each gene’s expression level to compute sample-specific risk scores (Figure 3D), showing strong discriminatory power in the GSE135092 dataset. These six genes were subsequently prioritized for downstream molecular docking and in vivo validation, serving as candidate targets for exploring the therapeutic mechanism of ALA in AMD.

We then performed molecular docking experiments to evaluate the binding affinity and interaction patterns of the ALA molecule to the six core molecular targets (Figure 4A). ALA showed favorable binding affinity to all six proteins at physiological temperature (37°C). Key interaction features included: hydrogen bonding between ALA and LYS-54 in MAPK1 (PDB: 4QTA); π-π stacking with His-353 in AHCY (PDB: 1LI4); binding to the ubiquinone site in DHODH (PDB: 4IGH); hydrophobic interactions with LEU-168 in MAPK8 (PDB: 4QTD); localization near the heme-binding pocket of HMOX1 (PDB: 1N45); and polar contacts involving ARG-68, TYR-120, and VAL-71 in NOS2 (PDB: 5XN3). The binding free energies ranged from −5.26 to −6.17 kcal/mol, with dissociation constants (Kd) between 28.9 and 148.3 μM (Figure 4B; Table 4).

Given the central role of HMOX1 in ferroptosis and the strong binding affinity observed for MAPK1, these two complexes were selected for 100-nanosecond MD simulations to assess the stability and binding behavior in detail. For HMOX1, the root mean square deviation (RMSD) stabilized between 0.21 and 0.30 nm after 45 ns (Figure 4C), while the MAPK1 complex reached stability between 0.31 and 0.36 nm after 40 ns (Figure 4D). Radius of gyration (Rg) values for HMOX1 and MAPK1 remained stable within 2.25–2.50 nm and 2.00–2.50 nm, respectively, suggesting compact and stable protein structures during simulation. Root mean square fluctuation (RMSF) profiles showed similar trends, with fluctuations remaining below 0.7 nm in both complexes, suggesting high protein rigidity and dynamic stability. Hydrogen bond analysis revealed 0–3 stable H-bonds in the HMOX1 complex, whereas the MAPK1 system exhibited more variability, reflecting greater conformational dynamics. Solvent-accessible surface area (SASA) remained relatively stable throughout the simulation process for both complexes. Finally, free energy landscape (FEL) analysis showed energy minima located in low RMSD and low Rg regions, suggesting that both complexes adopt stable conformations during the simulation.

Moreover, we evaluated and decomposed the binding free energy of the two complexes using MM-PBSA analysis (Supplementary Figure S1). For HMOX1, the gas-phase binding energy (GGAS) was −40.97 kcal/mol, with favorable van der Waals (VDWAALS <0) and electrostatic (EEL <0) contributions (Supplementary Figure S1A). The solvation energy (GSOLV) was 20.43 kcal/mol, including a small favorable non-polar term (ESURF <0) and an unfavorable polar solvation term (EGB >0). The total binding free energy, i.e., the sum of GGAS and GSOLV, was −20.54 ± 5.01 kcal/mol. The key contributing residues included HIS25, ALA28, MET34, LEU138, GLY139, and PHE207 in HMOX1 (Supplementary Figure S1B). The MAPK1-ALA complex showed a similar energy contribution pattern, with a favorable GGAS and a partially unfavorable GSOLV, leading to a total a binding free energy of −17.30 ± 6.33 kcal/mol (Supplementary Figure S1C), with dominant interactions involving residues GLU71, LEU115, LEU121, HIS125, LEU157, THR159, and CYS161 (Supplementary Figure S1D).

To validate the effects of ALA in vivo, we established a SI-induced retinal injury mouse model and evaluated the electrophysiological and visual functions of the mice (Figure 5A). The representative full-field ERG waves of the control, SI, and SI + ALA groups are shown in Figure 5B. Quantitative analyses further showed that, compared to the control group, the a- and b-wave amplitudes were significantly reduced in the SI group in a time-dependent manner, with the most significant deficits observed on day 15 (a-wave: 7.2% of control; b-wave: 6.4% of control; p < 0.001). In contrast, ALA treatment exhibited a significant protective effect (a-wave: 46.5% of control on day 15; b-wave: 39.8% of control on day 15; p < 0.001), representing 6.5-fold and 6.2-fold improvements over the SI group, respectively (Figures 5C,D). Furthermore, the OMR assay confirmed partial preservation of visual acuity in the SI + ALA group (Figure 5E). Visual thresholds in the SI + ALA group remained at 77.1% of the control levels, compared to only 19.4% in the SI group (p < 0.001; Figure 5F), representing a ∼4.0-fold improvement. These results clearly demonstrate that ALA treatment preserves retinal electrophysiological function and visual behavior in vivo.

Having confirmed the protective effects of ALA in SI-induced AMD mice, we further investigated its impact at the tissue and cellular levels. Histopathological analysis revealed that retinas from SI-treated mice were significantly thinner, with disorganized cellular architecture and abnormal deposit accumulation compared to the control group (Figures 6A,B). In contrast, ALA administration reduced the number of retinal deposits and significantly improved retinal thickness compared to the SI group (p < 0.001). To account for potential distortion from tissue edema or intercellular fluid, we quantified the number of nuclei in the outer nuclear layer (ONL) across defined intervals (200–1,600 μm from the optic nerve head). ALA treatment effectively rescued SI-induced ONL thinning (Figure 6C).

Prussian blue staining demonstrated substantially enhanced iron accumulation in SI-treated retinas compared to the controls. Significantly more iron-positive cells (blue arrows) could be observed in the inner retina in the SI group than in the control group (p = 0.006, Figures 6D,E). Meanwhile, due to SI-induced pigment epithelial atrophy, more depigmented RPE cells were observed in the SI group than in the SI + ALA group (black arrows, Figure 6D). In contrast, ALA treatment significantly reduced SI-induced iron deposition in the inner retina, with a 74.9% decrease in iron-positive cells compared to the SI group (p = 0.029; Figures 6D,E). Contrast-enhanced images further confirmed greater iron levels in the SI group compared to the control group, which was markedly reduced following ALA treatment (Supplementary Figure S2). Furthermore, DHE fluorescence analysis demonstrated significantly elevated ROS levels in SI-treated retinas, which were reduced by 22.2% following ALA treatment (p < 0.001; Figures 6F,G), suggesting substantial mitigation of oxidative stress. Additionally, we measured the retinal and serum levels of malondialdehyde (MDA) and GSH of the mice (Supplementary Figure S3). The retinal and serum MDA levels were significantly elevated in the SI-treated group, while GSH levels were decreased; notably, ALA treatment significantly reduced MDA levels and restored GSH levels, suggesting that ALA treatment could effectively relieve the oxidative stress in SI-treated mice. These pieces of evidence support the role of ALA as a potent iron chelator capable of reversing SI-induced iron deposition and pigmentary disruption.

Subsequently, we explored the expression levels of the six core genes identified previously. Using qRT-PCR, we found that Nos2 and Mapk8 were significantly upregulated following SI treatment, along with moderate but significant increases in Hmox1 and Mapk1 (p < 0.05, Figure 6H). ALA treatment further enhanced the expression level of Hmox1, while suppressing those of Nos2, Mapk8, and Mapk1. On the contrary, following SI treatment, the expression levels of metabolic regulators, Ahcy and Dhodh, were significantly downregulated by 55% and 48%, respectively (p < 0.05, Figure 6H), both of which were restored by ALA treatment (p < 0.05). Together, these findings suggest that ALA may alleviate SI-induced retinal degeneration by reducing iron accumulation, mitigating oxidative stress, and modulating ferroptosis-related gene expression.