Expanding the Boundaries of RNA Sequencing as a Diagnostic Tool for Rare Mendelian Disease

Abstract

Gene-panel and whole-exome analyses are now standard methodologies for mutation detection in Mendelian disease. However, the diagnostic yield achieved is at best 50%, leaving the genetic basis for disease unsolved in many individuals. New approaches are thus needed to narrow the diagnostic gap. Whole-genome sequencing is one potential strategy, but it currently has variant-interpretation challenges, particularly for non-coding changes. In this study we focus on transcriptome analysis, specifically total RNA sequencing (RNA-seq), by using monogenetic neuromuscular disorders as proof of principle. We examined a cohort of 25 exome and/or panel "negative" cases and provided genetic resolution in 36% (9/25). Causative mutations were identified in coding and non-coding exons, as well as in intronic regions, and the mutational pathomechanisms included transcriptional repression, exon skipping, and intron inclusion. We address a key barrier of transcriptome-based diagnostics: the need for source material with disease-representative expression patterns. We establish that blood-based RNA-seq is not adequate for neuromuscular diagnostics, whereas myotubes generated by transdifferentiation from an individual's fibroblasts accurately reflect the muscle transcriptome and faithfully reveal disease-causing mutations. Our work confirms that RNA-seq can greatly improve diagnostic yield in genetically unresolved cases of Mendelian disease, defines strengths and challenges of the technology, and demonstrates the suitability of cell models for RNA-based diagnostics. Our data set the stage for development of RNA-seq as a powerful clinical diagnostic tool that can be applied to the large population of individuals with undiagnosed, rare diseases and provide a framework for establishing minimally invasive strategies for doing so.

Keywords: Mendelian disease; RNA-seq; diagnostics; muscular dystrophy; myotubes; transcriptomics; transdifferentiation.

Figure 1

Overview of the Diagnostic Pipeline 70 samples were processed through our pipeline. Total RNA was extracted from muscle biopsies, fibroblasts, and t-myotubes, was polyA selected, and then sequenced at a depth of 50–100 paired-end reads. Our RNA-seq diagnostic algorithm compares undiagnosed individuals with our in-house database and with control transcriptome data obtained from GTEx. We worked first from a panel of genes known to be mutated in neuromuscular disorders (n = 132), and we focused our analysis in parallel on (1) novel splicing events (far left), (2) imbalances in allelic expression (middle left), (3) statistically significant differences in expression (middle right), and (4) rare sequence variants of clinical relevance (far right). Using this strategy, we were able to solve 36% of cases that were “negative” according to the gene panel and/or exome sequencing (n = 25 total unknown cases).

Figure 2

Sample Distribution and Gene-Expression Profile (A) A multi-dimensional scaling (MDS) plot of our cohort of 70 samples (primary fibroblasts, t-myotubes, and muscle) compared with tissue-matched sets from GTEx (blood, transformed fibroblasts, and muscle). There is significant overlap between the muscle samples from our cohort and GTEx. Samples from blood, muscle, and fibroblasts formed distinct clusters, but transdifferentiated myotubes cluster as a group between the fibroblasts and muscle samples. (B) A principal-component analysis (PCA) of only our muscle samples shows a clear clustering of samples on the basis of age (x axis; 1st PCA). There is increased variability in expression in the younger samples; this variability diminishes with increased age of the samples. (C) We identified 9,932 genes expressed at >1 RPKM. Of these genes, ∼92% were expressed at ≥1 RPKM in t-myotubes, 88% in fibroblasts, and 71% in blood. (D) Comparison of expression between our muscle samples and GTEx muscle for the five highest- and lowest-coverage genes. There is no significant difference between our samples and those from GTEx.

Figure 3

Analysis of Low-Frequency Novel Junctions (A) 177 novel junctions were detected in the transcriptomes of our muscle samples, and we divided them into high-frequency junctions (HFJs) and low-frequency novel junctions (LFNJs). The LFNJs were subdivided into new acceptor (n = 63), new donor (n = 18), and exon-skipping (n = 3) events. Pseudo exons (n = 4) are counted within new acceptor or donor events. (B) Example of a novel donor site detected in t-myotubes and muscle biopsy (family 4). The canonical donor splice site of exon 21 of RYR1 was reduced (551 reads), and two alternative, cryptic donor splice sites were activated. The first (chr19: 38,954,212) (new donor with 107 reads) extends exon 21 by 44 bp, causing a +1 frameshift and stop-gain at chr19: 38,955,338341. The second (chr19: 38,954,309) (new donor with 140 reads), extends exon 21 by 141 bp, causing a 1 nt frameshift with stop-gain at chr19: 38,954,395397. Overall, RYR1 expression is decreased 2.1-fold in affected cells compared to in controls. (C) Example of a pseudo exon (family 38). A novel exon was found in the muscle between exons 50 and 51 in DYSF (supported by 107 reads), creating a premature stop codon at chr2: 71,900,460–71,900,462. (D) Example of exon skipping (family 15). An exon-skipping event was detected in fibroblasts in LAMA2 (supported by four reads), causing a frameshift with stop-gain in exon 35 at chr6: 129,704,311–313. Overall, LAMA2 expression was decreased (see also Figure 4C).

Figure 4

Analysis of Variants, Expression Profiles, and Allele-Specific Expression (A) Detection of a 5′ UTR mutation with RNA-seq. WES for family 5 reported a heterozygous missense variant in GMPPB c.94C>T (p.Pro32Ser). We identified this variant in RNA from the fibroblasts (left arrow; labeled “Exon1 mutation”), along with a heterozygous 5′ UTR variant, g.349761246G>A (right arrow; labeled “5′UTR mutation”), that results in a new start codon. This variant was seen only in the transcriptome data and not by WES. It produces an in-frame insertion of 116 bp, potentially adding 29 new amino acids to the protein. Of note, both variants were absent from control fibroblasts. (B) The left panel depicts the n-fold difference in DMD RNA levels between GTEx control muscle and 4 cases of Duchenne muscular dystrophy (in family 9 = intronic splice gain, family 26 = mutation still elusive, family 34 = intronic splice gain, and family 35 = disruption of exon 1 ATG). The right panel shows the top 3 most downregulated genes (from a panel of 132 neuromuscular genes) for each of these cases. (C) Transcript expression ([C1], in RPKM) of LAMA2 from fibroblasts (in orange) from an individual with congenital muscular dystrophy and from the parents (father in green and mother in blue) (family 15). The proband had a greater-than-2-fold reduction in transcript expression as compared to his mother. (C2) He also had a maternally inherited pathogenic missense variant in LAMA2 (c.2548T>G [p.(Cys862Arg)]). This variant demonstrated mono-allelic expression, and it was found in 85% of the reads in the proband, 50% of the reads in the mother, and not found in the father.

Figure 5

Most Common Muscle-Disease-Causing Genes, Classified by Their Main Function, That Are Expressed at >1 RPKM We organized 123 of the 132 genes in our panel into seven distinct groups associated with muscle function, and we determined their expression in four different tissues (skeletal muscle, t-myotubes, fibroblasts, and blood). We calculated the percentage of genes in each category that were expressed at ≥1 RPKM in each of the tissues, and we saw that most of the muscle genes are poorly expressed in blood and fibroblasts and sufficiently expressed in t-myotubes (adapted from Dowling et al., 2017⁷⁹).