Effectiveness of whole exome sequencing in unsolved patients with a clinical suspicion of a mitochondrial disorder in Estonia.

OBJECTIVE: Reaching a genetic diagnosis of mitochondrial disorders (MDs) is challenging due to their broad phenotypic and genotypic heterogeneity. However, there is growing evidence that the use of whole exome sequencing (WES) for diagnosing patients with a clinical suspicion of an MD is effective (39-60%). We aimed to study the effectiveness of WES in clinical practice in Estonia, in patients with an unsolved, but suspected MD. We also show our first results of mtDNA analysis obtained from standard WES reads.
METHODS: Retrospective cases were selected from a database of 181 patients whose fibroblast cell cultures had been stored from 2003 to 2013. Prospective cases were selected during the period of 2014-2016 from patients referred to a clinical geneticist in whom an MD was suspected. We scored each patient according to the mitochondrial disease criteria (MDC) (Morava et al., 2006) after re-evaluation of their clinical data, and then performed WES analysis.
RESULTS: A total of 28 patients were selected to the study group. A disease-causing variant was found in 16 patients (57%) using WES. An MD was diagnosed in four patients (14%), with variants in the SLC25A4, POLG, SPATA5, and NDUFB11 genes. Other variants found were associated with a neuromuscular disease (SMN1, MYH2, and LMNA genes), neurodegenerative disorder (TSPOAP1, CACNA1A, ALS2, and SCN2A genes), multisystemic disease (EPG5, NKX1-2, ATRX, and ABCC6 genes), and one in an isolated cardiomyopathy causing gene (MYBPC3). The mtDNA point mutation was found in the MT-ATP6 gene of one patient upon mtDNA analysis.
CONCLUSIONS: The diagnostic yield of WES in our cohort was 57%, proving to be a very good effectiveness. However, MDs were found in only 14% of the patients. We suggest WES analysis as a first-tier method in clinical genetic practice for children with any multisystem, neurological, and/or neuromuscular problem, as nuclear DNA variants are more common in children with MDs; a large number of patients harbor disease-causing variants in genes other than the mitochondria-related ones, and the clinical presentation might not always point towards an MD. We have also successfully conducted analysis of mtDNA from standard WES reads, providing further evidence that this method could be routinely used in the future.

Introduction

Mitochondrial disorders (MDs) are one of the most frequent genetic disorders, with a prevalence of ~12.5 per 100,000 in adults [1] and ∼4.7 per 100,000 in children [2]. In Estonia, the live-birth prevalence for MDs (estimated during 2003–2009) was 1/20,764 live births [3]. However, diagnosing MDs is a very challenging task due to their extremely broad phenotypic and genotypic heterogeneity. They arise due to a primary mitochondrial dysfunction and can affect various functional systems; most often the central nervous system, neuromuscular system, and heart and liver, causing isolated or multisystemic syndromes. The underlying genetic cause can be a mitochondrial DNA (mtDNA) point mutation, single large-scale mtDNA deletion, or a variant in a nuclear gene, resulting in defective mtDNA maintenance, transcription, protein translation, or a defective ancillary process. Of the ∼1300 proteins in the mitoproteome encoded by nuclear DNA (nDNA), mutations have been reported in >250 genes [4]. Therefore, different scoring systems to evaluate the probability of MDs have been developed [[5], [6], [7], [8], [9]]. These scoring systems take into account the clinical features, biochemical changes, brain imaging anomalies, as well as pathomorphological findings in muscle, and enzymatic changes in muscle or fibroblasts.

Nowadays, as next generation sequencing technologies become more widely available for use in routine clinical practice, a vast number of studies have been carried out world-wide to establish the effectiveness of whole exome sequencing (WES) for the diagnosis of different groups of inherited diseases. Studies on the effectiveness of WES in MDs have shown a diagnostic yield of 39–60% [10,11]. In addition, these studies have shown that many patients suspected of having an MD, actually have a different disorder, usually neurological, which have great overlap in the phenotypes. Therefore, Morava and Brown, 2015 [12] have proposed a new diagnostic algorithm – sequencing mtDNA from blood and performing WES analysis, before taking a muscle/skin biopsy.

Our study aimed to evaluate the diagnostic effectiveness of WES in clinical practice, in Estonia for patients with an unsolved metabolic disease and a clinical suspicion of a mitochondrial disorder. We discuss these patients' phenotype, results of metabolic, enzymatic, and imaging studies, the pathomorphology of their muscle, and genotype. We also show our first results of mtDNA variant calling and interpretation from standard WES reads.

Methods

This study was approved by the Research Ethics Committee of the University of Tartu (approval date 17/11/2014 and number 242/M-10). Informed consent was obtained from all patients and/or their parents in the study group, except from those who were deceased.

Patient selection

The retrospective study group was selected from a database of 181 patients whose fibroblast cell cultures had been stored in the Department of Clinical Genetics, Tartu University Hospital, during the period of January 2003 – December 2013. These 181 patients had been consulted by a clinical geneticist and had a definite diagnosis or a suspicion of a metabolic disease. The fibroblasts were cultured from a skin biopsy to carry out necessary enzyme analyses, or for future diagnostic purposes. At the same time in most of the cases the DNA sample was also stored.To our study group, we included only those patients whose physician had a strong suspicion of a MD, but previous genetic investigations (targeted mtDNA and single nuclear gene analysis) resulted in no findings. WES analysis was not performed on any of these patients because it was not available in Estonia prior to 2013.

The prospective study group was selected during January 2014 – March 2016. We included patients who were referred to a clinical geneticist in the Department of Clinical Genetics, Tartu University Hospital, had an onset of the disease in childhood, and the clinical geneticist had a clinical suspicion of an MD.

Next, we re-examined all the medical history charts of each included patient to obtain their clinical signs and symptoms, results of metabolic and imaging studies, and other medical investigations, including muscle biopsy. The retrospective patients who were alive were invited for a follow-up genetic consultation. From the patients who were deceased, we had formerly taken a blood sample for DNA extraction and their clinical history was well described. We scored each patient according to the mitochondrial disease criteria (MDC) developed by Wolf and Smeitink, [8] and Morava et al. [9]. WES analysis was performed on each patient.

Whole exome sequencing

The DNA of all the patients was extracted from blood lymphocytes. The WES analysis was carried out using different methods due to availability and cost-effectiveness at various time points. The WES enrichment kits used were SureSelect XT Human All Exon v5 (Patients 2, 3, 5, 6, 8, 10, 11, 12, 17, 18, and 19) or v4 (Patient 1) enrichment kit (Agilent Technologies, Santa Clara, CA), Nextera Rapid Capture Exome 37 Mb kit (Illumina Inc., San Diego, CA, USA) (Patients 4, 7, 9, 14, 15, 16, 20, 21, 22, 23, 24, 25, 26, 27, and 28), and the TruSeq Exome Enrichment Kit 62 Mb (Illumina Inc., San Diego, CA, USA) (Patient 13). All patients were sequenced as proband only, except Patients 7 and 28, in whom offspring-parent trios were sequenced, and Patient 13 who had two affected sibs that were sequenced. After the sequencing, the reads were aligned to hg19 reference genome using Burrows-Wheeler Aligner (BWA) [13]. Further data processing, variant calling, and annotation was performed following Genome Analysis Toolkit (GATK) best practice guidelines [14] using Picard, GATK [15,16], Annovar [17], and SnpSift [18]. All reported variants were confirmed, and familial segregation analysis was done using Sanger sequencing.

mtDNA analysis from standard WES reads

Every patients' reads (except patient 1) generated during standard WES (as described in the section 2.2), but mapped to mtDNA, were specifically analyzed to detect mtDNA variants, and to assess the feasibility of investigating mtDNA from standard WES reads. In order to study mtDNA, SAMtools [19] were used to extract reads mapped to chrM. These reads were aligned to the b37 reference genome (GRCh37 including rCRS mitochondrial sequence) using BWA MEM [20], and subsequently sorted and marked for duplicates by Picard. The variants were called by combining SAMtools [19] mpileup with VarScan [21], and then annotated using Annovar [17]. The analysis of guided variants along with haplogroup calling was guided by Mitomap and Mitomaster [22]. The mtDNA coverage was calculated using the GATK DepthOfCoverage tool.

Results

Out of the 181 patients in the fibroblast database, we selected 21 patients for the retrospective group, of which four patients were excluded due to consent withdrawal, loss of clinical symptoms, no contacts with the family and/or no stored DNA. The prospective group consisted of 11 patients, who were all included. In total, we had 28 patients that we performed WES analysis on.

The main clinical, metabolic, imaging, and myopathological characteristics of all the patients are described in the Supplementary data (Table S1). The onset of disease ranged from prenatal to seven years of age, with onset at birth or early infancy in the majority of patients. The current ages ranged from 1 y to 29 y. Nine patients had deceased between the ages of 3 m to 26 y. Motor developmental delay was present in 16, and intellectual disability in 11 patients. Developmental regression was seen in three patients. Muscle weakness was evident in 16 patients. Eleven patients had muscular hypotonia and decreased/absent deep tendon reflexes (DTRs), whereas 14 patients presented with spasticity and/or elevated DTRs/positive Babinski sign. Other neurological manifestations, such as dystonia, choreoathetosis, ataxia, tremor, nystagmus, and strabismus occurred in 18 patients. Epilepsy was diagnosed in nine patients. The heart, vision, and gastrointestinal system were involved in seven patients, the skeletal system in five, liver in four, hearing in three, thyroid in two, and the peripheral nerves in one patient. Lactate in serum and/or cerebrospinal fluid was elevated in 13 patients, and tricarboxylic acid (TCA) cycle intermediates in urine were elevated in eight patients. Brain atrophy was detected in ten patients, white matter injury in six, and hypo/aplasia of the corpus callosum was found n in three patients upon MRI.

A muscle biopsy was taken from 18 patients. Unfortunately, the SDH/COX enzyme reaction analysis, which is one of the most important techniques to suggest a mitochondriopathy, was not performed in the majority of cases; however, Gomori trichrome stain, NADH-TR enzyme reaction, and electron microscopy were used to analyze the mitochondria. H&E, PAS, ORO staining, nonspecific esterase, ATPase reactions, and immunohistochemistry for alpha-4-laminin, macrophages, dystrophin, merosin, alpha-, beta- and gamma-sarcoglycans were also used. We saw ragged red fibers (RRF), peripheral increased staining on Gomori trichrome or NADH-TR enzyme reaction, or mitochondrial hyperplasia on electron microscopy in five patients. Mitochondrial ultrastructural changes were found in four patients. SDH/COX reaction was performed in two patients with normal results. In addition, nine patients had accumulations of lipid droplets in many fibers. Myofibrillar disorganization was present in two patients, and denervated fibers were seen in one patient. Patient 4 had a characteristic finding of spinal muscular atrophy (SMA) with large and small group atrophy of both fiber types, but nearly absent COX reaction (which was repeated and compared with a control). Patient 9 demonstrated vast fiber size variability with atrophic, regenerating, and vacuolated fibers, and massive inclusions with heterogenous structure, seen using electron microscopy. These features are very similar to the ones described by Hedberg-Oldfors et al., [23] in patients with Vici syndrome, a new subtype of vacuolar myopathies. Two patients did not have any pathological features on muscle biopsy.

The MDC scores and grouping (definite, probable, possible, unlikely), and the outcomes of the WES analyses are illustrated in Tables 1A and 2. We confirmed disease-causing gene variants in 16 of our patients (57%). An MD was diagnosed in four patients (14%), with variants in the SLC25A4, POLG, SPATA5, and NDUFB11 genes. Other variants found were associated with a neuromuscular disease (SMN1, MYH2, and LMNA genes), a neurodegenerative disorder (TSPOAP1, CACNA1A, and ALS2 and SCN2A genes), a multisystemic disease (EPG5, NKX1–2, ATRX, and ABCC6 genes), and one compound heterozygous variant in the MYBPC3 gene; an isolated cardiomyopathy-causing gene that is typically autosomal dominant, although autosomal recessive cases have been reported with a more severe phenotype [24], as was the case in our Patient 28. Patient 3 had no findings on WES, but a mtDNA point mutation was found in the MT-ATP6 gene upon subsequent mtDNA sequencing analysis from frozen muscle tissue (Table 1B).

Table 1

The pathogenicity, function, and disease causality of the variants found during our study.

Pt a	MDC	Gene (transcript)	Nucleotide and Amino Acid changes	Pathogenicity classified by the ACMG variant interpretation guidelines [25]	Function of the gene and protein, disease causality	Referencea
A. Retro- and prospective patients included in the study for assessing effectiveness of WES
1	11	SLC25A4(NM_001151.3)	Heterozygous variant:c.239G > A, p.(Arg80His)	Pathogenic variant:- Absent in population databases (moderate evidence); - Align GVGD, SIFT, MutationTaster predicts a deleterious effect, conserved region (supporting evidence); - Associated with loss of steady-state ACC1 levels in skeletal muscle and severely impaired ATP transport activity in vitro; decreased mtDNA copy number (~34%), and RCC I and IV activity in muscle, in our patient (strong evidence); - de novo (moderate evidence)	Encodes an ADP/ATP carrier isoform 1 (ACC1) protein that imports ADP into the mitochondrion and exports ATP into the intermembrane space. Its dysfunction causes insufficient adenine nucleotide availability for dATP synthesis and imbalanced dNTP pools, leading to mtDNA depletion. Its defects are associated with mitochondrial DNA depletion syndrome type 12 (OMIM 615418), adPEO (OMIM 609283), and a third distinct phenotypic group described by us	Thompson et al. [27]
2	10	POLG (NM_002693.2)	Compound heterozygous variants:1) c.1377G > A, p.(Met459Ile) 2) c.3154G > C, p.(Gly1052Arg)	1) Likely pathogenic variant:- Absent in population databases (moderate evidence); - PolyPhen-2, CADD predicts a deleterious effect, conserved region (supporting evidence); - Is located in the polymerase domain (http://polg.bmb.msu.edu), RRF on muscle histology in our patient (moderate evidence); - Detected in trans with a pathogenic variant (moderate evidence); - Patient's phenotype highly specific for gene (supporting evidence). 2) Likely pathogenic variant:- Absent in population databases (moderate evidence); - Is located in the polymerase domain, pathogenic partitioning loop cluster (http://polg.bmb.msu.edu) (moderate evidence); - Detected in trans with a pathogenic variant (moderate evidence); - Patient's phenotype highly specific for gene (supporting evidence)	Encodes the catalytic subunit of DNA polymerase gamma, which is essential for mtDNA replication and repair. Its deficiency can cause mtDNA depletion and/or multiple deletions and is associated with Alpers syndrome (OMIM 203700), MNGIE syndrome (OMIM 613662), mitochondrial recessive ataxia syndrome (OMIM 607459), adPEO (OMIM 157640), and arPEO (OMIM 258450).	Hikmat et al. [42]; http://polg.bmb.msu.edu
4	8	SMN1 (NM_000344.3)	Compound heterozygous variants:1) Deletion of exons 7 and 8 on one allele (detected by MLPA analysis) 2) c.410dup, p.(Asn137Lysfs*11)	3) Pathogenic variant:- the most common disease-causing variant; - 2–5% of patients are compound heterozygous for a deletion of at least exon 7, and an intragenic variant. 4) Pathogenic variant:- Absent in population databases (moderate evidence); - Predicted null variant (very strong evidence); - Detected in trans with a pathogenic variant (moderate evidence).	Encodes the SMN protein, which is directly and indirectly via spliceosomal snRNP biogenesis (critical for snRNP assembly) involved in general cellular RNA processing. Its deficiency is associated with spinal muscular atrophy (OMIM 253300).	Lorson and Androphy, [43]
5	7	SPATA5 (NM_145207.2)	Compound heterozygous variants:1) c.250C > T, p.(Arg84*) 2) 2) c.989_991del, p.(Thr330del)	1) Pathogenic variant:- One heterozygous carrier in ExAC, but no homozygous carriers (moderate evidence); - Predicted null variant (very strong evidence); - Detected in trans with a pathogenic variant (moderate evidence). 2) Pathogenic variant:- Most common disease-causing variant reported in multiple patients; - Detected in trans with a pathogenic variant.	Encodes a member of the AAA (ATPase Associated with diverse Activities) protein subfamily. SPATA5 deficiency leads to imbalance of mitochondrial fusion-fission events, mitochondrial shortening, decreased ATP production at axonal endings, and impaired axogenesis in primary cortical neurons in vitro. It results in an AR syndrome with severe global developmental delay, severe speech impairment, hearing loss, abnormal EEG, and microcephaly (OMIM 616577).	Puusepp et al. [28]
7	6	NDUFB11 (NM_019056.6)	Hemizygous variant:c.328C > T, p.(Pro110Ser)	Likely pathogenic variant:- Absent in population databases (moderate evidence); - Align GVGD, SIFT, MutationTaster predicts a deleterious effect, conserved region (supporting evidence); - Decreased RCC I enzyme activity in our patient, located in the NADH domain (moderate evidence); - de novo (strong evidence).	Encodes a component of the RCC I. Its defects are associated with mitochondrial complex I deficiency (OMIM 252010), and linear skin defects with multiple congenital anomalies (OMIM 300952).	Carroll et al. [44]
8	6	MYH2 (NM_017534.5)	Heterozygous variant:c.5609 T > C, p.(Leu1870Pro)	Pathogenic variant:- Previously reported in a patient with a similar phenotype to our patient; - Located in the highly conserved LMM domain within the myosin tail, substitution by proline introduces a structural distortion of the protein that affects the formation of myosin homodimers; - de novo.	Encodes the myosin heavy chain isoform that is expressed in fast type 2A muscle fibers, and is essential for muscle contraction. Its defects are associated with AR and AD proximal myopathy and ophthalmoplegia syndrome (OMIM 605637), and an unusual phenotype of neonatal-onset congenital myopathy matching our patient's phenotype.	D'Amico et al. [45]
9	6	EPG5 (NM_020964.2)	Homozygous variant:c.6690delT, p.(Asn2230Lysfs*6)	Pathogenic variant:- Absent in population databases (moderate evidence); - Predicted null variant (very strong evidence); - Patient's phenotype and muscle morphotype highly specific for gene (supporting evidence).	Encodes a protein with a key role in the autophagy pathway. Its deficiency is associated with Vici syndrome (OMIM 242840).	Cullup et al. [46]
16	5	LMNA (NM_170707.2)	Heterozygous variant:c.207_209delGGT, p.(Val70del)	Likely pathogenic variant:- Absent in population databases (moderate evidence); - Protein length changing variant (moderate evidence); - de novo (moderate evidence).	Encodes the structural proteins, lamin A and lamin C, which underlie the inner nuclear membrane and determine nuclear shape and size. Its defects are associated with various diseases, including EDMD 2 (OMIM 181350), EDMD 3 (OMIM 616516), CMD (OMIM 613205), LGMD 1B (OMIM 159001), cardiomyopathies, lipodystrophies, peripheral neuropathies, and progeria syndromes.	Fisher et al. [47]
17	5	CACNA1A (NM_023035.2)	Compound heterozygous variants:1) c.4315T>A, p.(Trp1439Arg) 2) c.472_478delGCCTTCC, p.(Ala158Thrfs*6)	1) Likely pathogenic variant:- Absent in population databases (moderate evidence); - PolyPhen-2, SIFT, CADD predict a deleterious effect, conserved region (supporting evidence); - Located in the ion transport domain (moderate evidence); - Co-segregation with disease in multiple affected family members (supporting evidence); - Detected in trans with a pathogenic variant (moderate evidence). 2) Pathogenic variant:- Absent in population databases (moderate evidence); - Predicted null variant (very strong evidence); - Co-segregation with disease in multiple affected family members (supporting evidence); - Detected in trans with a pathogenic variant (moderate evidence).	Encodes the pore-forming alpha-1A subunit of the calcium channel CaV2.1, acting as an ion pore and a voltage sensor. Impairment of this calcium channel leads to synaptic dysfunction and profound neuronal loss throughout the cerebellum. Its defects are associated with early infantile epileptic encephalopathy (OMIM 617106), episodic ataxia type 2 (OMIM 108500), familial hemiplegic migrane (OMIM 141500), and spinocerebellar ataxia 6 (OMIM 183086).	Reinson et al. [26]
18	5	TSPOAP1 (BZRAP1) (NM_004758.3)	Homozygous variant:c.2449_2450delinsTG p.(Gln817*)	Likely pathogenic variant:- Absent in population databases (moderate evidence); - Predicted null variant (very strong evidence).	Encodes RIM-binding protein 1. RIMs and RIM-BPs are multidomain scaffolding proteins that bind directly or indirectly to nearly all other presynaptic active zone proteins and Ca2+ channels, and are essential to all active zone functions. So far RIMBP1 has not been associated with any Mendelian disorders, but we have described a novel AR dystonia syndrome.	Acuna et al. [48]
19	5	NKX2–1 (NM_003317.3)	Heterozygous variant:c.578A > C, p.(His193Pro)	Likely pathogenic variant:- Absent in population databases (moderate evidence); - SIFT, MutationTaster predicts a deleterious effect, conserved position (supporting evidence); - de novo (moderate evidence); - Patient's phenotype is highly specific for disease (supporting evidence).	Encodes a transcription factor that is expressed during early development of thyroid, lung, and forebrain regions, particularly the basal ganglia and hypothalamus. Defects in the protein are associated with choreoathethosis, congenital hypothyroidism with or without lung dysfunction syndrome (OMIM 610978), and benign hereditary chorea (OMIM 118700).	Thorwarth et al. [49]
22	4	ALS2 (NM_020919.3)	Compound heterozygous variants:1) c.3559_3562dup, p.(Gly1188Glufs*19), 2) c.4750C > T, p.(Gln1584*)	1) Pathogenic variant:- Absent in population databases (moderate evidence); - Predicted null variant (very strong evidence); - Detected in trans with a pathogenic variant (moderate evidence). 2) Pathogenic variant:- Absent in population databases (moderate evidence); - Predicted null variant (very strong evidence); - Detected in trans with a pathogenic variant (moderate evidence).	Encodes alsin protein, which is a member of the guanine nucleotide exchange factors for the small GTPase RAB5, and plays a role in intracellular endosomal trafficking. Its defects are associated with juvenile amyotrophic lateral sclerosis (OMIM 205100), juvenile primary lateral sclerosis (OMIM 606353), and infantile onset ascending spastic paralysis (OMIM 607225).	Hadano et al. [50]
23	4	ATRX (NM_000489.3)	Hemizygous variant:c.4466 T > A, p.(Leu1489Gln)	Likely pathogenic variant:- Absent in population databases (moderate evidence); - Align GVGD, SIFT, MutationTaster predicts a deleterious effect, conserved region (supporting evidence); - de novo (moderate evidence); - Patient's phenotype highly specific for gene (supporting evidence).	Encodes an ATP-dependent chromatin remodeling factor, which is thought to play a variety of key roles at tandem repeat sequences within the genome. Its defects are associated with X-linked alpha-thalassemia/mental retardation syndrome (OMIM 301040), mental retardation- hypotonic facies syndrome (OMIM 309580), and tumorigenesis.	Clynes et al. [51]
24	4	SCN2A (NM_021007.2)	Heterozygous variant:c.4405A>T, p.(Ile1469Phe)	Likely pathogenic variant:- Absent in population databases (moderate evidence); - Align GVGD, SIFT, MutationTaster predicts a deleterious effect, conserved region (supporting evidence); - Located in the ion transport domain (moderate evidence); - de novo (moderate evidence).	Encodes the alpha-subunit of the voltage-sensitive sodium channel NaV1.2, which is responsible for the generation and propagation of action potentials early in development in neurons, in the hippocampus and cortex. Its defects are associated with early infantile epileptic encephalopathy (OMIM 613721) and benign familial infantile seizures (OMIM 607745).	Liao et al. [52]
27	2	ABCC6 (NM_001171.5)	Compound heterozygous variants:1) c.496C > T, p.(Arg166Cys); 2) c.1567G > A, p.(Gly523Ser)	1) Pathogenic variant:- Previously reported in a patient with a typical pseudoxanthoma elasticum. 2) VUS:- ExAC frequency is 0.0015% in Europe - SIFT tolerated, MutationTaster polymorphism - Located in the ABC transporter domain (moderate evidence)	Encodes the ATP-binding cassette transmembrane transporter subfamily C6. Its defects are associated with pseudoxanthoma elasticum (OMIM 264800).	Germain [53]
28	1	MYBPC3 (NM_000256.3)	Compound heterozygous variants:1) c.1828G > A, p.(Asp610Asn) 2) c.1543_1545delAAC, p.(Asn515del)	1) VUS?- Previously reported as likely pathogenic in a patient with dilated cardiomyopathy, but later re-classified this variant as variant of unknown significance. 2) Pathogenic variant:- Previously reported in a patient with dilated cardiomyopathy.	Encodes cardiac myosin-binding protein C, which is arrayed transversely in sarcomere A-band, and binds myosin heavy chain in thick filaments and titin in elastic filaments. Phosphorylation of this protein appears to modulate contraction. Its defects are associated with dilated cardiomyopathy (OMIM 615396), hypertrophic cardiomyopathy (OMIM 115197), and left ventricular noncompaction (OMIM 615396).	Previs et al. [54]
B. mtDNA variant
3	8	MT-ATP6 (NC_012920.1)	mtDNA variant:m.9176 T > C, p.ATP6:(Leu217Pro)Heteroplasmy 96% in muscle, 95% blood, 96% urine	Pathogenic variant:- Previously reported as a disease-causing variant; - Causes a deficit (~30%) in mitochondrial ATP production in yeast.	Encodes the alpha subunit of the ATP synthase or complex V, which produces most of the ATP in human cells. Its defects are associated with neuropathy, ataxia, retinitis pigmentosa (OMIM 551500), Leigh syndrome (OMIM 256000), and mitochondrial infantile bilateral striatal necrosis (OMIM 500003).	Kucharczyk et al. [55]
C. Mitochondrial disease causing variants found on WES in patients with no clinical suspicion of a mitochondrial disorder
29	7	PDHA1 (NM_001173454.1)	Heterozygous variant:c.1014_1017dup, p.(Arg340Leufs*13)	Pathogenic variant:- Absent in population databases (moderate evidence); - Predicted null variant (very strong evidence); - de novo (moderate evidence).	Encodes pyruvate dehydrogenase (PDH) complex E1-alpha subunit, which forms the structural core of the PDH complex that catalyzes the irreversible conversion of pyruvate into acetyl-CoA in mitochondria. Its defects are associated with pyruvate dehydrogenase E1-alpha deficiency (OMIM 312170).	Brown et al. [56]
30	5	BCS1L (NM_004328.4)	Compound heterozygous variant:1) c.232A > G, p.(Ser78Gly); 2) c.245C > T, p.(Ser82Leu)	1) Pathogenic variant:- Previously reported as a disease-causing variant; - Well-established functional studies show a deleterious effect. 2) Likely pathogenic variant:- Homozygous variant absent in population databases (moderate evidence); - There was a mild decrease in RCC III activity in our patient, showing a dysfunction in its assembly, however, Rieske Fe/S protein WB was normal (supporting evidence); - Detected in trans with a pathogenic variant (moderate evidence); - Patients' phenotype highly specific for gene (supporting evidence).	Encodes a homolog of S. cerevisiae bcs1 protein, which is presumed to facilitate insertion of Rieske Fe/S protein into precursors to complex III during assembly of the respiratory chain. Its defect is the most common cause for complex III dysfunction (OMIM 124000).	Cruciat et al. [57]

ACMG – American College of Medical Genetics and Genomics, AD – autosomal dominant, adPEO – autosomal dominant progressive external ophthalmoplegia, AR – autosomal recessive, arPEO – autosomal recessive progressive external ophthalmoplegia, CMD – congenital muscular dystrophy, EDMD – Emery-Dreifuss muscular dystrophy, LGMD – limb-girdle muscular dystrophy, MDC – mitochondrial disease criteria, MNGIE – mitochondrial neurogastrointestinal encephalopathy, RCC – respiratory chain complex, RRF – ragged red fibers, VUS – variant of unknown significance, WB – Western blot.

References are of the functions of the genes and the effect of the specific gene variant where described.

Table 2

The disease-causing variants found in our study group according to the mitochondrial disease criteria (MDC) grouping.

MDC grouping [9]	Number of all cases/solved cases	Confirmed disease-causing variants
		Mitochondrial function related genes	Other genes
Definite (score 8–11)	4/4	SLC25A4, POLG, MT-ATP6	SMN1
Probable (score 5–7)	16/8	SPATA5, NDUFB11	MYH2, EPG5, TSPOAP1, NKX2–1, CACNA1A, LMNA
Possible (score 2–4)	7/4		ALS2, ATRX, ABCC6, SCN2A
Unlikely (score 1)	1/1		MYBPC3

Is located in the polymerase domain, pathogenic partitioning loop cluster () (moderate evidence);

We detected eight recurrent disease-causing variants, which have been published before, and 19 novel variants. A detailed description of the variants found (pathogenicity prediction according to the ACMG guidelines [25], function of the gene, and disease causality) is given in Table 1A-B. During our study, we have published four articles reporting interesting cases, which provided new information on the specific disease/gene. We have described a new inheritance pattern for CACNA1A associated disorders – autosomal recessive early onset epileptic encephalopathy, progressive cerebral, cerebellar and optic nerve atrophy with reduced lifespan [26]. In collaboration with the Wellcome Trust Centre for Mitochondrial Research, Institute of Neuroscience, Newcastle University, UK, we have shown that autosomal dominant SLC25A4 variants can also cause early-onset severe MD. The functional studies on SLC25A4 variants, c.239G>A, p.(Arg80His), and c.703C>G, p.(Arg235Gly) showed that recombinant mitochondrial ADP/ATP carrier isoform 1 (AAC1) mutant proteins are severely impaired in ADP/ATP transport, most likely affecting the substrate binding and mechanics of the carrier, respectively. This highly reduced capacity for transport probably affects mtDNA maintenance and in turn respiration, causing a severe energy crisis [27]. We have also carried out functional studies on rat cortical neurons to elucidate the pathomechanisms of SPATA5 gene deficiency. Our results show that SPATA5 deficient neurons have an imbalance in their mitochondrial fusion-fission rate, shorter mitochondria, reduced ATP production, and shorter axons, which means that SPATA5 defects alter the dynamics of mitochondria and impair axogenesis. Therefore, SPATA5-related diseases can be indirectly categorized under MDs [28]. In addition, we have described a novel gene-phenotype association of biallelic mutations in the BZRAP1 gene, causing autosomal recessive dystonia syndrome. The delineation of the BZRAP1-associated syndrome and the results of functional studies are shown elsewhere (manuscript in preparation).

In total, we analyzed the mtDNA of 27 patients using exome sequencing data. The coverage of mtDNA was highly dependent on the enrichment kit used, thus all samples analyzed using the SureSelect enrichment kits had markedly lower coverage. For samples sequenced with the SureSelect kits, the 1× coverage ranged from 96.5% to 99.9% of the mtDNA nucleotides, with one outlier of 38.3%, but the 20× coverage ranged only from 0 to 9%. Regarding the Nextera kits, the sequence depth of 1× was achieved in 100% of the mtDNA positions, and the 20× coverage ranged from 63% to 100%, with 12/15 samples achieving 20× coverage of >95%. In two sibs sequenced using the TruSeq kits, the 1× coverage was 100% and the 20× coverage was 51% and 76%. We detected one disease-causing variant in the MT-ATP6 gene in Patient 3, although the sample was sequenced after SureSelect enrichment, so the coverage was lower (10 reads out of 12 supported the variant). In case 11, a subsequent mtDNA depletion test from muscle biopsy showed a moderately reduced mtDNA/nDNA ratio of 40% of the average ratio observed in controls. This can be classified as a “possible” mtDNA depletion syndrome.

Discussion

To evaluate the effectiveness of WES in Estonia, we analyzed 28 patients with a suspected MD who had been referred to a clinical geneticist in the Department of Clinical Genetics, Tartu University Hospital during the years 2003–2016. As the Department of Clinical Genetics at Tartu University Hospital is the only clinical genetics department in Estonia, and our routine approach to patients with an unsolved suspected mitochondrial disease is to take a skin biopsy, the patient selection bias should be quite low. The limitations of our study are the small sample size and the fact that we only included the patients who had been referred to a clinical geneticist. In addition, it is important to point out, that in the retrospective study group, we preliminary excluded patients whose MD has formerly been confirmed by molecular genetic studies. Therefore, this may influence the final results.

Our results show that the effectiveness of WES in patients with suspected MD and an onset in early childhood is quite high – 57% (16/28). However, MDs were found in only 14% (4/28) of the patients. Nevertheless, a large proportion of the patients with clinical features suggesting an MD, had variants in other disease-causing genes, which mainly result in neuromuscular diseases or neurodegenerative disorders (Table 1A). However, it is not known to us whether the mitochondrial dysfunction in those cases is due to the primary neurological defect or a secondary finding. Our results coincide with other studies on the effectiveness of WES in MD patients [10,11].

We also observed that the higher the MDC score, the greater the probability of an MD (Table 2). All of the cases with a definite MDC scoring were confirmed as having an MD, except for the patient with spinal muscular atrophy type I (SMA I). However, the MDC scoring of the latter case was very high because the COX enzyme reaction was nearly absent in her muscle biopsy. This resulted in an MDC score of four points in total, but if we had not performed a muscle biopsy, the patient would have still scored only four points, falling into the possible MD group. Apart from the absent COX reaction, this patient had a very typical phenotype and muscle morphotype for SMA I, but with an atypical genotype (exon 7 and 8 deletion on one allele, and a novel point mutation on the other allele). A recently performed systematic review of oxidative metabolism in 24 genetically confirmed SMA patients revealed several COX-negative fibers in muscle, in all of the cases, but it was more pronounced in SMA types I and II. In addition, RCC enzyme activities were also reduced, and mtDNA depletion was detected in all the patents studied [29]. However, a total absence of COX has not been described before in SMA patients. In addition, if we had not performed a muscle biopsy in the other three patients with the definite scoring, all of them would have had a lower score, falling into the probable MD group (scores 7, 6, and 6, respectively).

In the group with a probable MD scoring, the results vary widely. Two patients were confirmed to have an MD, however, here it is important to point out that muscle biopsies were not performed in these cases. This shows that MD patients tend not to gain a definite scoring based on clinical data and metabolic/imaging studies only. Two patients had variants in myopathy/muscular dystrophy-causing genes, two in synaptic dysfunction-causing genes, one in an autophagy-related gene, and one in a transcription factor encoding gene, with the latter two causing a multisystemic disease. In the group with a possible MD scoring, no patients were confirmed to have a MD. The four confirmed cases had variants in genes with different functions and heterogeneous phenotypes.

Our results show that the higher the MDC scores, the more probable it is to find not only an MD diagnosis, but an overall genetic diagnosis by WES analysis. Interestingly, out of 501 large gene panel sequencing analyses carried out in a clinical setting, in our Department of Clinical Genetics during 2015–2016, two patients not suspected of having an MD, had a disease-causing variant in mitochondrial function related genes (Table 1C) [30]. Their MDC scores, which were calculated afterwards, were 7 and 4.

Among the patients with a confirmed MD, the clinical symptoms varied from motor neuron involvement (SMA I, MYH2- and LMNA-myopathy), and other neurological signs (synaptic dysfunction), to multisystemic involvement. However, three of the seven MD patients presented with sensorineural hearing loss (SNHL), whereas the other patients did not. However, this does not allow us to delineate SNHL as being more indicative of MDs, as our study group was quite small. Serum lactate was elevated in five of the seven MD patients, compared to just two of the 12 patients with other diseases, suggesting that lactate remains a good marker for MDs, although not 100% specific. However, alanine was elevated only in one MD patient and urine organic acid analysis was abnormal in three of the seven MD patients, whereas TCA cycle intermediates were also elevated in patients with MYH2 and LMNA variants. Brain MRI imaging was not very indicative, with only brain atrophy being observed in most of the patients, including MD patients. All these details show that we cannot clearly differentiate MDs from other diseases on clinical grounds and metabolic/imaging studies alone. Verity et al. [31] have also shown that children do not have many of the symptoms that are characteristic of MD when they first present, and the typical features may only appear later in life. Therefore, using WES as a first-tier genetic analysis for patients with neuromuscular, neurological, and/or multisystem signs and symptoms should be a standard protocol. This is especially true in children because regarding MDs, nuclear gene defects are overrepresented in pediatric cases [4,32].

Muscle biopsy is an excellent choice for detecting mtDNA mutations, deletions, and depletion [10]. Serial cryosection analysis and respiratory chain enzyme activity analysis are also ideal for guiding genetic studies, because muscle is affected in many mitochondrial disorders [4]. Indeed, muscle cryosection analysis was quite indicative in our cohort, with all five MD patients in whom a muscle biopsy was performed showed RRF, marked mitochondrial subsarcolemmal aggregates, and/or pathological mitochondrial ultrastructure, whereas the other patients did not. However, this is not always the case; according to the literature, children with MDs may not always present with RRF or abnormally structured mitochondria [33,34]. On the other hand, a finding of mitochondriopathy in muscle could mislead us, as there are many reports of genetically confirmed neuromuscular disorders with secondary mitochondriopathy [35]. Respiratory chain enzyme deficiencies can also be secondary to other various cellular processes [36]. However, the respiratory chain enzyme activities were normal in five of our seven MD patients, indicating that normal enzyme analysis does not exclude an MD diagnosis. Also, because muscle biopsy is an invasive and quite traumatizing procedure for the child and the muscle, treating physicians, as well as parents would like to postpone it whenever possible. There is currently a scientific debate on whether WES analysis should precede muscle biopsy [34,[37], [38], [39]]. Based on our current study, we suggest that WES is performed prior to the muscle biopsy; however, muscle analysis might be needed afterwards for confirmation of the variants found, or to obtain additional information about the patient's phenotype.

A less invasive procedure to study respiratory chain enzyme activity is the analysis of cultured fibroblast cells derived from a skin biopsy. However, there might be cases with decreased activity in muscle, but not in fibroblasts [34], as was the case in two of our patients, 1 and 30.

For patients with no findings with WES, the next step should be mtDNA analysis. Presently, laboratories are already working on methods to analyze nDNA and mtDNA concomitantly from blood lymphocytes using commercially available exome sequencing kits [40]. Dinwiddie et al. [41] have shown that using exome sequencing without specifically enriching for mitochondrial sequences also provides sufficient coverage for mtDNA, and both homoplasmic and heteroplasmic mitochondrial variants may be detected. We have also developed a bioinformatics pipeline for analyzing mtDNA from the data of whole exome sequencing. We experienced large differences between samples, depending on which enrichment kit was used for the WES library preparation. The SureSelect kits (Agilent) had poor coverage, whereas the Nextera and TruSeq (Illumina) had quite good coverage. However, the SureSelect enrichment kit was used in the patient with a known mtDNA mutation, which we also detected using the WES, despite the low coverage. Therefore, this is a promising technique to analyze mtDNA, without the invasive procedure of muscle or skin biopsy. However, we should bear in mind that with this technique we analyze blood lymphocytes, which may not contain the mtDNA variants, especially in patients with mitochondrial myopathy. Another bottleneck is mitochondrial heteroplasmy – a low level of the variant may not be detected [34]. In addition, mtDNA deletions and mtDNA depletion may be missed. Therefore, a negative mtDNA analysis from the WES data does not exclude the possibility of a mtDNA mutation, which can only be achieved with certainty by examining DNA from clinically affected tissue. On the other hand, if a mutation is detected by WES, this likely prevents the need for further (invasive) procedures. In addition, based on experience from case 11, the lower coverage of mtDNA from the WES data may indicate mtDNA depletion, which can be further investigated.

Conclusions

In conclusion, the effectiveness of WES in Estonia for patients with a clinical suspicion of an MD is 57%. We suggest WES analysis as a first-tier method in clinical genetic practice for children with any multisystem, neurological, and/or neuromuscular problems. As nDNA variants are more common in children with MDs, a large number of patients harbor disease-causing variants in genes other than mitochondria-related ones, and the clinical presentation might not always point towards an MD. We have also successfully conducted the analysis of mtDNA from standard WES reads, providing further evidence that this method could be routinely used in the future.

The following is the supplementary data related to this article.

Supplementary Table S1

The clinical features, results of metabolic and imaging studies and muscle pathomorphology of patients included in our study group.