Patient-derived cornea organoid model to study metabolomic characterization of rare disease: aniridia-associated keratopathy
Aniridia (OMIM #106,210) is defined as a rare pan-ocular disease linked to mutations in the paired-box gene 6 ( PAX6 ), typically associated with complete or partial iris hypoplasia and foveal hypoplasia . This rare condition is associated with haploinsufficiency in the PAX6, caused by heterozygous nonsense mutations in over 90% of the gene copies . PAX6 is a crucial gene involved in ectodermal differentiation and eye development, and it is highly conserved across many species. Haploinsufficiency in PAX6 causes reduced levels of the PAX6 protein, leading to abnormal eye development . Aniridia, a disease marked by high penetrance and variable phenotype, even within the same family, is inherited in an autosomal dominant manner . Studies have reported significant phenotypic variability among individuals with the same mutation, highlighting the complexity of the disease . Understanding the detailed molecular and pathophysiological mechanisms of this condition is crucial for advancing treatment options. Generating Drosophila , zebrafish, quail, Xenopus , and rodents have been used to model aniridia by Pax6 deletions after the discovery of PAX6 homology across different species, as well as the pathways affected by PAX6 haploinsufficiency. Mouse models for aniridia harbors Pax6 deletions showed the dosage sensitivity of the Pax6 gene has been demonstrated and consequently, therapies that target gene dosage of Pax6 have been explored . However, patient-derived limbal epithelial stem cells (LESCs) were also extracted from aniridia patients as in vitro models with reduced PAX6 protein levels . Patient-derived induced human induced pluripotent stem cells (hiPSCs) are valuable models for studying ophthalmologic diseases, modeling eye development, and drug treatments. Only a few hiPSC lines with PAX6 mutations have been generated . 3D optic cups and 2D limbal epithelial stem cell (LESC) models were generated from aniridia patients with PAX6 nonsense mutation . More complex models of human eye tissues can recapitulate eye needs to be generated from aniridia patient-derived hiPSCs as 3D organoid models. Organoids are self-organized into sophisticated, multicellular assemblies that closely resemble their in vivo counterparts, making them invaluable for investigating human development, disease pathogenesis, and therapeutic responses . Recent advancements have enabled the modeling of various diseases such as cancer, neurodegenerative disorders, and infectious diseases, providing crucial insights into their underlying mechanisms and potential treatment strategies  (Yin et al., 2021). Several studies indicated the organoid models from hiPSCs for ocular development and disease models, including cornea, retina, and lacrimal gland, . Recent studies have highlighted the potential of cornea organoids as models for studying corneal diseases and evaluating therapeutic strategies . Cornea organoids, derived from iPSCs, closely mimic the structure and function of the human cornea, making them ideal for such investigations . However, metabolomic analyses using techniques like gas chromatography-mass spectrometry (GC–MS) have been instrumental in these studies, providing insights into the metabolic changes associated with corneal diseases and their treatments as well as promising biomarkers for clinical applications . In this study, we aimed to generate cornea organoids from AAK patients-derived hiPSCs (AAK-hiPSCs), to characterize their histologic and metabolomic profile compared to healthy control. 

Metabolomic profiles of undifferentiated aniridia patient-derived and healthy control-derived hiPSCs Patient-derived hiPSCs are valuable models for studying the pathophysiology underlying AAK in humans. In this study, we used three aniridia patient-derived hiPSCs lines with three different PAX6 heterozygous nonsense mutations and one wild-type (WT) hiPSC. hiPSCs were differentiated into cornea organoids and characterized their histologic and metabolomic profiles (Fig. ). WT and patients-derived hiPSCs showed typical pluripotent stem cell morphologies with defined borders for pluripotent stem cells (Fig. A). Their pluripotency was also confirmed through pluripotency markers expression by OCT3/4 and NANOG immunofluorescence staining (Fig. B) and SOX2 and KLF4 mRNA levels by qPCR (Fig. C). Overall, the morphological characteristics, immunofluorescence staining, and qPCR analysis collectively demonstrate that the hiPSCs utilized in this study exhibit robust pluripotency. This confirms their suitability for modeling AAK and furthering our understanding of the disease mechanisms in a human context. We next evaluated the metabolic profiles of four distinct hiPSC lines to further evaluate insight into the metabolic flux. Untargeted metabolomic analysis was conducted on undifferentiated cell samples using GC–MS. A total of 267 mass features were identified across all samples, with 90 of these being confirmed using reference databases (see Supplementary Excel file). The metabolomic profiles were subjected to multivariate statistical analysis, including principal component analysis (PCA), partial least squares discriminant analysis (PLS-DA), and a heatmap (Fig. A-D). The two-dimensional and three-dimensional PCA results in the absence of systemic errors or outliers in the dataset and demonstrated a clear contrast between the WT and AAK groups. The WT group (light blue) is observed to cluster separately from the AAK groups, indicating the presence of distinct metabolic profiles. Among the AAK groups, AAK-3 (green) exhibits some degree of overlap with WT, whereas AAK-1 (red) and AAK-2 (blue) are more distinctly separated. (Fig. A, B). The divergence between the WT and AAK groups is more pronounced in comparison to the PCA plot, which serves to illustrate the capacity of PLS-DA to differentiate the groups based on their metabolic profiles. A distinct separation between the metabolic profiles of the WT and AAK samples was evident in the PLS-DA models. The quality of the resulting PLS-DA model was assessed through internal cross-validation. The R2 and Q2 values were 0.99 and 0.77, respectively, indicating the presence of notable metabolic variations in the metabolome of the cell lines. (Fig. C). The metabolites that exhibited the greatest decrease in concentration between the WT and AAK samples were alanine, trans-3-hydroxy-L-proline, serine, glycine, succinic acid, leucine, 4-acetamidobutyric acid, glycerol 1-phosphate, and DL-glyceraldehyde. The concentrations of glyceraldehyde, cholesterol, threose, pantothenic acid, N-acetyl-L-aspartic acid, lysine, and asparagine are higher in WT samples, while other metabolites are lower significantly. This is in contrast to AAK samples, where the opposite was observed. The results of the hierarchical clustering analysis indicate that the samples derived from the WT cluster together, while those derived from AAK-1, AAK-2, and AAK-3 form separate clusters. (Fig. D). Generation of cornea organoids from aniridia patient-derived and healthy control-derived hiPSCs To generate corneal organoids from hiPSCs, differentiation was initiated when the cells reached 70–80% confluency. For the first two days, the cells were exposed to differentiation medium (DM) to promote spontaneous differentiation. Following this, the cells were cultured in retinal differentiation medium (RDM) for differentiation into ocular cells. After 30 days of differentiation, the adherent cells were gently detached with a cell scraper and transferred to an ultra-low attachment six-well plate suspension culture. Over the next 45 days, the cells were maintained in corneal differentiation medium (CDM) to support the development and maturation of corneal-like transparent (Fig. A). Morphologically, the cornea organoids appeared as spherical structures with translucent or transparent layers (Fig. B), suggesting they are fluid-filled. Histological analysis using hematoxylin and eosin staining of cornea organoid sections displayed a stratified epithelium and stroma-like lumen within the semi-transparent or transparent-like regions (Fig. C), indicating a cellular organization that closely resembles that of the natural cornea. Further characterization of the cornea organoids was performed using immunofluorescence staining for Collagen IV, a corneal stroma protein , and CK3, A corneal epithelial marker CK3 (Fig. D). Corneal organoids developed from hiPSCs showed significant positive expression of CK12, a biomarker of mature corneal epithelium (Fig. E) . This finding suggests the successful differentiation of corneal epithelial cells and that our organoid model reflects the corneal epithelial phenotype. Keratocan, a proteoglycan specific to stromal cells , and CK5, which is expressed in the basal layer of the corneal epithelium , were immunofluorescence-stained in corneal organoids (Fig. F). In our organoid model, CK5 expression facilitates precise epithelial cell feature mimicry. The possibility that organoids could mimic the corneal stroma is suggested by the expression of keratocan in their stromal-like layers. Pax6, the master regulator of eye development , is also expressed in corneal epithelial cells. p63 also regulates eye development and epithelial stem cell maintenance and is highly expressed in limbal stem/progenitor cells . Immunofluorescence staining for Pax6 and p63 showed that Pax6 expression was observed in corneal organoids developed from healthy (WT) hiPSCs, while no Pax6 expression was observed in organoids developed from Pax6 mutant hiPSCs. p63 expression was successfully expressed in the developed corneal organoids (Fig. G). To demonstrate the continued proliferation and development in corneal organoids, immunofluorescence staining of the cell proliferation marker Ki67 was performed (Fig. I) Nuclear localization of Ki67 in corneal organoids developed from healthy (WT) and Pax6 mutant hiPSC sources showed positive results. The combined positive results of limbal stem cell-associated protein p63 and Ki67 also support that the basic features of the corneal epithelium are mimicked. Tuj1 (β-III tubulin) positive expression was demonstrated in the developed corneal organoids (Fig. H). Tuj1 supports the presence of neuroectodermal cells that simulate the development of corneal innervation. Tuj1 is known as a marker of mature and developing neuronal cells and is critical for the demonstration of a phenotype that supports healthy innervation processes in corneal tissue . The distribution of Tuj1 positive cells by Tuj1 immunofluorescence staining supports the model representing the early stages of corneal innervation. These findings demonstrate that corneal organoids from both patients-derived and healthy control-derived hiPSCs can successfully mimic the structural organization of the natural cornea. Metabolomic profiles of cornea organoids from aniridia patient-derived and healthy control-derived hiPSCs Following the initiation of the differentiation process at days 30, 44 and 60, the metabolomic profiles of organoids are elucidated through the application of untargeted GC–MS metabolomic analysis. A total of 540 mass features were identified across all samples, with 134 of these being confirmed using reference databases (see Supplementary Excel file). A multivariate statistical analysis was conducted on the metabolomic profiles using partial least squares discriminant analysis (PLS-DA) heatmaps, as illustrated in Fig. A-F. The metabolomic profiles of AAK-1 and AAK-2 organoids were similar, while those of AAK-3 and WT were distinct at day 30 (Fig. A). At day 44, the groups were combined, and a similar metabolomic profile was observed between them, with the exception of AAK-1 (Fig. B). Finally, at day 60, AAK-1 and AAK-2 organoids exhibited a metabolomic profile that clustered together, while AAK-3 and WT organoids displayed a distinct profile (Fig. C). The heatmap for day 30 illustrates the comparative concentrations of diverse metabolites across the various experimental groups. Metabolites such as D-Ala/D-Ala, L-proline, and capric acid displayed distinctive patterns between the WT and AAK groups. The results of the hierarchical clustering analysis indicated that the samples from the WT group clustered together, while the samples from the AAK-1, AAK-2, and AAK-3 groups formed separate clusters (Fig. D). The heatmap for day 44 demonstrated a more pronounced differentiation in metabolite concentrations in comparison to day 30. Notably, metabolites such as D-Ala/D-Ala, L-proline, and xylitol exhibit pronounced differences between the WT and AAK groups. The clustering pattern remained consistent, with the WT samples forming a distinct cluster from the AAK groups (Fig. E). The day 60 heatmap demonstrated the existence of discrete metabolic profiles for the WT and AAK groups. Metabolites such as beta-glycerolphosphate, tartaric acid, and capric acid exhibited notable distinctions. The results of hierarchical clustering indicated the presence of discrete clusters for the WT and AAK groups, which is consistent with the findings of previous heatmaps (Fig. F). We further performed the comparison between the groups Metabolite Set Enrichment Analysis (MSEA) in metabolic pathways. AAK-1 and AAK-2 revealed a notable enrichment in pathways including alanine, aspartate, and glutamate metabolism; glyoxylate and dicarboxylate metabolism; glycolysis/gluconeogenesis; and starch and sucrose metabolism. Among these, alanine, aspartate, and glutamate metabolism exhibited the highest degree of enrichment, indicating substantial metabolic differences between AAK-1 and AAK-2 (Fig. A). Similarly, the comparison between AAK-1 and AAK-3 revealed significant enrichment in alanine, aspartate, and glutamate metabolism, along with glyoxylate and dicarboxylate metabolism, glycolysis/gluconeogenesis, and the citrate cycle (TCA cycle). The alanine, aspartate, and glutamate metabolism pathways exhibited the greatest degree of enrichment, indicating that it is subject to significant metabolic alterations between AAK-1 and AAK-3 (Fig. B). The comparison between AAK-2 and AAK-3 also revealed significant enrichment in alanine, aspartate, and glutamate metabolism, along with glyoxylate and dicarboxylate metabolism, glycolysis/gluconeogenesis, and starch and sucrose metabolism. This consistent enrichment highlights the existence of significant metabolic differences between AAK-2 and AAK-3 at D30 (Fig. C). A comparison between AAK-1 and AAK-2 at D44 revealed a notable enrichment in pathways including valine, leucine, and isoleucine biosynthesis; alanine, aspartate, and glutamate metabolism; glyoxylate and dicarboxylate metabolism; and the TCA cycle. The valine, leucine, and isoleucine biosynthesis pathways exhibited a pronounced increase in enrichment, underscoring significant metabolic disruptions between AAK-1 and AAK-2 (Fig. D). A similar pattern of enrichment was observed in the comparison between AAK-1 and AAK-3, with an increase in the valine, leucine, and isoleucine biosynthesis pathway, alanine, aspartate, and glutamate metabolism, glyoxylate and dicarboxylate metabolism, and the TCA cycle. This consistent pattern of enrichment underscored the existence of substantial metabolic alterations between AAK-1 and AAK-3 at D44 (Fig. E). In the comparison between AAK-2 and AAK-3, there was a notable enrichment in pathways including valine, leucine, and isoleucine biosynthesis; alanine, aspartate, and glutamate metabolism; glyoxylate and dicarboxylate metabolism; and the TCA cycle. These findings indicate that there are persistent metabolic differences between AAK-2 and AAK-3 at this time point (Fig. F). At D60, the comparison between AAK-1 and AAK-2 continued to demonstrate enrichment in valine, leucine, and isoleucine biosynthesis; alanine, aspartate, and glutamate metabolism; glyoxylate and dicarboxylate metabolism; and the TCA cycle. This consistent pattern of enrichment suggests the presence of ongoing metabolic differences between AAK-1 and AAK-2 (Fig. G). The comparison between AAK-1 and AAK-3 also revealed persistent enrichment in these pathways, indicating significant metabolic disruptions between these strains at D60 (Fig. H). Finally, the comparison between AAK-2 and AAK-3 demonstrated sustained enrichment in valine, leucine, and isoleucine biosynthesis, alanine, aspartate, and glutamate metabolism, glyoxylate and dicarboxylate metabolism, and the TCA cycle. These findings underscore the persistence of metabolic disparities between AAK-2 and AAK-3 at D60 (Fig. I). Overall, the metabolite enrichment analysis across diverse AAK organoids and time points consistently demonstrated substantial alterations in pivotal metabolic pathways, notably those pertaining to amino acid metabolism and primary energy production. These outcomes offer invaluable insights into the metabolic consequences of different AAK organoids. MSEA revealed the existence of enriched pathways associated with upregulated metabolites, including amino acid biosynthesis (e.g., L-serine, L-methionine, glutamic acid), carbohydrate metabolism (e.g., L-sorbose, cellobiose), nucleotide metabolism (e.g., uracil, thymine), and lipid metabolism (e.g., heptadecanoic acid, methyl linolenate). (Fig. A, C, E) In contrast, the downregulated metabolites were found to be enriched in pathways related to amino acid degradation (e.g., sarcosine, urea), lipid metabolism (e.g., myristic acid, methyl palmitoleate), and energy production (e.g., fumaric acid, pyruvic acid). The results demonstrate extensive metabolic alterations across the AAK mutants, with both shared and distinct metabolomic profiles indicating the specific functions of AAK-1, AAK-2, and AAK-3 in regulating amino acid metabolism, carbohydrate utilization, and energy balance (Fig. B, D, F). The metabolomic analysis revealed significant alterations in both upregulated and downregulated metabolites across the WT and AAK groups. Sixteen metabolites were commonly upregulated across all groups (WT vs. AAK-1, WT vs. AAK-2, and WT vs. AAK-3). These included diethyl oxalpropionate, L-methionine, uracil, aspartic acid, tryptophane, xanthosine, thymine, and phenyl-beta-glucopyranoside. These include ide, phenylalanine, sedoheptulose anhydride monohydrate, alpha-D-glucosamine 1-phosphate, 1-methylnicotinamide, cysteine, gluconic acid, L-pyroglutamic acid, and citric acid, which collectively reflect shared disruptions in amino acid biosynthesis, nucleotide metabolism, and carbohydrate pathways. A comparison of the WT vs. AAK-1 and WT vs. AAK-2 groups revealed that nine metabolites were upregulated. The presence of L-sorbose, 2-keto-L-gulonic acid, trans-4-hydroxy-L-proline, L-glutamic acid, maltose, cellobiose, tyrosine, glucoheptonic acid, and 5-aminovaleric acid indicates a common increase in amino acid and sugar metabolism. A comparison of WT vs. AAK-1 and WT vs. AAK-3 revealed three shared upregulated metabolites: sucrose, malonic acid, and 2-amino-1-phenylethanol. This suggests alterations in energy intermediates. A comparison of WT vs. AAK-2 and WT vs. AAK-3 revealed nine metabolites, including L-glutamine, L-serine, lysine, L-proline, putrescine, cycloleucine, 2-amino-2-methyl-1,3-propanediol, glycine, and glycerol 1-phosphate, that were consistently upregulated. This finding reflects an increase in amino acid biosynthesis and nitrogen metabolism. In the set of unique upregulated metabolites, WT vs. AAK-1 exhibited six: pipecolic acid, 2-ketoisocaproic acid, DL-3-aminoisobutyric acid, galacturonic acid, heptadecanoic acid, and lactamide. The WT vs. AAK-2 comparison exhibited a unique upregulation of isoleucine, L-cysteine, L-leucine, and talose, which suggests the regulation of branched-chain amino acids. In the WT vs. AAK-3 comparison, L-alanine, resveratrol, serine, D-lyxosylamine, L-homoserine, glycolic acid, D-Ala-D-Ala, lactobionic acid, and tartaric acid were uniquely upregulated, indicating the involvement of pathways associated with antioxidants and sugar metabolism (Fig. G). Four metabolites (sarcosine, serine, carbamic acid ethyl ester (urethane), and glutamic acid) exhibited consistent downregulation across WT vs. AAK-1, WT vs. AAK-2, and WT vs. AAK-3, indicating disruptions in amino acid degradation and nitrogen metabolism. In comparison between WT and AAK-1, urea was the sole metabolite that exhibited a unique downregulation, indicating a potential reduction in urea cycle activity. A comparison between WT vs AAK-1 and WT vs AAK-2 revealed that eleven metabolites were shared between the two groups and exhibited a downregulation: 3-methyl-2-oxobutanoic acid, D-lyxose, myristic acid, acetol, L-valine, threose, pyruvic acid, N-ethylglycine, DL-3-aminoisobutyric acid, pipecolic acid, and ribitol. This finding reflects impaired lipid catabolism and sugar processing. The unique downregulated metabolites in WT vs. AAK-1 included glycerol, fumaric acid, glyceric acid, phosphoric acid, acetylisatin, cis-4-hydroxycyclohexanecarboxylic acid, 2-piperidone, allo-inositol, N-methyl-DL-glutamic acid, and D-threitol, indicating disturbances in phosphate and energy metabolism. In the WT vs. AAK-2 comparison, twelve unique downregulated metabolites were identified. The presence of D-mannitol, hypotaurine, methyl linolenate, L-norleucine, DL-isoleucine, methyl palmitoleate, DL-glyceraldehyde, threonine, tyrosine, phenylacetaldehyde, iminodiacetic acid, and D-allose suggests disruptions in lipid and sugar pathways. In WT vs. AAK-3, the unique downregulated metabolites included L-homoserine, 1-methylhydantoin, glycolic acid, putrescine, methyl-beta-D-galactopyranoside, and resveratrol, which further highlight disruptions in polyamine biosynthesis, antioxidant pathways, and sugar metabolism (Fig. H). 

Patient-derived hiPSCs are valuable models for studying the pathophysiology underlying AAK in humans. In this study, we used three aniridia patient-derived hiPSCs lines with three different PAX6 heterozygous nonsense mutations and one wild-type (WT) hiPSC. hiPSCs were differentiated into cornea organoids and characterized their histologic and metabolomic profiles (Fig. ). WT and patients-derived hiPSCs showed typical pluripotent stem cell morphologies with defined borders for pluripotent stem cells (Fig. A). Their pluripotency was also confirmed through pluripotency markers expression by OCT3/4 and NANOG immunofluorescence staining (Fig. B) and SOX2 and KLF4 mRNA levels by qPCR (Fig. C). Overall, the morphological characteristics, immunofluorescence staining, and qPCR analysis collectively demonstrate that the hiPSCs utilized in this study exhibit robust pluripotency. This confirms their suitability for modeling AAK and furthering our understanding of the disease mechanisms in a human context. We next evaluated the metabolic profiles of four distinct hiPSC lines to further evaluate insight into the metabolic flux. Untargeted metabolomic analysis was conducted on undifferentiated cell samples using GC–MS. A total of 267 mass features were identified across all samples, with 90 of these being confirmed using reference databases (see Supplementary Excel file). The metabolomic profiles were subjected to multivariate statistical analysis, including principal component analysis (PCA), partial least squares discriminant analysis (PLS-DA), and a heatmap (Fig. A-D). The two-dimensional and three-dimensional PCA results in the absence of systemic errors or outliers in the dataset and demonstrated a clear contrast between the WT and AAK groups. The WT group (light blue) is observed to cluster separately from the AAK groups, indicating the presence of distinct metabolic profiles. Among the AAK groups, AAK-3 (green) exhibits some degree of overlap with WT, whereas AAK-1 (red) and AAK-2 (blue) are more distinctly separated. (Fig. A, B). The divergence between the WT and AAK groups is more pronounced in comparison to the PCA plot, which serves to illustrate the capacity of PLS-DA to differentiate the groups based on their metabolic profiles. A distinct separation between the metabolic profiles of the WT and AAK samples was evident in the PLS-DA models. The quality of the resulting PLS-DA model was assessed through internal cross-validation. The R2 and Q2 values were 0.99 and 0.77, respectively, indicating the presence of notable metabolic variations in the metabolome of the cell lines. (Fig. C). The metabolites that exhibited the greatest decrease in concentration between the WT and AAK samples were alanine, trans-3-hydroxy-L-proline, serine, glycine, succinic acid, leucine, 4-acetamidobutyric acid, glycerol 1-phosphate, and DL-glyceraldehyde. The concentrations of glyceraldehyde, cholesterol, threose, pantothenic acid, N-acetyl-L-aspartic acid, lysine, and asparagine are higher in WT samples, while other metabolites are lower significantly. This is in contrast to AAK samples, where the opposite was observed. The results of the hierarchical clustering analysis indicate that the samples derived from the WT cluster together, while those derived from AAK-1, AAK-2, and AAK-3 form separate clusters. (Fig. D). 

To generate corneal organoids from hiPSCs, differentiation was initiated when the cells reached 70–80% confluency. For the first two days, the cells were exposed to differentiation medium (DM) to promote spontaneous differentiation. Following this, the cells were cultured in retinal differentiation medium (RDM) for differentiation into ocular cells. After 30 days of differentiation, the adherent cells were gently detached with a cell scraper and transferred to an ultra-low attachment six-well plate suspension culture. Over the next 45 days, the cells were maintained in corneal differentiation medium (CDM) to support the development and maturation of corneal-like transparent (Fig. A). Morphologically, the cornea organoids appeared as spherical structures with translucent or transparent layers (Fig. B), suggesting they are fluid-filled. Histological analysis using hematoxylin and eosin staining of cornea organoid sections displayed a stratified epithelium and stroma-like lumen within the semi-transparent or transparent-like regions (Fig. C), indicating a cellular organization that closely resembles that of the natural cornea. Further characterization of the cornea organoids was performed using immunofluorescence staining for Collagen IV, a corneal stroma protein , and CK3, A corneal epithelial marker CK3 (Fig. D). Corneal organoids developed from hiPSCs showed significant positive expression of CK12, a biomarker of mature corneal epithelium (Fig. E) . This finding suggests the successful differentiation of corneal epithelial cells and that our organoid model reflects the corneal epithelial phenotype. Keratocan, a proteoglycan specific to stromal cells , and CK5, which is expressed in the basal layer of the corneal epithelium , were immunofluorescence-stained in corneal organoids (Fig. F). In our organoid model, CK5 expression facilitates precise epithelial cell feature mimicry. The possibility that organoids could mimic the corneal stroma is suggested by the expression of keratocan in their stromal-like layers. Pax6, the master regulator of eye development , is also expressed in corneal epithelial cells. p63 also regulates eye development and epithelial stem cell maintenance and is highly expressed in limbal stem/progenitor cells . Immunofluorescence staining for Pax6 and p63 showed that Pax6 expression was observed in corneal organoids developed from healthy (WT) hiPSCs, while no Pax6 expression was observed in organoids developed from Pax6 mutant hiPSCs. p63 expression was successfully expressed in the developed corneal organoids (Fig. G). To demonstrate the continued proliferation and development in corneal organoids, immunofluorescence staining of the cell proliferation marker Ki67 was performed (Fig. I) Nuclear localization of Ki67 in corneal organoids developed from healthy (WT) and Pax6 mutant hiPSC sources showed positive results. The combined positive results of limbal stem cell-associated protein p63 and Ki67 also support that the basic features of the corneal epithelium are mimicked. Tuj1 (β-III tubulin) positive expression was demonstrated in the developed corneal organoids (Fig. H). Tuj1 supports the presence of neuroectodermal cells that simulate the development of corneal innervation. Tuj1 is known as a marker of mature and developing neuronal cells and is critical for the demonstration of a phenotype that supports healthy innervation processes in corneal tissue . The distribution of Tuj1 positive cells by Tuj1 immunofluorescence staining supports the model representing the early stages of corneal innervation. These findings demonstrate that corneal organoids from both patients-derived and healthy control-derived hiPSCs can successfully mimic the structural organization of the natural cornea. 

Following the initiation of the differentiation process at days 30, 44 and 60, the metabolomic profiles of organoids are elucidated through the application of untargeted GC–MS metabolomic analysis. A total of 540 mass features were identified across all samples, with 134 of these being confirmed using reference databases (see Supplementary Excel file). A multivariate statistical analysis was conducted on the metabolomic profiles using partial least squares discriminant analysis (PLS-DA) heatmaps, as illustrated in Fig. A-F. The metabolomic profiles of AAK-1 and AAK-2 organoids were similar, while those of AAK-3 and WT were distinct at day 30 (Fig. A). At day 44, the groups were combined, and a similar metabolomic profile was observed between them, with the exception of AAK-1 (Fig. B). Finally, at day 60, AAK-1 and AAK-2 organoids exhibited a metabolomic profile that clustered together, while AAK-3 and WT organoids displayed a distinct profile (Fig. C). The heatmap for day 30 illustrates the comparative concentrations of diverse metabolites across the various experimental groups. Metabolites such as D-Ala/D-Ala, L-proline, and capric acid displayed distinctive patterns between the WT and AAK groups. The results of the hierarchical clustering analysis indicated that the samples from the WT group clustered together, while the samples from the AAK-1, AAK-2, and AAK-3 groups formed separate clusters (Fig. D). The heatmap for day 44 demonstrated a more pronounced differentiation in metabolite concentrations in comparison to day 30. Notably, metabolites such as D-Ala/D-Ala, L-proline, and xylitol exhibit pronounced differences between the WT and AAK groups. The clustering pattern remained consistent, with the WT samples forming a distinct cluster from the AAK groups (Fig. E). The day 60 heatmap demonstrated the existence of discrete metabolic profiles for the WT and AAK groups. Metabolites such as beta-glycerolphosphate, tartaric acid, and capric acid exhibited notable distinctions. The results of hierarchical clustering indicated the presence of discrete clusters for the WT and AAK groups, which is consistent with the findings of previous heatmaps (Fig. F). We further performed the comparison between the groups Metabolite Set Enrichment Analysis (MSEA) in metabolic pathways. AAK-1 and AAK-2 revealed a notable enrichment in pathways including alanine, aspartate, and glutamate metabolism; glyoxylate and dicarboxylate metabolism; glycolysis/gluconeogenesis; and starch and sucrose metabolism. Among these, alanine, aspartate, and glutamate metabolism exhibited the highest degree of enrichment, indicating substantial metabolic differences between AAK-1 and AAK-2 (Fig. A). Similarly, the comparison between AAK-1 and AAK-3 revealed significant enrichment in alanine, aspartate, and glutamate metabolism, along with glyoxylate and dicarboxylate metabolism, glycolysis/gluconeogenesis, and the citrate cycle (TCA cycle). The alanine, aspartate, and glutamate metabolism pathways exhibited the greatest degree of enrichment, indicating that it is subject to significant metabolic alterations between AAK-1 and AAK-3 (Fig. B). The comparison between AAK-2 and AAK-3 also revealed significant enrichment in alanine, aspartate, and glutamate metabolism, along with glyoxylate and dicarboxylate metabolism, glycolysis/gluconeogenesis, and starch and sucrose metabolism. This consistent enrichment highlights the existence of significant metabolic differences between AAK-2 and AAK-3 at D30 (Fig. C). A comparison between AAK-1 and AAK-2 at D44 revealed a notable enrichment in pathways including valine, leucine, and isoleucine biosynthesis; alanine, aspartate, and glutamate metabolism; glyoxylate and dicarboxylate metabolism; and the TCA cycle. The valine, leucine, and isoleucine biosynthesis pathways exhibited a pronounced increase in enrichment, underscoring significant metabolic disruptions between AAK-1 and AAK-2 (Fig. D). A similar pattern of enrichment was observed in the comparison between AAK-1 and AAK-3, with an increase in the valine, leucine, and isoleucine biosynthesis pathway, alanine, aspartate, and glutamate metabolism, glyoxylate and dicarboxylate metabolism, and the TCA cycle. This consistent pattern of enrichment underscored the existence of substantial metabolic alterations between AAK-1 and AAK-3 at D44 (Fig. E). In the comparison between AAK-2 and AAK-3, there was a notable enrichment in pathways including valine, leucine, and isoleucine biosynthesis; alanine, aspartate, and glutamate metabolism; glyoxylate and dicarboxylate metabolism; and the TCA cycle. These findings indicate that there are persistent metabolic differences between AAK-2 and AAK-3 at this time point (Fig. F). At D60, the comparison between AAK-1 and AAK-2 continued to demonstrate enrichment in valine, leucine, and isoleucine biosynthesis; alanine, aspartate, and glutamate metabolism; glyoxylate and dicarboxylate metabolism; and the TCA cycle. This consistent pattern of enrichment suggests the presence of ongoing metabolic differences between AAK-1 and AAK-2 (Fig. G). The comparison between AAK-1 and AAK-3 also revealed persistent enrichment in these pathways, indicating significant metabolic disruptions between these strains at D60 (Fig. H). Finally, the comparison between AAK-2 and AAK-3 demonstrated sustained enrichment in valine, leucine, and isoleucine biosynthesis, alanine, aspartate, and glutamate metabolism, glyoxylate and dicarboxylate metabolism, and the TCA cycle. These findings underscore the persistence of metabolic disparities between AAK-2 and AAK-3 at D60 (Fig. I). Overall, the metabolite enrichment analysis across diverse AAK organoids and time points consistently demonstrated substantial alterations in pivotal metabolic pathways, notably those pertaining to amino acid metabolism and primary energy production. These outcomes offer invaluable insights into the metabolic consequences of different AAK organoids. MSEA revealed the existence of enriched pathways associated with upregulated metabolites, including amino acid biosynthesis (e.g., L-serine, L-methionine, glutamic acid), carbohydrate metabolism (e.g., L-sorbose, cellobiose), nucleotide metabolism (e.g., uracil, thymine), and lipid metabolism (e.g., heptadecanoic acid, methyl linolenate). (Fig. A, C, E) In contrast, the downregulated metabolites were found to be enriched in pathways related to amino acid degradation (e.g., sarcosine, urea), lipid metabolism (e.g., myristic acid, methyl palmitoleate), and energy production (e.g., fumaric acid, pyruvic acid). The results demonstrate extensive metabolic alterations across the AAK mutants, with both shared and distinct metabolomic profiles indicating the specific functions of AAK-1, AAK-2, and AAK-3 in regulating amino acid metabolism, carbohydrate utilization, and energy balance (Fig. B, D, F). The metabolomic analysis revealed significant alterations in both upregulated and downregulated metabolites across the WT and AAK groups. Sixteen metabolites were commonly upregulated across all groups (WT vs. AAK-1, WT vs. AAK-2, and WT vs. AAK-3). These included diethyl oxalpropionate, L-methionine, uracil, aspartic acid, tryptophane, xanthosine, thymine, and phenyl-beta-glucopyranoside. These include ide, phenylalanine, sedoheptulose anhydride monohydrate, alpha-D-glucosamine 1-phosphate, 1-methylnicotinamide, cysteine, gluconic acid, L-pyroglutamic acid, and citric acid, which collectively reflect shared disruptions in amino acid biosynthesis, nucleotide metabolism, and carbohydrate pathways. A comparison of the WT vs. AAK-1 and WT vs. AAK-2 groups revealed that nine metabolites were upregulated. The presence of L-sorbose, 2-keto-L-gulonic acid, trans-4-hydroxy-L-proline, L-glutamic acid, maltose, cellobiose, tyrosine, glucoheptonic acid, and 5-aminovaleric acid indicates a common increase in amino acid and sugar metabolism. A comparison of WT vs. AAK-1 and WT vs. AAK-3 revealed three shared upregulated metabolites: sucrose, malonic acid, and 2-amino-1-phenylethanol. This suggests alterations in energy intermediates. A comparison of WT vs. AAK-2 and WT vs. AAK-3 revealed nine metabolites, including L-glutamine, L-serine, lysine, L-proline, putrescine, cycloleucine, 2-amino-2-methyl-1,3-propanediol, glycine, and glycerol 1-phosphate, that were consistently upregulated. This finding reflects an increase in amino acid biosynthesis and nitrogen metabolism. In the set of unique upregulated metabolites, WT vs. AAK-1 exhibited six: pipecolic acid, 2-ketoisocaproic acid, DL-3-aminoisobutyric acid, galacturonic acid, heptadecanoic acid, and lactamide. The WT vs. AAK-2 comparison exhibited a unique upregulation of isoleucine, L-cysteine, L-leucine, and talose, which suggests the regulation of branched-chain amino acids. In the WT vs. AAK-3 comparison, L-alanine, resveratrol, serine, D-lyxosylamine, L-homoserine, glycolic acid, D-Ala-D-Ala, lactobionic acid, and tartaric acid were uniquely upregulated, indicating the involvement of pathways associated with antioxidants and sugar metabolism (Fig. G). Four metabolites (sarcosine, serine, carbamic acid ethyl ester (urethane), and glutamic acid) exhibited consistent downregulation across WT vs. AAK-1, WT vs. AAK-2, and WT vs. AAK-3, indicating disruptions in amino acid degradation and nitrogen metabolism. In comparison between WT and AAK-1, urea was the sole metabolite that exhibited a unique downregulation, indicating a potential reduction in urea cycle activity. A comparison between WT vs AAK-1 and WT vs AAK-2 revealed that eleven metabolites were shared between the two groups and exhibited a downregulation: 3-methyl-2-oxobutanoic acid, D-lyxose, myristic acid, acetol, L-valine, threose, pyruvic acid, N-ethylglycine, DL-3-aminoisobutyric acid, pipecolic acid, and ribitol. This finding reflects impaired lipid catabolism and sugar processing. The unique downregulated metabolites in WT vs. AAK-1 included glycerol, fumaric acid, glyceric acid, phosphoric acid, acetylisatin, cis-4-hydroxycyclohexanecarboxylic acid, 2-piperidone, allo-inositol, N-methyl-DL-glutamic acid, and D-threitol, indicating disturbances in phosphate and energy metabolism. In the WT vs. AAK-2 comparison, twelve unique downregulated metabolites were identified. The presence of D-mannitol, hypotaurine, methyl linolenate, L-norleucine, DL-isoleucine, methyl palmitoleate, DL-glyceraldehyde, threonine, tyrosine, phenylacetaldehyde, iminodiacetic acid, and D-allose suggests disruptions in lipid and sugar pathways. In WT vs. AAK-3, the unique downregulated metabolites included L-homoserine, 1-methylhydantoin, glycolic acid, putrescine, methyl-beta-D-galactopyranoside, and resveratrol, which further highlight disruptions in polyamine biosynthesis, antioxidant pathways, and sugar metabolism (Fig. H). 

In this study, we successfully generated corneal organoids from patient-derived hiPSCs with PAX6 heterozygous nonsense mutations associated with AAK. These organoids closely recapitulated the structure and function of the human cornea and were utilized to investigate the metabolomic profiles of both healthy (WT) and AAK corneal organoids. Our findings highlighted important distinctions between WT and AAK groups, providing new insights into the molecular and metabolic pathways that may underlie AAK pathophysiology. One of the key outcomes of this study was the ability to replicate the structural characteristics of corneal tissue in vitro using corneal organoids derived from patient-specific iPSCs. Histological analysis revealed that these organoids consisted of stratified epithelial layers and stroma-like regions, closely resembling the architecture of the human cornea. The presence of corneal-specific markers such as Collagen IV and CK3 further confirmed the successful differentiation and maturation of these organoids. This demonstrates the potential utility of iPSC-derived corneal organoids in modeling human corneal diseases, as previously reported by studies on retinal and other ocular organoids . The results shown are consistent with the results found by Ilmarien et al. . Corneal epithelial cells and limbal epithelial stem cells may not be the primary cause of the AAK disease phenotype. The metabolites and secretions secreted by these corneal epithelial and limbal epithelial cells may affect cells such as keratocytes, microglia and immune cells in the tissue environment, which may explain the initial defect in the patient, inflammation and neurodegeneration. Another main outcome of this study was the generation of corneal organoids from patient-derived iPSCs that resemble the structural and cellular properties of the human cornea. The successful expression of markers specific to different cells of the cornea, such as Pax6, p63, Keratocan, CK5, CK12, and CK3, provides evidence for the appropriate differentiation and maturation of these organoids. These findings are consistent with previous publications on corneal and multiocular organoids developed from iPSCs . Single-cell transcriptomic analyses reveal the complex cellular diversity of these organoids and the heterogeneity present in the organoids . In these studies, the presence of basal and superficial epithelial cells, stromal keratocytes, endothelial cells, and limbal stem/progenitor cell populations was identified in corneal organoids developed from iPSCs. The findings obtained in our study support the similar dynamics of corneal development as shown in these single-cell transcriptomic analyses. While Pax6 expressions were observed in corneal organoids developed from healthy (WT) hiPSCs, which is consistent with the literature, it is noteworthy that Pax6 expressions were not observed in corneal organoids developed from Pax6 mutant hiPSC sources. Pax6 is a transcription factor known to play a critical role in eye development and is important in corneal epithelial differentiation. Haploinsufficiency in individuals with heterozygous mutations leads to a high decrease or complete loss of Pax6 expression, which is consistent with the findings observed in corneal organoids we developed from patient-derived iPSCs . These organoids, beyond their similarity to native tissues, also provide unique platforms for understanding the mechanisms underlying patient-specific diseases. In a study by Lima Cunha et al. , restoration of functional Pax6 in ocular models developed from iPSCs demonstrates the feasibility of using patient-generated organoids for therapeutic discovery . Corneal organoids developed from healthy (WT) iPSCs have shown successful results for modeling, while organoids developed from Pax6 mutant iPSCs show variable expression indicating the complexity of rare genetic diseases. This highlights the value of organoids developed from patient-specific iPSCs for understanding the molecular basis of genetic diseases and developing targeted therapeutic strategies. Metabolomics studies examine the flow and chemical mobility of metabolites in biological samples, providing a full picture of its composition . Global metabolomics studies can provide insights into biological changes induced by disease or environmental interactions, as well as information on the potential pathways involved, which can aid in the identification of potential biomarkers to better disease diagnosis and therapy . The integration of metabolomics with organoid research provides a comprehensive analysis of metabolic pathways and the identification of biomarkers pertinent to disease and treatment outcomes . Our metabolomic studies consistently demonstrated distinct clustering between the WT and AAK mutant groups, indicating significant metabolic differences. Similarly, the PLS-DA plots across different time points and media further reinforce the notion that AAK mutations induced significant metabolic alterations. These clustering patterns align with prior research indicating that mutations in pivotal metabolic genes can result in pronounced shifts in metabolite profiles, thereby driving distinct clustering patterns in multivariate analyses . The heat maps provided comprehensive insights into the distinctive metabolic alterations between the WT and AAK groups. It is noteworthy that metabolites such as D-Ala/D-Ala, L-proline, capric acid, and xylitol exhibited distinctive patterns between the two groups. For instance, the heatmap for day 30 demonstrated notable discrepancies in D-Ala/D-Ala and L-proline levels, indicating the potential for disruptions in amino acid metabolism. These discrepancies become more pronounced by day 44, thereby underscoring the progressive nature of metabolic alterations over time. By day 60, the heatmap demonstrated evident disparities in metabolites such as beta-glycerolphosphate, tartaric acid, and capric acid, thereby further exemplifying disparate metabolic states between the WT and AAK groups. The application of hierarchical clustering to these metabolites consistently produced separate clusters for WT and AAK, mutations can alter the metabolic regulations . The observed metabolic alterations indicated that AAK mutations have a profound impact on fundamental metabolic pathways. The observed alterations in carnitine and putrescine levels indicate disruptions in lipid metabolism and polyamine biosynthesis, respectively. These pathways are essential for maintaining cellular energy homeostasis and responding to metabolic stress . The distinct clustering of AAK-1 and AAK-2 away from WT, with AAK-3 positioned intermediately, indicates varying degrees of metabolic disruption among the AAK mutants, with AAK-3 having a less severe impact on metabolism compared to AAK-1 and AAK-2. MSEA presents a comprehensive metabolite enrichment analysis, highlighting key metabolic pathways altered across different experimental conditions within three datasets: RDM D30, RDM D44, and CDM D60. Each panel (A-I) compares metabolic shifts between pairs of conditions (AAK-1 vs AAK-2, AAK-1 vs AAK-3, and AAK-2 vs AAK-3), consistently revealing significant enrichment in specific amino acid metabolism pathways. During the differentiation process of organoids from D30 to D60 their energy metabolism shifted from glycolysis to oxidative phosphorylation, as suggested by previous studies . In the RDM D30 dataset "Alanine, aspartate, and glutamate metabolism" is consistently identified as the most significantly enriched pathway, particularly in AAK-1 vs. AAK-2 and AAK-1 vs. AAK-3. This pathway is of great importance for maintaining nitrogen balance, neurotransmitter cycling, and energy production through the TCA cycle. The enrichment of "Glycine, serine, and threonine metabolism," which is essential for nucleotide synthesis and one-carbon metabolism, further underscores the pivotal role these pathways play in maintaining cellular metabolic homeostasis under varying conditions . Furthermore, pathways such as "Cysteine and methionine metabolism," which are crucial for antioxidant defense and glutathione production, are conspicuous in AAK-1 vs AAK-3, indicating adaptive metabolic responses to oxidative stress . The consistent presence of amino acid biosynthesis pathways, such as "Phenylalanine, tyrosine, and tryptophan biosynthesis," in AAK-2 vs. AAK-3 (Fig. C) serves to underscore their fundamental role in cellular differentiation and response to environmental changes . In the RDM D44 dataset, the analysis reveals a notable enrichment of "Valine, leucine, and isoleucine biosynthesis" (BCAA biosynthesis) in all comparisons. BCAAs are essential for protein synthesis, energy production, and metabolic regulation, indicating that these processes are vital for maintaining cellular energy homeostasis and anabolic activity under these experimental conditions . This finding is particularly relevant in light of the growing recognition of the role of BCAA metabolism in the pathogenesis of metabolic diseases and cancer, where it provides support for anabolic growth and survival. The significant enrichment of pathways such as "Phenylalanine metabolism" and "Ubiquinone and another terpenoid-quinone biosynthesis" in AAK-1 vs. AAK-3 indicates shifts in mitochondrial function and oxidative stress responses, which are further evidenced by the enrichment of pathways such as BCAA biosynthesis. This highlights the adaptability of cellular metabolism in response to varying conditions . The CDM D60 dataset lends further support to these findings, with BCAA biosynthesis consistently emerging as a key pathway. The repeated enrichment of amino acid metabolism pathways across all datasets and comparisons underscores their central role in cellular adaptation. This likely reflects a conserved metabolic strategy to optimize anabolic and catabolic processes under different experimental stimuli . Overall, the consistent identification of amino acid metabolism pathways, particularly those involved in BCAA biosynthesis and the metabolism of alanine, aspartate, glutamate, glycine, serine, and threonine, emphasizes their critical role in cellular adaptation and stress responses . 

The findings presented in this study provide compelling evidence that AAK PAX6 mutations induce significant alterations in metabolic pathways, particularly those involved in amino acid metabolism. The distinct clustering patterns observed in PCA and PLS-DA analyses, coupled with the detailed metabolite profiling through heatmaps and MSEA, highlight the profound impact of AAK mutations on cellular metabolism. The consistent enrichment of pathways related to BCAA biosynthesis and the metabolism of alanine, aspartate, and glutamate across different datasets underscores their central role in maintaining cellular homeostasis and adaptation to environmental changes. These findings not only deepen our understanding of the metabolic consequences of AAK mutations but also point to potential therapeutic targets for conditions where these metabolic pathways are dysregulated. As a limitation of our metabolomic studies on cornea organoids, further validation experiments confirming the expected roles of the amino acids retrieved in each group have to be conducted. In addition, future research should focus on elucidating the mechanistic underpinnings of these metabolic shifts and exploring their implications in disease contexts. Furthermore, integrated multi-omics could be employed to elucidate the effect of AAK mutations on the organoids. 

Culture of hiPSCs Patient-derived hiPSC lines (AAK1, AAK2, AAK3) generated by Ilmarien et al., 2024 were used in this study (Ilmarien et al., 2024). WT Control hiPSC line (WT) generated by  was used in this study . All cell lines were maintained under feeder-free conditions on six-well plates coated with hESC-qualified Matrigel (Corning), with mTeSR1 medium (STEMCELL Technologies). Feeder-free hiPSCs were passaged every 5–7 days at a ratio of up to 1:10 using ReleSR (STEMCELL Technologies). Differentiation of hiPSCs to cornea organoids Four hiPSC lines were initially cultured in six-well plates coated with 1% hESC qualified Matrigel (STEMCELL Technologies) and when reached 70%−80% confluency, the growth medium was switched with differentiation medium (DM) containing DMEM/F12 (Thermo Fisher Scientific), 8% knockout serum replacement (KOSR, Thermo Fisher Scientific), 1 × non-essential amino acids (NEAA, Lonza), 1 × Glutamax (Thermo Fisher Scientific) and 1 × Pen-Strep (Thermo Fisher Scientific) for 2 days to induce spontaneous differentiation. Subsequently, the medium was replaced to retinal differentiation medium (RDM) containing DM plus 2% B27 supplement (Thermo Fisher Scientific) to induce eye field specification. After 4 weeks of culture ocular phenotype cells were scraped off using cell scrapers and suspended in low adhesion plates (Corning) in RDM for another 2 weeks. Cells were starting to give rise to 3D, transparent-like corneal structures development in the low adhesion plates in suspension culture. To maturation of these transparent-like corneal structures, the medium was switched with corneal differentiation medium (CDM) containing DM plus 1% N2 MAX supplement (R&D Systems), 1 × ITS (Thermo Fisher Scientific), 5 ng/ml FGF2 and 10 ng/ml EGF (Thermo Fisher Scientific). Cornea organoids were maintained in suspension culture for up to 90 days with daily medium refreshment. RNA isolation and quantitative PCR Undifferentiated hiPSCs total RNA extraction was performed with TRIzol reagent, and samples were kept at −80 °C until further processing. Following the manufacturer’s instructions, RNA was isolated using the Macherey–Nagel Nucleospin RNA isolation kit. The OneScript Plus cDNA Synthesis Kit (Applied Biological Materials) was used to synthesize cDNA. The GoTaq Master Mix (Promega) was used for RT-PCR. On the Roche LightCycler 480 device, qPCR was run for two minutes at 95 °C, then for forty cycles (ten seconds at 95 °C, ten seconds at 60 °C, and thirty seconds at 72 °C). There was also a cooling step (40 °C for 30 s) and a melting step (95 °C for 10 s and 65 °C for 30 s). The gene expression levels were normalized to the endogenous control glyceraldehyde 3-phosphate dehydrogenase (GAPDH). The fold change relative to the control was computed using the 2 − Δ(ΔCt) technique. The primers used to determine mRNA levels are presented in the Supplement Table. Immunofluorescence staining hiPSC-derived cornea organoids were fixed with 4% paraformaldehyde (Sigma-Aldrich) for 20 min at room temperature and embedded in cryomatrix (OCT, Fisher Healthcare). 5 μm thick serial sections were blocked/permeabilized with staining solution (1% Bovine Serum Albumin (BSA; Sigma-Aldrich) (w/v) and 0.3% Triton-X100 (neoFroxx) (v/v) in PBS) for 1 h at room temperature. Cells were then treated with primary antibodies at 4 °C overnight and afterward with secondary antibodies for 1 h at room temperature. Nuclei were stained with 4,6-diamidino-2-phenylindole (0.5 μg/ml) (DAPI; Neofroxx). iPSCs were treated with the staining solution for blocking and permeabilization after fixation and staining proceeded as previously explained. Samples were visualized either with Fluorescence Microscopy (Olympus IX71) or Confocal Microscopy (Zeiss LSM880). Antibodies used are listed in the Supplement Table. Sample preparation for metabolomics analysis Four distinct iPSC lines (AAK-1, AAK-2, AAK-3, and WT) in their undifferentiated state were subjected to washing with a 0.9% NaCl solution (Sigma-Aldrich) and subsequently fixed using a methanol:water (9:1, v/v) solution (Sigma-Aldrich). The samples were then flash-frozen in liquid nitrogen and scraped into 2 ml centrifuge tubes. Following centrifugation at 15,000 rpm for 15 min at + 4 °C, the supernatants were carefully transferred into 2 ml Eppendorf tubes and stored at − 80 °C until further analysis. Additionally, culture media from AAK organoid samples were collected on days 30, 44, and 60, and stored at − 80 °C. Metabolomic analysis by GC–MS Samples stored at − 80 °C were thawed at room temperature and evaporated to complete dryness in a vacuum centrifuge (Savant ISS110 SpeedVac Concentrator, Thermo Scientific). The dried samples were methoxylated with 20 µL of methoxyamine hydrochloride (20 mg/ml in pyridine) by incubation in an oven at 30 °C for 90 min. Subsequently, 50 µL of N-methyl-N-trimethylsilyl trifluoroacetamide containing 1% trimethylchlorosilane (Sigma-Aldrich) was added, and the samples were further incubated at 37 °C for 30 min. The derivatized samples were then transferred into silylated GC–MS vials and analyzed using a GC–MS system (Shimadzu-QP2010 Ultra) equipped with a DB-5MS stationary phase column (30 m + 10 m DuraGuard × 0.25 mm i.d. and 0.25-μm film thickness). The solvent delay was set to 5.90 min. The oven temperature was initially maintained at 60 °C for 1 min, followed by a temperature increase at a gradient of 10 °C/min until reaching 325 °C, which was held for 10 min before cooling. The MSD transfer line temperature was maintained at 290 °C, and the column flow rate was held constant at 1 ml He/min. The mass range was set to 50–650 Da, with a total run time of 37.5 min. Data deconvolution, peak alignment, normalization, and data matrix creation were conducted using MS-DIAL (version 4.8) software. Metabolite identification for GC–MS was performed utilizing a commercially available retention index library (Fiehn Retention Index Library) with a 70% or higher identification cut-off score. The data matrix generated by MS-DIAL was exported to an Excel file. Metabolite traits with more than 50% missing values were excluded from the data matrix. Missing values were imputed with half of the smallest concentration value within the metabolite group. The final data matrix was imported into MetaboAnalyst (version 6.0) software for enrichment analysis and multivariate analyses, including principal component analysis (PCA), partial least squares discriminant analysis (PLS-DA), and heatmap analysis. Statistical analysis All experiments were represented for three biological replicates. Data are represented as the mean ± SEM or mean ± SD. The GraphPad Prism 8.0 (GraphPad Software, San Diego, USA) program was used for statistical analysis. Image processing was performed using ImageJ and Zen 3.1 software. Metabolomics analyses were carried out by Metaboanalyst 6.0. 

Patient-derived hiPSC lines (AAK1, AAK2, AAK3) generated by Ilmarien et al., 2024 were used in this study (Ilmarien et al., 2024). WT Control hiPSC line (WT) generated by  was used in this study . All cell lines were maintained under feeder-free conditions on six-well plates coated with hESC-qualified Matrigel (Corning), with mTeSR1 medium (STEMCELL Technologies). Feeder-free hiPSCs were passaged every 5–7 days at a ratio of up to 1:10 using ReleSR (STEMCELL Technologies). 

Four hiPSC lines were initially cultured in six-well plates coated with 1% hESC qualified Matrigel (STEMCELL Technologies) and when reached 70%−80% confluency, the growth medium was switched with differentiation medium (DM) containing DMEM/F12 (Thermo Fisher Scientific), 8% knockout serum replacement (KOSR, Thermo Fisher Scientific), 1 × non-essential amino acids (NEAA, Lonza), 1 × Glutamax (Thermo Fisher Scientific) and 1 × Pen-Strep (Thermo Fisher Scientific) for 2 days to induce spontaneous differentiation. Subsequently, the medium was replaced to retinal differentiation medium (RDM) containing DM plus 2% B27 supplement (Thermo Fisher Scientific) to induce eye field specification. After 4 weeks of culture ocular phenotype cells were scraped off using cell scrapers and suspended in low adhesion plates (Corning) in RDM for another 2 weeks. Cells were starting to give rise to 3D, transparent-like corneal structures development in the low adhesion plates in suspension culture. To maturation of these transparent-like corneal structures, the medium was switched with corneal differentiation medium (CDM) containing DM plus 1% N2 MAX supplement (R&D Systems), 1 × ITS (Thermo Fisher Scientific), 5 ng/ml FGF2 and 10 ng/ml EGF (Thermo Fisher Scientific). Cornea organoids were maintained in suspension culture for up to 90 days with daily medium refreshment. 

Undifferentiated hiPSCs total RNA extraction was performed with TRIzol reagent, and samples were kept at −80 °C until further processing. Following the manufacturer’s instructions, RNA was isolated using the Macherey–Nagel Nucleospin RNA isolation kit. The OneScript Plus cDNA Synthesis Kit (Applied Biological Materials) was used to synthesize cDNA. The GoTaq Master Mix (Promega) was used for RT-PCR. On the Roche LightCycler 480 device, qPCR was run for two minutes at 95 °C, then for forty cycles (ten seconds at 95 °C, ten seconds at 60 °C, and thirty seconds at 72 °C). There was also a cooling step (40 °C for 30 s) and a melting step (95 °C for 10 s and 65 °C for 30 s). The gene expression levels were normalized to the endogenous control glyceraldehyde 3-phosphate dehydrogenase (GAPDH). The fold change relative to the control was computed using the 2 − Δ(ΔCt) technique. The primers used to determine mRNA levels are presented in the Supplement Table. 

hiPSC-derived cornea organoids were fixed with 4% paraformaldehyde (Sigma-Aldrich) for 20 min at room temperature and embedded in cryomatrix (OCT, Fisher Healthcare). 5 μm thick serial sections were blocked/permeabilized with staining solution (1% Bovine Serum Albumin (BSA; Sigma-Aldrich) (w/v) and 0.3% Triton-X100 (neoFroxx) (v/v) in PBS) for 1 h at room temperature. Cells were then treated with primary antibodies at 4 °C overnight and afterward with secondary antibodies for 1 h at room temperature. Nuclei were stained with 4,6-diamidino-2-phenylindole (0.5 μg/ml) (DAPI; Neofroxx). iPSCs were treated with the staining solution for blocking and permeabilization after fixation and staining proceeded as previously explained. Samples were visualized either with Fluorescence Microscopy (Olympus IX71) or Confocal Microscopy (Zeiss LSM880). Antibodies used are listed in the Supplement Table. 

Four distinct iPSC lines (AAK-1, AAK-2, AAK-3, and WT) in their undifferentiated state were subjected to washing with a 0.9% NaCl solution (Sigma-Aldrich) and subsequently fixed using a methanol:water (9:1, v/v) solution (Sigma-Aldrich). The samples were then flash-frozen in liquid nitrogen and scraped into 2 ml centrifuge tubes. Following centrifugation at 15,000 rpm for 15 min at + 4 °C, the supernatants were carefully transferred into 2 ml Eppendorf tubes and stored at − 80 °C until further analysis. Additionally, culture media from AAK organoid samples were collected on days 30, 44, and 60, and stored at − 80 °C. 

Samples stored at − 80 °C were thawed at room temperature and evaporated to complete dryness in a vacuum centrifuge (Savant ISS110 SpeedVac Concentrator, Thermo Scientific). The dried samples were methoxylated with 20 µL of methoxyamine hydrochloride (20 mg/ml in pyridine) by incubation in an oven at 30 °C for 90 min. Subsequently, 50 µL of N-methyl-N-trimethylsilyl trifluoroacetamide containing 1% trimethylchlorosilane (Sigma-Aldrich) was added, and the samples were further incubated at 37 °C for 30 min. The derivatized samples were then transferred into silylated GC–MS vials and analyzed using a GC–MS system (Shimadzu-QP2010 Ultra) equipped with a DB-5MS stationary phase column (30 m + 10 m DuraGuard × 0.25 mm i.d. and 0.25-μm film thickness). The solvent delay was set to 5.90 min. The oven temperature was initially maintained at 60 °C for 1 min, followed by a temperature increase at a gradient of 10 °C/min until reaching 325 °C, which was held for 10 min before cooling. The MSD transfer line temperature was maintained at 290 °C, and the column flow rate was held constant at 1 ml He/min. The mass range was set to 50–650 Da, with a total run time of 37.5 min. Data deconvolution, peak alignment, normalization, and data matrix creation were conducted using MS-DIAL (version 4.8) software. Metabolite identification for GC–MS was performed utilizing a commercially available retention index library (Fiehn Retention Index Library) with a 70% or higher identification cut-off score. The data matrix generated by MS-DIAL was exported to an Excel file. Metabolite traits with more than 50% missing values were excluded from the data matrix. Missing values were imputed with half of the smallest concentration value within the metabolite group. The final data matrix was imported into MetaboAnalyst (version 6.0) software for enrichment analysis and multivariate analyses, including principal component analysis (PCA), partial least squares discriminant analysis (PLS-DA), and heatmap analysis. 

All experiments were represented for three biological replicates. Data are represented as the mean ± SEM or mean ± SD. The GraphPad Prism 8.0 (GraphPad Software, San Diego, USA) program was used for statistical analysis. Image processing was performed using ImageJ and Zen 3.1 software. Metabolomics analyses were carried out by Metaboanalyst 6.0. 

Supplementary Material 1. Supplementary Material 2. Supplementary Material 3. Supplementary Material 4. 

