INTRODUCTION

Fetal akinesia (FA) describes a clinical syndromic entity characterized by reduced or absent fetal movements leading to multiple phenotypic abnormalities. These abnormalities may include intrauterine growth restriction (IUGR), craniofacial dysmorphic features, limb pterygia, pulmonary hypoplasia, and arthrogryposis. These multiple contractures are commonly known as arthrogryposis multiplex congenita (AMC) or, when associated with pterygia, multiple pterygium syndrome. It was initially described as the Pena–Shokeir phenotype, and is further subdivided into 20 clinical subtypes that are currently linked to at least 37 different genetically defined entities.1,2 “Arthrogryposis” has been commonly used as a descriptive term for congenital contractures of the joints in at least two different body parts, in contrast to the most common isolated congenital contracture, the congenital clubfoot or pes equinovarus.3 And arthrogryposis may occur isolated or as part of a broad spectrum of syndromic conditions. The phenotypical spectrum includes both the clinical entities of AMC and the fetal akinesia deformation sequence (FADS). Therefore, we propose to use fetal akinesia (FA) as an overarching term covering the entire phenotypical spectrum from a mild AMC phenotype with a reasonable quality of life for the patient to a severe FADS phenotype with a prenatally lethal outcome.

Even though there are 166 genes linked to an arthrogrypotic or fetal akinetic phenotype (Supplementary Table S4), currently many patients with FA have no causal genetic diagnosis. Therefore counseling of affected patients and their families with regard to prognoses and recurrence risks remains challenging.2 Furthermore, this lack of knowledge is a major obstacle for the development of molecular therapies for these patients. There is an unmet clinical need to decipher the genetic basis and defective pathways in AMC and FADS.

Here, we present novel genomic insights gained from implementing next-generation sequencing (NGS) analysis in a cohort of patients presenting with FA, with a phenotypic severity ranging from mild musculoskeletal defects to prenatally lethal phenotypes.

MATERIALS AND METHODS

Inclusion criteria

Approval for the research performed in this publication has been granted by the ethics board of the University of Cologne (sign 19–1269). Written consent for publication from all included patients and/or legal guardians has been acquired and archived. The patients included in this cohort were recruited primarily with regard to the presence of multiple joint contractures manifested at or even before birth as the key diagnostic criterion for the fetal akinesia spectrum. Patients were referred to the authors when no causative pathogenic variant had been determined by candidate gene sequencing, array comparative genomic hybridization (array CGH), or karyogram. Exclusion criteria were a nongenetic etiology such as maternal antibodies or a prenatal infection. Prenatal fetal abnormalities, dysmorphic facial features, neurological defects, cognitive disabilities, limb deformations, and dysfunction involving other organs besides the musculoskeletal system were clinically phenotyped in depth if possible. Wherever accessible, the clinical findings from follow-up examinations were included in the clinical description of the patient. Clinical findings were recorded using a clinical core data form for the relevant clinical data supplied by the primary specialized health-care provider of the patient or by the fetal pathologist and the parents in cases of prenatally lethal FA. Four of the novel variants (cases 2, 11, 12, 14; see Tables 1 and 2) were previously published by the authors.4,5,6

Table 1 Clinical overview of the cohort
Full size table
Table 2 List of primary candidate pathogenic variants
Full size table

Next-generation sequencing

DNA used for NGS was extracted using a standard extraction protocol (details in Supplementary Materials file 1). The Mendeliome panel was run initially with Illumina TruSight One covering 4813 genes and later the expanded TruSight One base panel with 6709 genes used by the authors in cooperation with Illumina during the development of TSOne V2.0.7,8 Exome sequencing (ES) was performed initially using the NimbleGen SeqCap EZ Human Genome Library v2.0 and later the Agilent SureSelect V6 panel, which was chosen due to better overall performance and target coverage.9 The average coverage of all currently known FA genes has been analyzed for the different NGS kits used in this study (Supplementary Figs. S9–13). If single ES did not lead to the diagnosis, a trio exome (9 trios) was performed provided parental consent was available.

The data analysis of the generated raw NGS data was described previously7,10 using up-to-date versions of the algorithms and programs implemented within the Varbank pipeline (https://varbank.ccg.uni-koeln.de/varbank2/).

The NGS data were filtered with the aid of the Varbank pipeline to successively remove low-quality variant artifacts, deep intronic variants, variants that were predicted to be benign and variants in genes without a matching disease association that was already documented in public control databases. The remaining variants were then manually curated by gene function, known disease association, and severity of the predicted effect on the protein product to generate a list of potential candidates. Detailed information about the NGS data analysis and filtering strategy is provided in the “Methods” section of Supplementary Materials file 1 and Supplementary Fig. 14.

Dideoxy sequencing validation and segregation analysis

The validation of the variants determined by NGS was performed by dideoxy sequencing and segregation analysis was performed whenever parental DNA was available (Supplementary Materials file 2 for detailed pedigrees).

Variant classification

The remaining candidate variants were graded using both the American College of Medical Genetics and Genomics (ACMG) classification system11 as well as the proposed European Society of Human Genetics (ESHG) classification system (https://www.eshg.org/index.php?id=949). All variants achieving a pathogenic or likely pathogenic ACMG score were considered to be definitely solved; all ACMG variants of uncertain significance (VUS) achieving an ESHG rating of C–D (mildly pathogenic or susceptibility variant) were considered to be potentially solved with a variant of interest.

Systems biology analysis

The protein–protein interactions between previously published FA genes (Supplementary Table S4) and our new candidates (Supplementary Tables S1 and S2) have been visualized using Cytoscape StringApp,12 which uses STRINGdb data.

An over-representation analysis (ORA) has been performed for Gene Ontology (GO) biological process (BP), molecular function (MF), and cellular component (CC) classes via R-clusterProfiler library.13 We also used the Spearman rank order correlation to identify statistical significant differences between categorical classifications for phenotypical distribution, etiological classification, NGS success, and molecular function distribution (further details in Supplementary Materials file 1, “Methods” section).

RESULTS

Fetal akinesia phenotype classification

We recruited 51 patients from 47 families, including four affected sibling pairs, who were diagnosed with FA. Of 51 individuals in the cohort, 26 were female and 25 were male. In 21 cases the pregnancy was either terminated prematurely, the child was stillborn, or the child died shortly after birth. The age range of the other patients are 0–28 years (median age 1 year).

The cohort can be grouped into five categories based on their clinical phenotype (Table 1 and Supplementary Table S5). We propose to expand the clinical classification first introduced by Hall,3 and recently critically reviewed,14 as in our opinion the viability of the phenotype is an important feature. Category I consists of ten patients with a phenotype limited primarily to limb contractures. Category II contains patients with limb contractures and other systemic malformations and is split into IIa, indicating a viable phenotype with 14 patients, and IIb, indicating the seven patients who died either prenatally or within their first year of life. Similarly, category III contains patients with limb contractures and neurological abnormalities, and again IIIa contains the 12 patients with viable phenotypes whereas IIIb contains the eight lethal phenotypes (Fig. 1a, b).

Fig. 1
figure 1

Fetal akinesia (FA) cohort overview. (a) Phenotypical spectrum of the cohort ranging from prenatally lethal FA to a mild arthrogrypotic phenotype. We introduce a modified version of the classification first proposed by Hall.3 Category I consists of patients with a phenotype limited primarily to limb contractures. Category II contains patients with limb contractures and other systemic malformations, and is split into IIa, indicating a viable phenotype, and IIb, indicating patients who died either prenatally or within their first year of life. Similarly, category III contains patients with limb contractures and neurological abnormalities, and again IIIa contains the viable phenotypes whereas IIIb contains the lethal phenotypes. (b) The cohort has been sorted into the five phenotypical categories based on the modified Hall classification described above. (c) The cohort has been sorted into the following categories with regard to the success of next-generation sequencing (NGS) in identifying potential pathogenic variants: group 1—identified variants or copy-number variants (CNVs) in known AMC/FADS genes; group 2—identified variants in known disease genes; group 3—identified putative variants in genes previously not linked to a disease; group 4—currently unsolved cases. (d) Etiological distribution grouped into primarily neurogenic or myogenic pathomechanisms, pathomechanisms involving the neuromuscular junction, limb malformation and syndromic malformation, as well as unclassified in three cases where no prior disease association is known. (e) Breakdown of the molecular function of the genes carrying pathogenic variants. (f) This network map shows the gene interactions generated from STRINGdb. Gray nodes represent the previously published FA-related genes that are not present in our study; yellow nodes represent the previously published FA genes that are also found in our study; red nodes represent the new primary and secondary candidate genes that have been found in our study. Size differences illustrate the variant counts for a particular gene in this study. Blue lines between nodes represent the combined confidence score of interaction between 0.4 and 0.9 (medium confidence–high confidence). Green lines represent the combined confidence score above 0.9 (highest confidence).

Full size image

Molecular findings

The pathogenic variants in 45 of 51 cases can be divided into the following three categories (Fig. 1c):

In 27 cases we found the underlying pathogenic variants in known FA-associated genes, including a sib pair with a homozygous deletion in RAPSN (cases 16 and 17; Supplementary Fig. S15b). Two of the variants in ACTA1 (case 1) and PIEZO2 (case 15) respectively were already known according to ClinVar (ClinVar IDs 381639 and 137629). The remaining 31 variants are to our knowledge novel pathogenic variants.

In a further 15 cases, the underlying pathogenic variants are in known disease genes that until now had not been linked to FA. These novel FA candidate genes containing 18 pathogenic variants are as follows: ADSSL1, ASAH1, ASPM, ATP2B3, EARS2, FBLN1, PRG4, PRICKLE1, ROR2, SETBP1, SCN5A, SCN8A, and ZEB2, as well TNNT1 with a homozygous copy-number variant (CNV) deletion in a sibling pair (cases 40 and 41; Supplementary Fig. S15a). The pathogenic variant in SETBP1 (case 37) is a known pathogenic variant (ClinVar ID 1031) for Schinzel–Giedion syndrome (OMIM 269150), which fits the patient's phenotype. A homozygous variant in PRG4 (case 34) and the homozygous CNV deletion in TNNT1 mark the first time these genes have been linked to FA beyond in silico predictions.15

Furthermore, in three cases we have found likely pathogenic variants in putative novel FA candidate genes without clear and convincing prior FA-disease link until now: GCN1, IQSEC3, and RYR3. A de novo heterozygous variant in GCN1 and compound heterozygous variants in VPS13D were found in case 43. While GCN1 appears to play an important role in translational regulation during stress response, the scientific literature provides no clear disease association beyond a tenuous link to schizophrenia.16,17 VPS13D, on the other hand, is involved in mitochondrial regulation and has been linked to ataxia and intrauterine fetal death.6,18 However, the variant in GCN1 is a confirmed de novo variant that is absent from controls with a Combined Annotation Dependent Depletion (CADD) score of 22, the gene has a Residual Variation Intolerance Score (RVIS) of −3.77, is expressed in skeletal muscle and cerebellum brain tissue according to GTEx, and according to the Gene Network 2.019 is coregulated with the FA gene DYNC1H1 (p value 6.7 × 10−42), whereas one of the compound heterozygous variants in VPS13D has a rather high allele frequency. Therefore, we consider the variant in GCN1 as the primary candidate while the compound heterozygous variants in VPS13D might be modifiers or contribute to an oligogenic phenotype. The VPS13D variants are therefore listed as secondary variants. In case 44 we found compound heterozygous variants in IQSEC3, a gene known to play a role in the regulation of synapse formation20 and recently linked to a neurodevelopmental phenotype in two cases.21 One of the two variants causes an amino acid change, is absent from controls, and reaches a CADD score of 26.5, while the other one is predicted to trigger a cryptic splice site activation. According to GTex, IQSEC3 is almost exclusively expressed in the brain, and particularly in the cerebellum. It has a RVIS score of −0.84 and is also coregulated with another novel FA candidate gene, ATP2B3 (p value 8.8 × 10−21), according to the Gene Network 2.0. Therefore, we conclude that given its function in neuronal development IQSEC3 is a likely candidate as a novel FA gene. A highly conserved heterozygous variant in RYR3 has been found in case 45 with a high CADD score of 32. RYR3 is a Ca2+ ion channel for excitation–contraction coupling with a RVIS score of −5.87 and closely related to the known FA gene RYR1. A stringApp analysis of the interactome showed that it interacts closely with the FA genes RYR1 and STAC3 (see Supplementary Fig. S25). It has recently been linked to nemaline myopathy22 as well as a phenotype with global developmental delay.23 Given its expression in both brain and muscle tissue and the murine animal knockout model with a severe neuromuscular phenotype it is a promising novel FA disease gene candidate.24

This translates to a 73% (37/51) success rate of solved cases and in total 88% (45/51) of possibly solved cases including potentially pathogenic variants linked to the observed phenotype. Detailed overviews of the clinical features of each patient and of the proposed causative pathogenic variants can be found in Tables 1 and 2, respectively.

Disease association of candidate genes

With regard to known disease associations and anatomical–physiological locations of the primary defect, the disease candidate genes were grouped (Fig. 1d) as follows. Eleven genes are known to cause neurogenic disease phenotypes: ASAH1, ASPM, ATP2B3, CNTNAP1, EARS2, GLDN, LGI4, NALCN, PRICKLE1, SCN8A, and ZEB2. Another seven genes are linked to myogenic phenotypes: ACTA1, ADSSL1, GBE1, RYR1, SCN4A, SCN5A, and TNNT1. The third group includes four genes that are associated with dysfunction of the neuromuscular synaptic junction: CHRND, CHRNG, RAPSN, and UNC50. Additionally, ASCC1, PIEZO2, and SETBP1 are linked to complex developmental syndromic malformations including the nervous system. Another set of genes are known to cause limb malformations—this group includes FBLN, PRG4, and ROR2. For detailed information regarding these disease associations please refer to the PubMed ID (PMID) links provided in Supplementary Table S1. Recently, the first disease associations have been reported for IQSEC321 and RYR3.22,23,25 GCN1 had been associated with schizophrenia in genome-wide association studies (GWAS).16 All primary pathogenic variant candidates are listed in Table 2, and Supplementary Table S1 contains additional information for the primary pathogenic variant candidates. In addition to the candidate genes listed above, secondary variants besides the primary candidates survived our filtering strategy in several cases: ALDH5A1, DQX1, DYNC1H1, GFRA4, HKR1, KIAA1109, MAGI3, NAGA, PIEZO2, SPAG16, TMPO, and VPS13D (Supplementary Table S2). They might function as disease phenotype modifiers or form part of an oligogenic phenotype in the case of likely pathogenic variants in ALDH5A1, KIAA1109, or PIEZO2.

Gene functions

The potential candidate genes discovered in our cohort cover a broad spectrum of molecular functions based on the gene function descriptors of OMIM (Fig. 1e) as expected in the case of phenotypes such as AMC and FADS with their heterogeneous genetic etiology. The most common group among our candidates are the ion channel or ion pump genes such as ATP2B3, CHRNG, CHRND, NALCN, PIEZO2, RYR1, SCN4A, SCN5A, SCN8A. Functionally related are the receptor protein genes such as ROR2, ion channel modulators like RAPSN, as well as myelinization modulators like CNTNAP1 and LGI4. Another major group contains genes involved in the regulation of transcription and translation such as ASCC1, GCN1, EARS2, PRICKLE1, SETBP1, and ZEB2 as well as UNC50, a gene involved in protein trafficking. Motor protein genes like ACTA1, DYNC1H1, or TNNT1 are also potential disease candidates for some of our patients, while the extracellular matrix structure gene FBLN1 and the proteoglycan gene PRG4 contain other potentially pathogenic variants. In other cases, pathogenic variants are found in genes involved in cell cycle regulation (ASPM) or cell development (GLDN). The last major group of disease candidate genes contains several genes involved in cytosolic, mitochondrial, or lysosomal metabolic function such as ADSSL1, ASAH1, GBE1, and IQSEC3. Figure 1f illustrates the protein–protein interaction network based on STRINGdb and contains all proteins linked to an FA phenotype as well as the primary and secondary candidates found in this study.

Correlation analysis

A correlation analysis between the phenotype of the patients, the etiology or the pathomechanism of the disease-associated genes, and their molecular function (based on the data collated in Supplementary Table S5) revealed (Supplementary Fig. S22) that the more severe phenotypes (IIb and IIIb) could all be solved whereas the milder phenotypes (I, IIa, and IIIa) contain several cases without any promising candidate genes that could be discovered via NGS. While the molecular function of genes where the pathomechanism leads to primary muscle diseases is concentrated on motor proteins, ion channels, and cell metabolism proteins, the genes with a neurogenic pathomechanism appear more diverse in their molecular function.

DISCUSSION

The aim of this study was to determine the underlying genetic cause for FA by NGS, a genetically extremely heterogeneous syndromic disease entity. This heterogeneity comprises hundreds of genes containing causative pathological variants (Supplementary Table S4).

The results from this study suggest a major role for pathogenic variants in ion channel genes and genes coding for ion channel modulators in the pathogenesis of FA (Fig. 1e). RYR1 alone contained the causative pathogenic variants in five unrelated cases and one sibling pair of a total of 37 solved cases. Our study proves the importance of pathogenic variants in RYR1 as a potential cause for FA in a large cohort together with an earlier description by Alkhunaizi.26

Some of our findings are particularly noteworthy. TNNT1 was first linked to human diseases by Johnston et al.,27 who proved its link to a subtype of nemaline myopathy (OMIM 605355) occurring in old order Amish families. The first case reported outside the old order Amish population by van der Pol et al.28 shares remarkably similar clinical symptoms such as pectus carinatum, kyphoscoliosis, multiple joint contractures, respiratory insufficiency, and dysphagia with our patients (40 and 41) from one afflicted family of Eastern European (Romania/the former Yugoslavia) descent (Fig. 2). The reported results of the muscle biopsy generally mirror our findings with highly divergent muscle fiber calibers and increased endo- and perimysial connective tissue. Our findings differ in that we see trilaminary fibers (Fig. 2r) that have not yet been reported in TNNT1 cases whereas no classical nemaline rods were present in the examined tissue samples.29 To our knowledge, this is the first time that a pathogenic gene variant has been discovered for myopathy with trilaminary fibers.30 While in all previous cases the causative pathogenic variant was homozygous or there were compound heterozygous small mutations, here we report a homozygous CNV deletion in the two affected siblings as a likely cause for the phenotype, as reported very recently in another case with different myopathic phenotype.31

Fig. 2
figure 2

Clinical presentation of cases 26, 30, and 40. (a) Fetal akinesia (FA) patient 26 with a de novo SCN4A variant with respiratory support at 24 months. (b) Close-up view of the right hand of patient 26 with an ulnar deviation of the hand at 24 months. (c) Close-up view of the right foot of patient 26 with characteristic joint contractures at 24 months. (d) Patient 26 at 2 months with prominent joint contractures. (e) Close-up view of the facial features of patient 26 at 12 months. (f) Patient 26 shortly after birth with prominent joint contractures. (g) Congenital contractures of the left foot observed shortly after the birth of patient 26. (h) NADH stain of muscle tissue biopsy of patient 26 showing darkly stained type 1 fibers, myopathic presentation. (i) Slow myosin immunostain of type 1 muscle fibers (dark) of patient 26 showing an increased variation of fiber size. No fiber grouping is observable, myopathic presentation. (j) Patient 30: aborted fetus in the 18th gestational week with microcephaly and arthrogryposis due to compound heterozygous ASPM pathogenic variants. (k) Lateral view left/right of the microcephalic head of patient 30. (l) Lateral view left/right and axial view of the fetal brain of patient 30 with hypoplastic forebrain and lack of gyration. (m) Histological overview of a fetal brain tissue slice of patient 30. (n) Increased magnification of fetal brain tissue of patient 30 with reduced size of the cortical plate (C) and neuroblast stretch between cortex and lateral ventricle (V) resembling subcortical band heterotopia. (o) Patient 40 with prominent pectus carinatum, facial features, and multiple joint contractures including jaw contractures at the age of 2 due to a homozygous copy-number variant (CNV) deletion in TNNT1. (p) Hematoxylin and eosin (H&E) stain of muscle tissue biopsy of patient 40 with several hypertrophic and numerous small atrophic fibers and pathological caliber variance, myopathic presentation. (q) Lateral view of patient 40 at the age of 2 years with flexed hips due to contractures. (r) NADH stain of muscle tissue biopsy of patient 40 with trilaminary fibers, in the upper right corner, is an enhanced picture detail showing a trilaminary fiber.

Full size image

ADSSL1 has only recently been discovered as a myopathy gene.32,33 Here we report the first case of a pathogenic homozygous frameshift ADSSL1 variant in a patient of Turkish descent (patient 28). Previously, ADSSL1 pathogenic variants (OMIM 617030) have been reported in the Korean population. Unlike the eight Korean patients, our unrelated patient presented with congenital joint contractures and a more severe neurological phenotype, while all previously described patients with recessive ADSSL1 pathogenic variants suffered from adolescent-onset myopathy and did not show any contractures. A possible explanation for the earlier onset and more severe phenotype could be the fact that all previously described mutations resulted in amino acid exchanges except for one compound heterozygous frameshift variant. This frameshift variant truncates the sequence of the protein after a substrate binding domain, unlike the frameshift variant we present here. Another factor could be the presence of a frameshift variant in ALDH5A1, which might modify the phenotype of the patient as a multilocus phenotype.34

Our study emphasizes that neuronal migration defects may lead to fetal akinesia. We present the histopathological examination of the brain of a fetus (patient 30) that suffered from an ASPM nonsense variant at—to our knowledge—the earliest developmental stage published (Fig. 2). The neuronal migration was delayed and resulted in a neuroblast layer in the subventricular cortical zone still present at gestational week 18 and consequently a hypoplastic cortex—in particular in the frontal cortex (Fig. 2).35 The results of Passemard et al.36 agree with our findings that ASPM defects in humans primarily affect cortical development, which in turn leads to reduced intrauterine movement and consequently an FA phenotype.

With a success rate of about 88%, we could determine the underlying genetic defect for most of the patients in the cohort presented here. The initially employed targeted panel approach using the Mendeliome proved unsatisfactory because only 41% (21 cases) of all cases could be solved this way. For the remaining initially unsolved 30 cases, the ES analysis had a success rate of 53% (16 cases). This includes the four cases that could only be solved by systematic CNV analysis, which significantly increased the power of pathogenic variant detection in ES data sets.

Both of the currently employed ES kits provided a satisfactory average coverage of 20× or above and allowed the analysis of all 166 known FA genes, whereas the initially used Mendeliome enrichment kit did not sufficiently cover 45 known FADS disease genes. Even the improved version of the Mendeliome enrichment kit failed to cover 20 genes.

Another issue of recent interest due to the increasing availability of ES data is multilocus or oligogenic pathogenic variations and their combined effect on the phenotype of individual patients. The work of Posey et al.34 and Karaca et al.37 indicated that pathogenic multilocus variants are a valid alternative to classical phenotypic expansion of a known disease gene when it comes to explaining divergent phenotypes in consanguineous families. Here, we demonstrate the validity of this hypothesis also in nonconsanguineous cases. The potential multilocus candidates that survived the filtering process are provided in Supplementary Table S2.

The genes containing presumed pathogenic variants in our study vary in both molecular function and disease phenotype association, but most can be broadly grouped into five groups. The two major groups contain myogenic and neurogenic genes, while the three smaller groups consist of genes linked to neuromuscular junction dysfunction, limb malformations, and syndromic malformations involving multiple organ systems (Fig. 1d, Supplementary Table S5).

A total of 192 genes were included in the systems biology approach. One hundred fifty of those were previously published with an FA association but not found in our study; 16 were found in both gene lists and 26 were found only in our study (Supplementary Tables S1, S2, S4). The STRINGdb analysis (Fig. 1f) performed for this study showed that most of the novel genes proposed to be involved in the development of an FA phenotype are interacting directly with genes known to cause AMC/FADS. This supports the idea that pathogenic variants in these novel genes adversely affect the same pathways as pathogenic variants in known FA disease genes. An example of this would be SCN8A, which interacts with SCN4A and CNTNAP1, or ASPM, which interacts with KIF14, CENPJ, and CEP55. Moreover, even within the group of the known FA disease genes, there are multiple outliers with no direct interaction with any other known FA disease gene, indicating that there might be gaps in our knowledge regarding the protein networks and metabolic pathways involved in the FA phenotypical spectrum. This of course means that even if no direct interaction with a known disease gene can be found via this method, a novel gene carrying a pathogenic variant should not be excluded based solely on this lack of interaction partners, as it might belong to the outlier group with no direct interaction with the large majority of known FA disease genes.

The Gene Ontology enrichment analysis (Supplementary Figs. S16–21) showed that our novel candidate genes reveal a similar enrichment pattern as the known FA genes used as a comparison group with regard to the biological process affected, the molecular function of the protein, and cellular compartment localization, favouring ion channels, primary muscle/myogenic and excitation-coupling genes/pathways. The notable outcome was that the over-representation analysis (ORA) performed with known FA genes had a higher statistical significance compared with the ORA performed with our candidates for all GO classes (known FA genes GO:BP ORA p > 2 × 10−6; all candidates GO:BP ORA p > 0.01; known FA genes GO:CC ORA p > 0.001; all candidates GO:CC ORA p > 0.025). This was to be expected given the difference in the number of genes used in the analysis. It was also supported by the fact that the enriched genes are often common to both gene lists (Supplementary Figs. S16–21).

Overall, our findings support the usage of exome sequencing and stringent variant filtering combined with thorough clinical assessment of patients together with refiltering as an essential tool to uncover the genetic etiology of complex syndromic diseases. Applying a selected panel of disease genes—the Mendeliome—can certainly be helpful in identifying novel pathogenic variants in known disease genes as we solved 41% of all cases this way. However, by analyzing not only a particular set of enriched genes as in a gene panel but examining the entirety of the patients’ exome with a superior target coverage previously unnoticed pathogenic variants even in known disease genes can be detected. In addition, the effect of multiple pathogenic variants on an oligogenic phenotype might become evident, and CNVs with pathogenic effect are observable, as we have shown here (Supplementary Fig. S15).

The methodical analysis presented in this publication and subsequent efforts will hopefully enable clinicians to determine the underlying genetic cause of many FA cases that until now have been labeled as sporadic and allow for improved genetic counseling and prenatal diagnosis. This makes ES a tool that is most suited to uncover potential issues for couples planning further pregnancies, in particular given the recent developments of rapid exome sequencing and ultrarapid genome sequencing.38