Mitochondrial Genetic Heterogeneity in Leber’s Hereditary Optic Neuropathy: Original Study with Meta-Analysis

Leber’s hereditary optic neuropathy (LHON) is a mitochondrial disorder that causes loss of central vision. Three primary variants (m.3460G>A, m.11778G>A, and m.14484T>C) and about 16 secondary variants are responsible for LHON in the majority of the cases. We investigated the complete mitochondrial DNA (mtDNA) sequences of 189 LHON patients and found a total of 54 disease-linked pathogenic variants. The primary variants m.11778G>A and m.14484T>C were accountable for only 14.81% and 2.64% cases, respectively. Patients with these two variants also possessed additional disease-associated variants. Among 156 patients who lacked the three primary variants, 16.02% harboured other LHON-associated variants either alone or in combination with other disease-associated variants. Furthermore, we observed that none of the haplogroups were explicitly associated with LHON. We performed a meta-analysis of m.4216T>C and m.13708G>A and found a significant association of these two variants with the LHON phenotype. Based on this study, we recommend the use of complete mtDNA sequencing to diagnose LHON, as we found disease-associated variants throughout the mitochondrial genome.


Introduction
Leber's hereditary optic neuropathy (LHON; OMIM# 535000) was the first disease reported to be associated with defective mitochondrial DNA (mtDNA) [1,2]. LHON is one of the most commonly studied mitochondrial disorders, affecting at least 1 in 14,000 males in Northern England and 1 in 25,000 in the British population [2][3][4]. The condition preferentially affects young adult males and is characterised by acute to subacute, painless central vision loss, which can be unilateral or bilateral [2,[4][5][6].
In >95% of cases, three major point variants, m.3460G>A, m.11778G>A, and m.14484T>C, account for the disease phenotype [7]. These three point variants affect the genes which encode for different mitochondrial subunits (m.3460G>A: MT-ND1, m.11778G>A: MT-ND4, and m.14484T>C: MT-ND6) of complex I of the respiratory chain complex. However, not all individuals who harbour one or more these primary variants develop a LHON phenotype. Despite years of extensive genetic and functional studies on LHON, the exact aetiology and pathophysiology of the disease is not yet fully understood. Existing evidence indicates that additional factors like genetic factors (unknown nuclear genes and mtDNA haplogroups), environmental factors (smoking, alcohol intake, etc.), and epigenetic factors may control the penetrance and phenotypic expression of LHON disorder [2,7].
Previous studies from our lab have shown that the incidence of m.11778G>A among the Indian population is high as compared to that of m.14484T>C, and that both these variants are present in diverse haplogroups [8,9]. Several other independent studies across India have also reported that primary LHON variant m.11778G>A is more prevalent in the Indian population [9][10][11]. However, the role of other secondary variants is less explored in the Indian context.
Extensive research conducted across the world has indicated the importance of haplogroup background on LHON pathogenesis. However, the results are not consistent. For instance, independent studies on European cohorts suggest that haplogroup-J (J2, J1) and K confer increased penetrance of m.11778G>A, m.14484T>C, and m.3460G>A, respectively, while the H haplogroup protects patients carrying the m.11778G>A variant on the European continent [12]. Furthermore, studies on a Chinese cohort showed that M7b1 2 increases the risk of LHON, whereas M8a can impart a protective effect. However, no significant difference in the frequency of M7b1 2 and M8a between case and control subjects was observed, despite their evident effect on the clinical expression of the disease [13]. Later, Saikia et al. (2017), reported for the first time an association of the M haplogroup (p = 0.028) with LHON [14]. In contrast, our previous studies did not show any specific correlation between haplogroup background and LHON expression [8,9,15]. Since the results are ambiguous and inconclusive, it is imperative to investigate the effect of haplogroup with respect to LHON pathogenesis.
Most of the previous studies from India have attempted to address LHON pathogenesis in the context of primary LHON variants. However, in the past few years, several studies across the world have highlighted the importance of other secondary/intermediate variants in LHON pathogenesis. In addition, several pieces of evidence suggest the involvement of other factors like haplogroup background, environmental factors, etc., in LHON pathogenesis. Therefore, in the present study, we intended to understand the involvement of haplogroup(s) and other disease-associated variants in mtDNA, which could directly or indirectly lead to the LHON phenotype. We also performed meta-analysis of the most frequent LHON-associated variants identified in our study to decipher the quantitative association between these two variants and LHON pathogenesis.

Ethical Declaration and Patient Evaluation
Informed written consent was obtained from patients prior to sample collection. The study was approved by the Institutional Ethical Committee (IEC) of the CSIR Centre for Cellular and Molecular Biology (CCMB) and L V Prasad Eye Institute (LVPEI), Hyderabad. All the patients underwent appropriate ophthalmological examinations at LVPEI. An expert team of ophthalmologists categorised the patients based on clinical features of LHON, including sudden, painless acute or subacute vision loss, without any specific reason. All partial optic atrophy cases of unknown aetiology were subjected to comprehensive ophthalmic evaluation, followed by Magnetic resonance imaging (MRI), Neuromyelitis Optica (NMO), Myelin oligodendrocyte glycoprotein (MOG), and vitamin B12 level tests. After typical and atypical optic neuritis, nutritional optic neuropathy, or toxic optic neuropathy were ruled out, patients were subjected to LHON variant or mitochondrial optic neuropathy tests. Patients with other causes of optic neuropathy, such as glaucoma, optic neuritis, ischemic optic neuropathy, trauma, and exposure to drugs such as ethambutol, were also excluded from the study. However, patients who were initially diagnosed with optic neuritis and who did not respond to intravenous steroids were included after the results of genetic tests. We collected~5.0 mL peripheral blood samples in an EDTA vacutainer and stored them at 4 • C until further use.

DNA Isolation and Sequencing of Complete mtDNA
Total genomic DNA was isolated using the standard phenol-chloroform protocol described elsewhere [16]. The complete mitochondrial genome of each patient was amplified in 24 overlapping fragments [17]; the amplicons were purified using Exo-SAP-IT (USB, Affymetrix, Santa Clara, CA, USA) at 37 • C for 15 min followed by 15 min incubation at 80 • C, and then sequenced using the ABI BigDye Terminator cycle sequencing kit (Applied Biosystems, Foster City, CA, USA). Extended PCR products were precipitated using a 1:5 mixture of 5M sodium acetate and ethanol, followed by washing with 80% distilled alcohol [8]. After gentle drying, the purified products were dissolved in Hi-Di formamide and analysed in the ABI 3730 DNA Analyzer (Applied Biosystems) [9]. The obtained DNA sequences of the patients were edited using Sequencing Analysis v5.2 (Applied Biosystems) and assembled with a mitochondrial reference sequence (revised Cambridge Reference Sequence with GenBank number (NC_012920)) [18,19] using AutoAssembler sequence assembly software (Applied Biosystems, Foster City, CA, USA), and the variants were noted.

Meta-Analysis of Primary LHON Variants
A literature search was done for all three primary LHON variants (m.3460G>A, m.11778G>A, and m.14484T>C) using the keyword "LHON" and respective primary variant number and/or details. The frequency data from each study were extracted. Three individual matrices were created for each of the three primary variants using variant frequency and their latitude and longitude. Furthermore, an interpolated frequency spectrum was constructed using the Kriging gridding method of the automap package in R and the geospatial frequency distribution was plotted on the world map using the ggplot2 package in R.

Search and Selection of Relevant Studies
A meta-analysis was performed for the two most frequent variants identified in our study (m.4216T>C and m.13708G>A) as per the protocol described in our previous study [22]. Studies on the m.4216T>C and m.13708G>A variants in LHON were selected after a systematic literature search through databases including Pubmed and Google Scholar using the keywords "4216 and LHON" and "13708 and LHON", respectively. Articles published in the English language and without restriction of publication date up to March 2021 were included in the study. Furthermore, studies were selected based on the following inclusion criteria: (1) clear information about the variants in both cases and controls, (2) high-resolution SNP identification, and (3) inclusion of patients as per the standard guidelines. Exclusion criteria were as follows: (1) studies lacking detailed genotypic data, (2) conference abstracts, and (3) reports lacking accurate description of the study.

Data Extraction
All the articles fulfilling the inclusion conditions were read thoroughly, and details like first author, year of publication, country of origin, total number of cases and controls, and frequency of variants were mined.

Statistical Analysis
Meta-analysis was executed using the meta package in R (v4.0.3). heterogeneity was calculated using the Q test [23]. The levels of heterogeneity were determined based on Higgins and Thompson classification, i.e., 25%, 50%, and 75% I 2 values correspond to a low, medium, and high level of heterogeneity, respectively [24]. The odds ratios were calculated using both fixed-effect (Mantel-Haenszel method) and random-effect (DerSimonian-Laird estimator) models. High-resolution forest plots were made using the forest function in the meta package in R (v1.3.10) to estimate the odds ratios from both fixed and random models. A p-value of <0.01 was considered to be statistically significant.
To investigate the publication bias in the studies involved in the meta-analysis, we generated funnel plots using the funnel function in meta package of R. The symmetry of these funnel plots suggests the presence or absence of publication bias. Furthermore, quantitative estimation of publication bias was done using Egger's regression intercept test using the metabias function in R. This test is a weighted linear regression of the standardized effect sizes on their corresponding precisions. The estimated regression intercept divided by standard error follows the t-distribution, which provides the p-value of the regression test. We considered a p-value of <0.01 to be statistically significant [25]. All statistical analysis for the study was done in R.

Study Samples
A total of 189 LHON patients were recruited in the study, including 172 males and 17 females. Of these 189 probands, two individuals were from one family, and the remaining 187 were sporadic cases. The average age of the patients was 25.51 years (SD ± 11.62 years, ranging from 1 year to 60 years). Clinical features of the patients with primary LHON variants are provided in Table 1.

Whole Mitochondrial Genome Screening
We sequenced the complete mitochondrial genomes of 189 LHON patients and identified two primary LHON variants, i.e., m.11778G>A and m.14484T>C, in our cohort. The inclusive frequency of these two primary LHON variants was 17.46 (33/189). Within the sample, 14.81% of patients (28/189) were positive for the m.11778G>A variant, whereas the m.14484T>C variant was present in 2.64% (5/189) patients (Table 2). In addition, out of 33 patients with primary LHON variants, 13 patients (11 with m.11778G>A and 2 with m.14484T>C) harboured additional variants associated with the spectrum of disorders (Table 3). Interestingly, five of these variants [m.4216T>C (p.Y304H), m.5460G>A (p.A331T), m.11253T>C (p.I165T), m.13708G>A (p.A458T), and m.15927G>A MT-tRNA-Threonine] were found to be associated with LHON, and might influence patient phenotype. Since these variations are associated with the spectrum of mitochondrial disorder, it is possible that these variants, either alone or in combination with primary LHON variants, can lead to clinical variability in the LHON patients.

Haplogroup Analysis of Patients with Primary LHON Variants
The mtDNA-based haplogroup analysis suggested that patients harbouring the m.11778G>A variant belonged to six different haplogroups (D, H, M, R, U, and W), and patients with m.14484T>C variants belonged to three haplogroups (I, M, and U) ( Table 2). About 60.60% (20/33) of patients with primary LHON variants were in the M haplogroup.

Haplogroup Analysis of Patients with Secondary LHON-Associated Variants
Patients with LHON-associated variants belonged to eight different haplogroups (D, F, H, M, R, T, U, and W) ( Table 4). We found that a total of 40% (10/25) of patients with secondary LHON-associated variants were in the M haplogroup, and these were further categorised into the following sub-haplogroups: M43b (8%, 2/25), and 1 each of M2a1b (4%,

Meta-Analysis of LHON-Associated Variants
The literature search resulted in 495 articles for m.4216T>C and 452 articles for m.13708G>A. After initial screening based on the title and presence of keywords in the articles, a total of 73 articles were reviewed for m.4216T>C and 71 for m.13708G>A. For the final meta-analysis, we included five articles related to m.4216T>C [26,52,[63][64][65] and eight articles related to m.13708G>A [49,61,[63][64][65][66][67][68]. Other reviewed articles (68 for m.4216T>C and 63 for m.13708G>A) were excluded from the study as a result of (a) the absence of relevant information, (b) lack of relevance of the study to this analysis, (c) lack of control data, or (d) unavailability of full articles. Including our current study, the m.4216T>C and m.13708G>A meta-analyses covered six and nine datasets, respectively,

Meta-Analysis of LHON-Associated Variants
The literature search resulted in 495 articles for m.4216T>C and 452 articles for m.13708G>A. After initial screening based on the title and presence of keywords in the articles, a total of 73 articles were reviewed for m.4216T>C and 71 for m.13708G>A. For the final meta-analysis, we included five articles related to m.4216T>C [26,52,[63][64][65] and eight articles related to m.13708G>A [49,61,[63][64][65][66][67][68]. Other reviewed articles (68 for m.4216T>C and 63 for m.13708G>A) were excluded from the study as a result of (a) the absence of relevant information, (b) lack of relevance of the study to this analysis, (c) lack of control data, or (d) unavailability of full articles. Including our current study, the m.4216T>C and m.13708G>A meta-analyses covered six and nine datasets, respectively, which comprised 658 cases and 707 controls for m.4216T>C, and 675 cases and 1551 controls for m.13708G>A.
The frequencies of m.4216T>C and m.13708G>A were compared between cases and controls. We identified significant amounts of heterogeneity in the studies included for our meta-analysis: (a) m.4216T>C; p = 0.0063, Q = 16.19, df = 5, tau 2 = 1.0139, tau = 1.0069, I 2 = 69.1%, H = 1.80; (b) m.13708G>A; p = 0.0006, Q = 27.59, df = 8, tau 2 = 0.6999; tau = 0.8366; I 2 = 71.0%; H = 1.86. The meta-analysis was done using both fixed-and random-effect models to give an unbiased estimate of effect size (odds ratio). However, given the significant heterogeneity we found in our analysis, the random-effect model gave the best estimate of the effect size. The pooled odds ratios calculated using fixed-and random-effect models were (m.  The publication bias was assessed using a funnel plot and Egger's test. The funnel plot (Figure 3) indicated no publication bias, which was quantitatively confirmed by Egger's test (p4216 = 0.6545; p13708 = 0.3448).

Discussion
In more than 95% of the LHON cases, three point variants (m.3460G>A, m.11778G>A, and m.14484T>C) in the mtDNA are associated with disease condition [7]. Unfortunately, a substantial number of LHON cases remain undiagnosed due to the absence of these primary variants. In addition, due to incomplete penetrance and sex bias, the diagnosis of LHON remains challenging. Therefore, it is apparent that there are additional genetic factors such as secondary variants in mtDNA, defects in nuclear DNA (nDNA), ethnic background, and/or epigenetic factors (including environmental factors and lifestyle) that may influence the expression and/or penetrance of the disease.
Previous studies from our lab showed that the frequencies of m.11778G>A and m.14484T>C were 29.22% and 4.2%, respectively and found significant sex bias [8,9]. In alignment with our previous studies, we pooled 189 LHON samples from different parts of India and screened whole mtDNA to find the existing and novel mtDNA variants in our cohort. We first aimed to determine the frequencies of three primary LHON variants (m.3460G>A, m.11778G>A, and m.14484T>C) and to further corroborate the effect of haplogroup on the clinical expression of disease. We observed the presence of only two primary variants, i.e., m.11778G>A and m.14484T>C in 14.81% and 2.64% of the patients, respectively; this was much lower than in our previous reports i.e., 29.22% (m.11778G>A)

Discussion
In more than 95% of the LHON cases, three point variants (m.3460G>A, m.11778G>A, and m.14484T>C) in the mtDNA are associated with disease condition [7]. Unfortunately, a substantial number of LHON cases remain undiagnosed due to the absence of these primary variants. In addition, due to incomplete penetrance and sex bias, the diagnosis of LHON remains challenging. Therefore, it is apparent that there are additional genetic factors such as secondary variants in mtDNA, defects in nuclear DNA (nDNA), ethnic background, and/or epigenetic factors (including environmental factors and lifestyle) that may influence the expression and/or penetrance of the disease.
Previous studies from our lab showed that the frequencies of m.11778G>A and m.14484T>C were 29.22% and 4.2%, respectively and found significant sex bias [8,9]. In alignment with our previous studies, we pooled 189 LHON samples from different parts of India and screened whole mtDNA to find the existing and novel mtDNA variants in our cohort. We first aimed to determine the frequencies of three primary LHON variants (m.3460G>A, m.11778G>A, and m.14484T>C) and to further corroborate the effect of haplogroup on the clinical expression of disease. We observed the presence of only two primary variants, i.e., m.11778G>A and m.14484T>C in 14.81% and 2.64% of the patients, respectively; this was much lower than in our previous reports i.e., 29.22% (m.11778G>A) [9] and 4.2% (m.14484T>C) [8], as well as in another Indian study (3.33%, m.14484T>C) [11]. Detection of low frequencies of primary LHON variants in India could be due to the unique genetic architecture of the Indian population [69,70]. Furthermore, similarly to our previous study [9], we observed that the m.11778G>A variant was present at a higher frequency than the m.14484T>C variant. This suggests that the m.11778G>A variant is primarily responsible for LHON pathogenesis in the Indian population. Both these variants were absent in our control samples (n = 192). A few studies have reported the co-existence of m.11778G>A and m.14484T>C primary variants in a single patient [10,71]; however, in our study, co-occurrence of these two primary variants were not observed. Similarly to our previous studies [8,9], we did not identify m.3460G>A variation in the current cohort. Interestingly, in the present study, we observed that only males harboured both the primary variants, which correlates to the sex bias reported to be a critical determining factor in vision loss [9].
Furthermore, in 11 out of 28 patients with the m.11778G>A variant, we observed additional variants, including two synonymous, nine nonsynonymous, three tRNA, and four rRNA variants associated with the spectrum of disease. Interestingly, six of these variants were associated with LHON. In the patients harbouring LHON-specific primary and other LHON-associated variants, one patient (P2) harboured two LHON-associated variants, i.e., m.13708G>A and m.15927G>A, along with m.11778G>A variants and displayed early onset of disease (age: 10 years). Another patient (P16) harboured m.13708G>A along with m.11778G>A and presented with vision loss along with tiredness and fatigue, pain in the upper and lower limbs, and occasional paraesthesia. Therefore, the varying clinical presentation in these patients (P2 and P16) could be due to synergistic effects of the other LHON-associated variants. In three patients (P23, P24, and P27), other LHON-associated variants (m.4216T>C, m.11253T>C, and m.5460G>A, respectively) were observed along with m.11778G>A; however, they presented only LHON-specific phenotypes. Similarly to our previous reports, we did not find secondary LHON-associated variants in patients harbouring the m.14484T>C variant [8]. However, in two patients (P32 and P33), we identified two variants, m.2361G>A and m.921T>C, respectively, both associated with left ventricular noncompaction cardiomyopathy (LVNC), although neither of these patients showed any heart-related complications. Nonetheless, the identified variations are associated with the spectrum of mitochondrial disorders, and the presence of these variations with primary LHON variants could impact the clinical variability of LHON patients.
Haplogroup analysis based on complete mtDNA sequencing suggested that individuals harbouring two primary variants belong to diverse haplogroup backgrounds. Individuals with the m.11778G>A variant belonged to six major haplogroups (D, H, M, R, U, and W), while individuals with the m.14484T>C variant belonged to three different haplogroups (I, M, and U). Interestingly, except in one case, there was no overlap in the haplogroup of patients harbouring these two primary variants; this suggests that these variants exist in different haplogroup backgrounds which may alter the expression of LHON. Furthermore, among the different haplogroups identified in our patients harbouring primary LHON variants, the most predominant was the M haplogroup (60.60%), followed by R (15.15%), U (9.09%), W (6.06%), D (3.03%), H (3.03%), and I (3.03%). The predominance of the M haplogroup could be due to the ubiquitous existence of the M haplogroup, accounting for >70% of mtDNA lineage in the Indian population [72]. The presence of such diverse haplogroups among the patients carrying both the primary variants suggests that these patients are maternally unrelated, and no particular haplogroup is associated with the LHON phenotype. Similarly to our previous report, we did not find the J haplogroup (which has been reported to increase the disease expression in European cohorts) in the patients carrying primary LHON variants [9]. However, we identified the H haplogroup in one patient (P7), which is reported to reduce the chance of visual failure in a patient with the m.11778G>A variant [12] and could be associated with the late onset of disease in this patient.
Of the 156 patients lacking primary LHON variants, we identified secondary LHONassociated variants in 16.02% (25/156) individuals. A total of 17 unique LHON-associated variants were identified either alone or in combination with other disease-associated variants. The two most predominant LHON-associated variants identified in our cohort were m.4216T>C in the MT-ND1 gene and m.13708G>A in the MT-ND5 gene, each present in 16% (4/25) of the patients. In addition, three variants, m.3316G>A, m.11696G>A, and m.9966G>A, were identified as the second most common variants, and each of them was present in 12% (3/25) of the patients. As in other studies [52,73], a maximal number of variants (53.4%) was observed in MT-ND gene, suggesting the significant involvement of the MT-ND gene in disease pathogenies. Since the LHON-associated variants were linked with several diseases, including LHON, and some of these variants are haplogroup-defining, it is likely that these variants might have a synergistic effect on LHON pathogenesis. However, the frequency of these LHON-associated variations varies from one population to another and they are also reported in control samples in varying frequencies. Therefore, it is critical to consider and evaluate the role of these LHON-associated variations in LHON pathogenesis. Haplogroup analysis in patients with secondary LHON-associated variants suggested that patients belonged to eight different haplogroups (D, F, H, M, R, T, U, and W), with the majority (40%) of them belonging to the M haplogroup, followed by the R (20%), T (8%), F (8%), U (8%), W (8%), D (4%), and H (4%) haplogroups. The frequencies of some of these haplogroups in the Indian population are M (>70%) [72], H (16%), R (1.86%), T (7.46%), U (38.6%), and W (8.7%) [74]. Since the frequencies of identified haplogroups in this study and previous studies from the Indian population are comparable, it is unlikely that the haplogroup of an individual might play a significant role in the expression of LHON; rather, this can be interpreted as the LHON-specific variants having arisen in different haplogroup backgrounds.
In addition, we also identified other disease-associated variants in 22 patients. These patients did not harbour any known primary or other LHON-associated variants. Interestingly, approximately 50% of the identified variants in these patients were associated with hearing loss. We also identified variants associated with CPEO, MELAS, myopathy, intellectual disability, heart diseases, etc. We have seen similar genetic heterogeneity in other diseases, including cardiovascular diseases, neuromuscular diseases, and several mitochondrial disorders. The genotype-phenotype did not correlate in these patients; however, due to the clinical and genetic heterogeneity in mitochondrial disorders, the involvements of these variants must not be ignored and should be further investigated.
In this study, using the mtDNA-sequencing approach, we identified variants (including primary variants, LHON-associated variants, and other disease-associated variants) in 42.32% (80/189) individuals; the remaining 57.67% (109/189) of patients lacked any pathogenic variants in the mitochondria. It is possible that the patients without mtDNA variants either were not affected with a mitochondrial disorder or had abnormalities specific to the optic nerve, resembling mitochondrial myopathy [75]. Another reason could be a very low level of heteroplasmy which could go undetected by conventional PCR approaches [27,39,73]. Furthermore, considering the varied level of penetrance of LHON, studying the effects of other factors like nuclear gene defects, mtDNA copy number, haplogroup heterogeneity, the existence of higher mutation load (both nonsynonymous and synonymous variants), and environmental factors could be interesting and might help to better our understanding of the disease.
Meta-analysis was performed for the two most frequent variants identified in our study samples, i.e., m.4216T>C and m.13708G>A. Based on the pooled odds ratio, we found a significant association of m.4216T>C and m.13708G>A with the LHON phenotype. Therefore, it is possible that these two variants might have a synergistic effect on LHON expression in the presence of other variants or under the influence of epigenetic factors.
In summary, the results of our current study support the data previously published by our lab. We have reported primary LHON variants, several secondary LHON-associated variants, and other disease-associated variants that could be responsible for the disease phenotype. Although these variants are known to be strongly associated with disease phenotypes, the evidence regarding their contribution to LHON or other diseases is limited and therefore requires further investigation. In addition, our study highlights the importance of screening the whole mitochondrial genome for a comprehensive diagnosis of LHON and always examining the role of other important factors like variants in the nuclear genes and/or epigenetic factors for a better understanding of the disease mechanism.
Author Contributions: K.T. conceived and designed the study and provided the reagents. A.P., V.Y.V., and R.K. evaluated LHON patients and provided clinical materials. R.K.J. performed all the experiments. R.K.J., C.D., Q.H., and G.G. analysed the data. R.K.J. and K.T. wrote the manuscript with the inputs from all the co-authors. All authors have read and agreed to the published version of the manuscript. Data Availability Statement: Data exist in CCMB repository and would be made available to the researchers, upon request.