Fluorescent-activated cell sorting of muscle cells and RNA seq. myo-3::gfp-expressing body wall muscle cells from synchronized early L1 larvae were isolated by FACS as described in Material and Methods and in [41, 42] (Figure 1). Muscles were isolated from the wild-type (N2) strain and etr-1(lq61) mutants. Three biological replicates for each genotype were isolated. RNA was also isolated from triplicate samples of non-dissociated L1 larvae for the all-cell control group.
A total of twelve RNA-seq libraries were constructed, representing three biological replicates of both N2 and etr-1(lq61) with both sorted muscle cells and whole L1 larval stage cells (see Materials and Methods). Paired-end 150-bp reads were generated from each of the twelve samples using the Illumina Nextseq550 platform FASTQ files can be accessed in the Sequence Read Archive, Project number PRJNA733501
A muscle cell transcriptome defined by RNA-seq.
We used DEseq2 to identify genes with significant differential expression in wild-type muscle cells compared to all L1 stage wild-type cells (Supplemental File 1). There were 3,718 protein-coding genes with significantly higher expression in muscle cells compared to all cells (log2-fold change ≥ 0.5849 (1.5x); q ≤ 0.05), including many
canonical muscle structure and function genes previously shown to be expressed in muscle (e.g. unc-15/paramyosin, unc-54/myosin, unc-95/paxillin, and the myofilament structure pat genes [43]. DEseq2 also identified etr-1 as being more highly expressed in muscle (Supplemental File 1; log2fold change 2.21; q = 9.93e−37). Differentiall exon usage using DEXseq showed that multiple etr-1 exons showed significantly increased expression in muscle compared to all cells as a whole (Figure 2A). The etr-1 locus is extensively alternatively spliced [20, 44], and isoforms with exon 8 are required in muscles for Q neuroblast descendant migration [20]. Exon 8 expression was increased in muscles compared to all cells as a whole (Figure 2A), suggesting that isoforms with exon 8 might be more abundant in muscles compared to all cells as a whole. Expression of 8,763 protein coding genes was significantly reduced in muscle cells compared to all cells (Supplemental File 1). Differential expression of non-coding RNAs and pseudogenes is presented in Supplemental File 1.
Genes with exon usage affected by etr-1(lq61) encode molecules involved in myofilament lattice structure and attachment, and muscle physiology.
CELF family proteins are known to regulate splicing [12, 45], and we endeavored to determine the effects of ETR-1/CELF on the muscle transcriptome, including splicing. We compared exon representation across the genome in wild-type and etr-1(lq61) mutant muscle cells using the Bioconductor package DEXseq [46] (see Materials and Methods). Across the genome, there were 242 protein-coding genes and seven non-coding RNA genes with at least one exon significantly differentially represented in etr-1(lq61) muscle compared to wild-type (q ≤ 0.05) (Supplemental File 2).
The etr-1(lq61) mutation is a premature stop codon in alternatively-spliced exon 8 [20]. Exon 8 was significantly underrepresented in etr-1(lq61) muscle compared to wild-type (Figure 2B), suggesting that transcripts containing exon 8 are reduced in etr-1(lq61) muscle cells, as predicted. Exon 13 was also significantly underrepresented in etr-1(lq61). Possibly, transcripts with exon 8 might preferentially
contain exon 13. Alternatively, ETR-1 containing exon 8 might be involved in the regulation of processing of etr-1 exon 13.
We used the Database for Annotation, Visualization and Integrated Discovery (DAVID) [47, 48] to perform a gene ontology term (GO term) analysis on this gene set that showed differential exon representation in etr-1(lq61) mutant muscle compared to wild-type muscle, including both over-and underrepresented exons in wild-type compared to etr-1(lq61) (see Materials and Methods). We analyzed GO terms for the three categories: biological process (BP), cellular component (CC) and molecular function (MF) (Supplemental File 3). The six most significantly enriched GO terms in each category are shown in Figure 3. These include GO terms associated with myofilament lattice formation and function (e.g. striated muscle myosin thick filament assembly, locomotion, M band, striated muscle thin filament assembly, striated muscle dense body, I band, actin filament, and actin filament binding). Also included are muscle physiology GO terms (e.g. voltage gated ion channel activity, calcium ion binding, voltage gated potassium channel activity, and kinase activity). These are consistent with the previously-reported effects of etr-1 RNAi knockdown on muscle development and attachment [19]. GO terms of apoptotic process and reproduction are also enriched, consistent with a known role of etr-1 in germline development and engulfment of germ cell apoptotic corpses [44].
Genes with transcript accumulation affected by etr-1(lq61) encode molecules involved in translation and ribosome function.
The CELF-family proteins control transcript stability [12, 13, 45]. We used stringTie and DEseq2 to identify transcripts with differential accumulation in wild-type versus etr-1(lq61) muscle cells (see Materials and Methods). We identified transcripts that were differentially represented with a log2 fold change ≥ 1 (2x) and a false discovery rate less than 0.05 (q ≤ 0.05), to increase stringency given the large number of genes returned in this analysis (Supplemental File 4).
We identified 1180 transcripts representing 971 loci with altered accumulation in etr-1(lq61) (Figure 4A), including coding and non-coding RNAs (Supplemental File 4). 506 loci had transcripts that were overrepresented, and 414 loci had transcripts that were underrepresented in etr-1(lq61) muscle cells (Figure 4B). There were 51 loci in which some transcripts were overrepresented and some underrepresented (Figure 4B and Table 1).
Gene Ontology enrichment analysis was conducted on four separate groups of these muscle genes with transcripts affected by etr-1(lq61) (Supplemental File 5): all genes; genes with transcripts that were only underrepresented; genes transcripts that were only overrepresented; and genes that had transcripts both over- and underrepresented. The six most significant GO terms for each group are shown in Figure 5. Considering all genes, GO terms associated with translation and ribosomal function were apparent, as well as myofilament structure, muscle physiology, reproduction, and embryonic and larval growth (Figure 5A). Genes with transcripts underrepresented in etr-1(lq61) were described by GO terms representing translation and ribosomal function (10 of the 18 top GO terms) (Figure 5B). Genes with overrepresented transcripts were described by GO terms representing a broad cross section of cellular function, but translation and the ribosome were not among these (Figure 5C). Genes with both over- and underrepresented transcripts were described by GO terms representing myofilament lattice and muscle physiology, and other cellular functions (Figure 5D). In sum, this GO term analysis suggests that etr-1(lq61) influences a broad spectrum of cellular events in muscle, including myofilament lattice and muscle physiology, as well as translation and ribosomal function. Notably, genes involved in translation and ribosomal function are strongly represented among those with transcripts reduced in etr-1(lq61). A similar reduction of expression of genes involved in translation and ribosomal function was described after siRNA knock-down of CELF1 in chicken cardiomyocytes [17].
Genes identified by both differential exon usage and transcript accumulation
etr-1(lq61) affected transcript expression of 971 genes, and alternative exon representation of 244 genes. There were 102 genes shared between the 244
alternatively-spliced genes and the 971 genes that had differentially-expressed transcripts (Table 3), a significant association (p < 0.0001). These might represent genes with transcripts for which ETR-1 controls both splicing and transcript accumulation. Alternately, the differential use of exons could influence the transcript accumulation DEseq2 algorithm, leading to under- or overrepresentation of transcripts by alternate exon usage. In any event, identification in both analyses suggests that ETR-1 might have strong effects on the transcripts of these genes.
Genes affecting AQR and PQR migration
This analysis suggests that ETR-1 regulates multiple aspects of muscle cell function, most notably myofilament lattice, muscle physiology, and translation and ribosomal function. Previous studies indicated that etr-1(lq61) had a muscle-derived, non-autonomous effect on migration of AQR and PQR neurons [20]. Thus, ETR-1 might regulate a secreted signal from the muscles that directs AQR and PQR migration. Genetic analysis suggests this signal could act with or in parallel to Wnt signaling [20], which directs AQR and PQR migration [38].
We used feeding RNAi (see Materials and Methods) and mutants to knock down a subset of genes with known roles in cell migration, and genes that had transcripts both overrepresented and underrepresented in etr-1(lq61) muscles versus wild-type muscles (Table 1). AQR and PQR position in these animals was scored (see Materials and Methods) (Table 3). In unc-52, unc-71, and skn-1, between 2% and 5% of AQR and PQR neurons migrated in the wrong direction, as evidenced by AQR in position 5, and PQR in positions 1, 2, and 3. UNC-52 is the basement membrane heparan sulfate proteoglycan Perlecan [49] and will be discussed in more detail below. UNC-71 is an ADAM metalloprotease that has been shown to act in anterior Q descendant migration [50, 51]. SKN-1 is an ortholog of the human NFE2L1 transcriptional regulator that controls a wide variety of developmental events including muscle differentiation [52, 53]. Some genes displayed defects in the ability of AQR and PQR to migrate (>10%), but not directional defects, as evidenced by AQR in positions 2, 3, and 4 and PQR in position 4. Many genes displayed few (<10%) or no defects in the ability of AQR and PQR to migrate. Thus, some genes with transcripts that are regulated by ETR-1 in muscle had instructional roles in directing AQR and PQR migration, and some had permissive roles in the ability of AQR and PQR to migrate.
The heparan sulfate proteoglycan UNC-52 affects AQR and PQR migration
unc-52 encodes the basement membrane heparan sulfate proteoglycan Perlecan and is involved in myofilament lattice attachment to the basement membrane [49]. unc-52 is extensively alternatively spliced [54, 55], including in the epidermis by CCAR-1 [56] and MEC-8 [57, 58] with consequences on hemidesmosome formation, muscle attachment and mechanosensory neuron function.
We found that etr-1(lq61) affected unc-52 transcript expression in muscle cells. unc-52 was identified in both exon representation by DEXseq and transcript accumulation by DEseq2 (Table 2). In etr-1(lq61) muscle cells, some 5’ exons were significantly overrepresented and 3’ exons significantly underrepresented compared to wild-type muscle cells. The far 5’ exons (see Figure 6A) predicted for unc-52 were not highly expressed in muscle in either background. We visualized unc-52 splice junctions using the Sashimi plot function in the Integrated Genome Viewer (see Materials and Methods) (Figure 6B). In wild-type animals, 3’ exons were well-represented (Figure 6B). In etr-1(lq61), 3’ exons were significantly underrepresented and 5’ exons significantly overrepresented. These data suggest that ETR-1 proteins that include exon 8 are required to produce the long isoforms of unc-52 containing the 3’ exons in muscles.
As described above, unc-52(RNAi) resulted in AQR and PQR directional migration defects (Table 1). For RNAi, we utilized the ZC101.2 Source Bioscience II-9A20 clone, which is located located in the 3’ region of unc-52, the region that is underrepresented in etr-1(lq61) muscles (Figure 6B). unc-52(RNAi) animals displayed the paralyzed, arrested at two-fold stage (Pat) phenotype similar to strong loss-of-function alleles of unc-52 (Figure 7A and B) [59]. Notably, etr-1(lq61) animals, which do not express the 3’ exons from this unc-52 region in muscle, do not show the Pat phenotype. The viability of etr-1(lq61) mutants could be
due to a muscle-specific effect of ETR-1 on unc-52, with unc-52 long isoforms with 3’ exons expressed and functional in other tissues (e.g. hypodermis) in etr-1(lq61) mutants. Alternatively, it is possible that the phenomenon of RNAi amplification [60] results in all unc-52 isoforms being affected by RNAi, in which case conclusions about the roles of the 3’ exons cannot be drawn.
Despite embryonic lethality and the Pat phenotype, AQR and PQR were visible in unc-52(RNAi) arrested Pat animals. AQR and PQR displayed defects in the ability of AQR and PQR to migrate (Figure 7A and B) as well as defects in direction of migration (Table 3). In sum, ETR-1 proteins that contain exon 8 are required for the accumulation of long isoform transcripts of unc-52 containing 3’ exons in muscles. Targeting unc-52 with RNAi resulted in AQR and PQR defects, suggesting a role in AQR and PQR migration.
RNAi knockdown of lev-11 results in AQR and PQR migration defects
lev-11 was also identified by both exon representation and transcript accumulation in etr-1(lq61) (Table 2). lev-11 encodes a tropomyosin [61], which is a known target of vertebrate CELF1 [16]. lev-11 encodes multiple isoforms regulated in a tissue-specific manner [62, 63].
lev-11 exons 8 and 15 were significantly overrepresented in etr-1(lq61) muscle compared to wild-type muscle (Figure 8A). IGV-Sashimi splice junction analysis revealed that lev-11 exons 8 and 15 were included in etr-1(lq61) muscle and largely excluded in wild-type muscle (Figure 8B). These results suggest that ETR-1 proteins with exon 8 controlled alternative splicing of lev-11 in muscle cells, most notably removing exons 8 and 15 from lev-11 transcripts.
lev-11 RNAi resulted in AQR not migrating the full distance anteriorly, a phenotype also observed in lev-11 mutants (Table 3 and Figure 7C and D). The region at which AQR stopped migrating in both lev-11 RNAi and lev-11 mutants was highly stereotyped, in position 4 just posterior to the normal position. These results suggest that lev-11 might affect a specific developmental sign post rather than the general ability of AQR to migrate. In sum, our functional analysis indicates that genes with transcripts regulated by ETR-1 in muscles include those that control AQR and PQR neuron migration.