m6A modification of U6 snRNA modulates usage of two major classes of pre-mRNA 5’ splice site
Abstract
Alternative splicing of messenger RNAs is associated with the evolution of developmentally complex eukaryotes. Splicing is mediated by the spliceosome, and docking of the pre-mRNA 5’ splice site into the spliceosome active site depends upon pairing with the conserved ACAGA sequence of U6 snRNA. In some species, including humans, the central adenosine of the ACAGA box is modified by N6 methylation, but the role of this m6A modification is poorly understood. Here, we show that m6A modified U6 snRNA determines the accuracy and efficiency of splicing. We reveal that the conserved methyltransferase, FIONA1, is required for Arabidopsis U6 snRNA m6A modification. Arabidopsis fio1 mutants show disrupted patterns of splicing that can be explained by the sequence composition of 5’ splice sites and cooperative roles for U5 and U6 snRNA in splice site selection. U6 snRNA m6A influences 3’ splice site usage. We generalise these findings to reveal two major classes of 5’ splice site in diverse eukaryotes, which display anti-correlated interaction potential with U5 snRNA loop 1 and the U6 snRNA ACAGA box. We conclude that U6 snRNA m6A modification contributes to the selection of degenerate 5’ splice sites crucial to alternative splicing.
Editor's evaluation
This is an important paper reporting that an adenosine methyltransferase in the model plant Arabidopsis functions to target a key RNA component of the spliceosome, as in fission yeast, and thereby contributes to intron recognition as well as flowering timing. By contrast, the authors report no major role for the methyltransferase in targeting mRNAs, as reported in several previous studies in Arabidopsis. Overall, the approaches are convincing, although a correlative analysis identifying an intronic sequence feature characteristic of methyltransferase-sensitive introns is not followed up with tests to establish causality of this feature.
https://doi.org/10.7554/eLife.78808.sa0eLife digest
All the information necessary to build the proteins that perform the biological processes required for life is encoded in the DNA of an organism. Making these proteins requires the DNA sequence of a gene to be transcribed into a ‘messenger RNA’ (mRNA), which is then processed into a final, mature form. This blueprint is then translated to assemble the corresponding protein.
When an mRNA is processed, segments of the sequence that do not code for protein are removed and the remaining coding sequences are joined together in the right order. An intricate molecular machine known as the spliceosome controls this mechanism by recognising the ‘splice sites’ where coding and non-coding sequences meet. Depending on external conditions, the spliceosome can ‘pick-and-mix’ the coding sequences to create different processed mRNAs (and therefore proteins) from a single gene. This alternative splicing mechanism is often used to regulate when certain biological processes take place based on environmental cues; for example, the splicing of genes which control the timing of plant flowering is sensitive to ambient temperatures.
To investigate this mechanism, Parker et al. focused on Arabidopsis thaliana, a plant that blooms later when temperatures are low. This precise timing partly relies on a gene whose mRNA is efficiently spliced in the cold, resulting in an active form of its protein that blocks blooming. Parker et al. grew and screened many A. thaliana plants to find individuals that could flower early in the cold, in which splicing of this gene was disrupted.
A mutant fitting these criteria was identified and subjected to further investigation, which revealed that it could not produce FIONA1. In non-mutant plants, this enzyme chemically modifies one of the components of the spliceosome, a small nuclear RNA known as U6. Parker et al found that there are two types of splice site – one more likely to interact with U6 and another that preferentially interacts with another small nuclear RNA, U5. When FIONA1 is inactive (such as in the mutant identified by Parker et al.), splice sites that tend to strongly interact with U5 are selected. However, when the enzyme is active, splice sites that tend to bind with the chemically modified U6 are used instead.
Further work by Parker et al. showed that these two types of splice sites (‘preferring’ either U5 or U6) are found in equal proportions in the genomes of many species, including humans. This suggests that Parker et al. have uncovered an essential feature of how genomes are organised and splicing is controlled.
Introduction
Split genes are a defining characteristic of eukaryotic genomes (Plaschka et al., 2019). During pre-mRNA transcription, intervening sequences (introns) are excised, and the flanking sequences (exons) are spliced together. In developmentally complex eukaryotes, such as humans or Arabidopsis, most protein-coding genes have introns, and the cis-elements controlling splicing are degenerate. Such sequence variation facilitates alternative splice site selection that, in turn, permits the regulation of mRNA expression and the production of functionally different protein isoforms (Lee and Rio, 2015; Nilsen and Graveley, 2010). Consistent with this, alternative splicing is the foremost genomic predictor of developmental complexity (Chen et al., 2014).
Pre-mRNA splicing is carried out by the spliceosome (Plaschka et al., 2019; Wilkinson et al., 2020). This large, dynamic molecular machine comprises more than 100 proteins and 5 uridylate-rich snRNAs (UsnRNAs). Splicing requires two sequential transesterification reactions. The first reaction, called branching, occurs when the 2’ hydroxyl of the conserved intron branchpoint adenosine performs a nucleophilic attack on the 5’ splice site (5’SS), resulting in a cleaved 5’ exon and a branched lariat-intron intermediate. The second reaction, exon ligation, occurs via nucleophilic attack of the 3’ hydroxyl of the 5’ exon on the 3’ splice site (3’SS). Typically, the 5’SS is first recognised by U1 snRNP, whilst U2 snRNP recognises the branchpoint sequence of pre-mRNA. The preassembled U4/U6.U5 tri-snRNP subsequently joins the spliceosome. The conserved U6 snRNA ACAGA sequence replaces U1 snRNP at the 5’SS and loop 1 of U5 snRNA binds 5’ exon sequences adjacent to the 5’SS. In humans, selection of the 5’SS occurs during this hand-off to U6 and U5 snRNAs, and this step is decoupled from formation of the active site, thus potentially allowing for plasticity of 5’SS selection (Charenton et al., 2019; Fica, 2020). Subsequently, U6 pairs with U2 snRNA and intramolecular U6 snRNA interactions facilitate the positioning of the catalytic metal ions to effect branching and exon ligation (Wilkinson et al., 2020).
U6 snRNA is the most highly conserved UsnRNA, reflecting its crucial roles at the active site of the spliceosome (Wilkinson et al., 2020). Different base modifications are found in each of the UsnRNAs, including U6 snRNA, but the role of these modifications is poorly understood (Morais et al., 2021). The essential ACAGA sequence of Saccharomyces cerevisiae U6 snRNA is unmodified. However, in other species including Schizosaccharomyces pombe (Gu et al., 1996), the plant Vicia faba (Kiss et al., 1987), and human (Shimba et al., 1995), the central adenosine in the corresponding sequence is modified by methylation at the N6 position: ACm6AGA. In S. cerevisiae, U6 snRNA ACAGA recognises the stringently conserved GUAUGU sequence at the 5’ end of introns, with the central adenosine making a Watson-Crick base pair with the almost invariant U at the +4 position of the 5’SS (5’SS U+4) (Neuvéglise et al., 2011; Wan et al., 2019). Conversely, in species with m6A-modified U6 snRNA, the 5’SS can be degenerate and the identity of the base at the +4 position varies but is usually enriched for A. Therefore, the pairing between the 5’SS and the U6 snRNA ACAGA box is unlikely to be driven primarily by canonical Watson-Crick base pairs in these species. The precise function of the U6 snRNA m6A modification is unknown.
In humans, the U6 snRNA ACAGA sequence is methylated by the conserved methyltransferase METTL16 (Aoyama et al., 2020; Pendleton et al., 2017; Warda et al., 2017). This modification depends upon a specific sequence and distinct structural features of U6 snRNA that are recognised by METTL16. Hairpin sequences that mimic these features of U6 snRNA are found within the 3’UTR of S-adenosylmethionine (SAM) synthetase MAT2A mRNA and are also methylated by METTL16 (Pendleton et al., 2017; Shima et al., 2017). SAM synthetase is the enzyme responsible for production of the methyl donor SAM, which is required for methylation reactions in the cell. Binding of METTL16 to MAT2A mRNA influences splicing and/or stability, regulating MAT2A expression and SAM levels (Pendleton et al., 2017; Shima et al., 2017). In Caenorhabditis elegans, the METTL16 orthologue, METT-10, methylates U6 snRNA and also influences SAM levels by targeting SAM synthetase (sams) genes (Mendel et al., 2021). METT-10 methylates the 3’SS of sams-3/4 intron 2, although in a different sequence and structural context to U6 snRNA (Mendel et al., 2021). In S. pombe, the METTL16 orthologue, MTL16, appears to target only U6 snRNA (Ishigami et al., 2021). S. pombe mutants defective in MTL16 function show less efficient splicing of some introns, resulting in increased levels of intron retention. The S. pombe introns sensitive to loss of MTL16 function are distinguished by having adenosine at the +4 position of the intron (5’SS A+4). However, the effect of 5’SS A+4 is partially suppressed in introns that have stronger predicted base pairing between U5 snRNA loop 1 and the 5’ exon (Ishigami et al., 2021).
Alternative splicing plays crucial roles in gene regulation, affecting diverse biological processes including the control of flowering in response to ambient temperature (Airoldi et al., 2015; Capovilla et al., 2017). Slight changes in ambient temperature can profoundly influence flowering time (Andrés and Coupland, 2012). Indeed, documented changes in the flowering times of many plant species have provided some of the best biological evidence of recent climate change (Fitter and Fitter, 2002). At cooler ambient temperatures, pre-mRNAs encoding repressors of Arabidopsis flowering, FLM (FLOWERING LOCUS M) and MAF2 (MADS AFFECTING FLOWERING 2), are efficiently spliced (Airoldi et al., 2015; Balasubramanian et al., 2006; Lee et al., 2013; Lutz et al., 2015; Posé et al., 2013; Scortecci et al., 2001). As a result, FLM and MAF2 proteins are produced at cooler temperatures and form higher-order protein complexes with other MADS box factors, such as MAF3 (MADS AFFECTING FLOWERING 3), SVP (SHORT VEGETATIVE PHASE), and FLC (FLOWERING LOCUS C), to directly repress the expression of the genes FT (FLOWERING LOCUS T) and SOC1 (SUPPRESSOR OF OVEREXPRESSION OF CO 1), which promote flowering (Gu et al., 2013). At elevated ambient temperatures, the splicing of MAF2 and FLM introns is adjusted, in different ways, effecting increased levels of non-productive transcripts (Airoldi et al., 2015; Balasubramanian et al., 2006; Lee et al., 2013; Lutz et al., 2015; Posé et al., 2013; Scortecci et al., 2001). Consequently, the abundance of active MAF2 and FLM proteins is compromised at elevated temperatures; repression of genes promoting flowering is relieved, and flowering is enabled (Airoldi et al., 2015; Capovilla et al., 2017; Posé et al., 2013; Rosloski et al., 2013). In the case of MAF2, this temperature-responsive splicing is characterised by a progressive increase in retention of intron 3 as ambient temperature increases (Airoldi et al., 2015). The thermosensory mechanism that underpins the temperature-responsive alternative splicing of MAF2 and FLM is unknown.
We conducted a two-step mutant screen for factors required for the efficient splicing of MAF2 pre-mRNA intron 3 at cooler temperatures. This screen identified an early flowering mutant disrupting the gene encoding the Arabidopsis METTL16 orthologue, FIONA1 (FIO1; Kim et al., 2008). We show that FIO1 is required for m6A modification of U6 snRNA. We detect widespread changes in pre-mRNA splicing in fio1 mutants that can be explained by the identity of the base at the 5’SS +4 position of affected introns and by cooperative roles for U5 and U6 snRNA in 5’SS selection. We also reveal that U6 snRNA interactions at 5’SSs affect 3’SS choice. Analysis of other annotated genomes reveal the existence of two major classes of 5’SS. Our findings suggest cooperative and compensatory roles of U5 and U6 snRNA may contribute to 5’SS selection in a range of eukaryotes.
Results
Identification of FIO1 in a two-step screen for early flowering mutants with increased retention of MAF2 intron 3
To identify factors responsible for promoting efficient splicing of MAF2 intron 3 at low ambient temperature, we performed a two-step genetic screen using ethyl methanesulfonate (EMS) as a mutagen (Figure 1A). In the first step, we screened 15,000 EMS lines for early flowering Arabidopsis mutants at 16°C. In the second step, the earliest flowering 100 individuals were re-screened for enhanced MAF2 intron 3 retention at 16°C, using RT-PCR. Two independent mutants (EMS-129 and EMS-213) showed early flowering and increased levels of MAF2 intron 3 retention at 16°C. We focus here on EMS-129 (Figure 1—figure supplement 1A). The causal mutation for the early flowering phenotype of EMS-129 was mapped by back-crossing to the M0 parental line to generate a segregating population (Figure 1—figure supplement 1B). Early flowering plants and aphenotypic sisters were grouped into separate pools for genomic DNA sequencing. EMS-induced G to A transitions were identified, and allele fractions were compared between the two pools (Javorka et al., 2019). A total of 19 homozygous G to A transitions in a 2.6 Mb region on chromosome 2 were associated with the early flowering phenotype (Figure 1B). Of these, only the SNP at position 9,041,454 of chromosome 2, which disrupts the 5’SS of AT2G21070 intron 2, was predicted to have a major impact on gene expression. An EMS mutation in AT2G21070 has previously been described and named fio1-1 (Kim et al., 2008). Crosses between EMS-129 and fio1-1 confirmed allelism because 100% of F1 progeny from EMS-129 lines crossed to fio1-1 flowered early (35 plants) compared to 0% of F1 progeny from EMS-129 lines crossed to Col-0 (15 plants). We isolated an independent transfer-DNA insertion in AT2G21070 (Alonso et al., 2003), which we refer to as fio1-3. Henceforth, we refer to EMS-129 as fio1-4. We summarise the fio1-1, fio1-3, and fio1-4 alleles in Figure 1C. fio1-1 and fio1-3 show increased levels of MAF2 intron 3 retention at low temperature (Figure 1D) and like maf2 mutants, flower early (Figure 1E–F). fio1 alleles share additional visible phenotypes, including reduced apical dominance and stature (Figure 1F). FIO1 encodes the Arabidopsis orthologue of the human methyltransferase, METTL16 (Kim et al., 2008; Pendleton et al., 2017).
-
Figure 1—source data 1
- https://cdn.elifesciences.org/articles/78808/elife-78808-fig1-data1-v2.xlsx
-
Figure 1—source data 2
- https://cdn.elifesciences.org/articles/78808/elife-78808-fig1-data2-v2.txt
-
Figure 1—source data 3
- https://cdn.elifesciences.org/articles/78808/elife-78808-fig1-data3-v2.txt
-
Figure 1—source data 4
- https://cdn.elifesciences.org/articles/78808/elife-78808-fig1-data4-v2.xlsx
-
Figure 1—source data 5
- https://cdn.elifesciences.org/articles/78808/elife-78808-fig1-data5-v2.xlsx
-
Figure 1—source data 6
- https://cdn.elifesciences.org/articles/78808/elife-78808-fig1-data6-v2.zip
FIO1-dependent methylation of poly(A)+ RNA is rare or absent
To reveal which RNAs were methylated by FIO1, we first used two orthogonal approaches to determine the impact of FIO1 on poly(A)+mRNA m6A modification: liquid chromatography tandem mass spectrometry (LC-MS/MS) and nanopore direct RNA sequencing (DRS). As a positive control for these experiments, we included a hypomorphic fip37-4 allele defective in the Arabidopsis orthologue of the METTL3 m6A writer complex component WTAP (Zhong et al., 2008).
LC-MS/MS analysis (Figure 2A) revealed an 81.8% decrease in m6A levels in fip37-4 compared to Col-0 (t test p=2.2 × 10−6). In contrast, a 6.7% decrease in m6A levels was observed in poly(A)+RNA purified from fio1-1 (t test p=0.40). A similarly modest decrease in mRNA m6A levels has recently been reported for other fio1 alleles (Wang et al., 2022).
-
Figure 2—source data 1
- https://cdn.elifesciences.org/articles/78808/elife-78808-fig2-data1-v2.xlsx
-
Figure 2—source data 2
- https://cdn.elifesciences.org/articles/78808/elife-78808-fig2-data2-v2.xlsx
-
Figure 2—source data 3
- https://cdn.elifesciences.org/articles/78808/elife-78808-fig2-data3-v2.xlsx
-
Figure 2—source data 4
- https://cdn.elifesciences.org/articles/78808/elife-78808-fig2-data4-v2.xlsx
We previously used nanopore DRS to map m6A sites dependent on the Arabidopsis METTL3-like writer complex component VIRILIZER (VIR), revealing m6A enriched in a DRm6ACH consensus in the 3’ terminal exon of protein coding mRNAs (Parker et al., 2020). We used a similar nanopore DRS experiment to map FIO1-dependent m6A. We sequenced poly(A)+RNA purified from four biological replicates each of wild-type (Col-0), fio1-1 and fip37-4, resulting in a total of 22.7 million mapped reads. The corresponding sequencing statistics are detailed in Supplementary file 1. We used the software tool Yanocomp to map differences in mRNA modifications (Parker et al., 2021a). Applying this approach, we identified 37,861 positions that had significantly different modification rates in a three-way comparison between Col-0, fip37-4, and fio1-1 (Figure 2B). Of these, 97.9% had significant changes in modification rate in a pairwise comparison between fip37-4 and Col-0. In contrast, only 7.6% had altered modification rates in fio1-1. Of these, 85.8% also had altered modification rates in fip37-4 (Figure 2B), with larger effect sizes (Figure 2—figure supplement 1A), indicating that they are METTL3-like writer complex-dependent m6A sites whose modification rate is indirectly affected by loss of FIO1. The FIP37-dependent modification sites, including those with small effect-size changes in fio1-1, were found in a DRACH consensus in 3’UTR regions (Figure 2—figure supplement 1B–C), consistent with established features of m6A modifications deposited by the Arabidopsis METTL3-like writer complex (Parker et al., 2020).
We used the antibody-based technique miCLIP to perform orthogonal validation of predicted modification sites (Parker et al., 2020): FIP37-dependent m6A sites were well supported, with 52.8% less than 5 nt from an miCLIP peak (Figure 2B). In contrast, of the 408 sites discovered only in fio1-1, just 21.3% (87) were less than 5 nt from an miCLIP peak, indicating that the majority are likely to be false positives. The 87 positions with miCLIP support were found in a DRACH consensus in 3’UTRs, suggesting that they are false negative FIP37-dependent m6A sites rather than genuine FIO1-dependent sites (Figure 2—figure supplement 1D–E).
We did not find FIO1-dependent methylation changes in transcripts of the four Arabidopsis homologues of SAM-synthetase (Figure 2—figure supplement 2A–D). Consequently, the METTL16-dependent regulation of SAM homeostasis identified in metazoans (Mendel et al., 2021; Pendleton et al., 2017; Shima et al., 2017) appears not to be conserved in Arabidopsis. Consistent with previous evidence that m6A in 3’UTRs affects cleavage and polyadenylation (Parker et al., 2020), we identified 1104 genes with altered poly(A) site choice in fip37-4, compared to only 49 in fio1-1 (Figure 2—figure supplement 1F). Overall, these findings indicate that FIO1-dependent m6A sites in poly(A)+mRNA are rare or absent.
m6A modification of U6 snRNA depends upon FIO1
Since FIO1 orthologues methylate U6 snRNA, we asked if FIO1 was required to modify Arabidopsis U6 snRNA. We performed two experiments using two independent anti-m6A antibodies to immunopurify methylated RNAs from Col-0 and fio1 alleles. For control, we used U2 snRNA, which is m6A methylated in humans, but not by METTL16 (Chen et al., 2020). Using RT-qPCR, we could detect equivalent levels of U6 snRNA in the input RNA purified from different genotypes, suggesting that the abundance of U6 snRNA is unaffected by loss of FIO1 function (Figure 2—figure supplement 3A–B). We detected enrichment of U6 and U2 snRNAs in RNA immunopurified with anti-m6A antibodies from Col-0 (Figure 2C–D). The enrichment of U2 snRNA in these experiments was unaffected by loss of FIO1 function. In contrast, we identified significant depletion of U6 snRNA in the anti-m6A immunopurified RNA from fio1-1, fio1-3, and fio1-4 alleles (Figure 2C–D). These data are consistent with recent reports (Wang et al., 2022; Xu et al., 2022) indicating that Arabidopsis U6 snRNA is m6A-modified at the conserved position of the ACm6AGA box and that this modification requires active FIO1.
Widespread disruption of pre-mRNA splicing in fio1 mutants
Having found that FIO1 is required for U6 snRNA m6A modification, we examined changes in gene expression and pre-mRNA splicing in fio1 mutants. Nanopore DRS is insightful for mapping the complexity of RNA processing and modification, but throughput currently limits the statistical power to detect splicing changes, and basecalling errors complicate the alignment of exon/intron boundaries (Parker et al., 2021b). Consequently, we used Illumina RNA-seq of poly(A)+RNA to analyse gene expression and splicing. Since FIO1 affects the splicing of MAF2, which is sensitive to changes in temperature, we included plants grown at 20°C and shifted to either 4, 12, or 28°C for 4 hr prior to harvesting and poly(A)+RNA purification. We sequenced six biological replicates each of wild-type Col-0 and fio1-3, generating a minimum of 46 million paired end reads, 150 bp in length, per replicate. On average, 97.9% of read pairs were mappable per replicate, resulting in a total of 3.6 billion mapped read pairs. A summary of the sequencing read statistics is given in Supplementary file 1. The RNA-seq data clearly reveal disruption of FIO1 gene expression in fio1-3 (Figure 3—figure supplement 1A).
To detect cryptic splice sites that might be activated in fio1 mutants, we built a condition-specific reference transcriptome derived from the Col-0/fio1-3 Illumina RNA-seq and Col-0/fio1-1 nanopore DRS data using the software tool Stringtie2 (Kovaka et al., 2019). We quantified the expression of these transcripts using Salmon (Patro et al., 2017), then used SUPPA2 to calculate event-level percentage changes in splicing, which also classifies different types of splicing event (Trincado et al., 2018). These percent spliced indices (PSIs) were used to fit linear models for each splicing event to identify changes in splicing dependent on temperature, genotype, and temperature×genotype interactions. Analysis of the Col-0 temperature-shift data served as a control. Consistent with previous studies (Calixto et al., 2018), we found that intron retention is the predominant class of alternative splicing event detected when Arabidopsis is subjected to different temperatures (Figure 3A). In contrast, when we repeated this analysis to classify alternative splicing events that differ between fio1-3 and Col-0, we found that a larger proportion of alternative splicing events were classified as alternative 5’SS usage – 34.4%, compared to only 18.6% of temperature-dependent events (Figure 3B, Figure 3—figure supplement 1B). In addition, we detected changes in the PSI of retained introns, exon skipping, and alternative 3’SS selection (Figure 3B, Figure 3—figure supplement 1B). There was a significant overlap between the alternative splicing events that were sensitive to loss of FIO1, and those that were sensitive to temperature (hypergeometric-test p<1 × 10−16). However, in 64.2% of fio1-sensitive events, loss of FIO1 did not alter splicing responses to temperature (Figure 3B). Of the remaining 2505 splicing events that did have altered temperature sensitivity in the absence of FIO1, 38.4% were alternative 5’SSs, and of these, 69.9% had greater sensitivity to loss of FIO1 at 28°C than at 4°C. This suggests that canonical 5’SS selection at elevated temperature requires FIO1 (Figure 3—figure supplement 1C).
Analysis of the RNA-seq data confirmed the increased retention of MAF2 intron 3 in fio1-3, consistent with the results of our two-step mutant screen (Figure 3C). Loss of FIO1 decreased the splicing efficiency of MAF2 intron 3 by approximately equivalent amounts at all temperatures (Figure 3D), implying that FIO1 is not required to generate MAF2 intron 3 temperature sensitivity.
We detected changes in gene expression consistent with the early flowering phenotype of fio1 mutant alleles (Figure 3—figure supplement 2A–G). For example, FT and SOC1 mRNA levels were increased in fio1-3 (Figure 3—figure supplement 2A–B). In contrast, the mRNA levels of the floral repressors FLM, MAF2, MAF3, MAF4, MAF5, and FLC were reduced (Figure 3—figure supplement 2C–H). In the case of FLC, both sense and antisense transcript levels were reduced in fio1-3 (Figure 3—figure supplement 2H–I), but there were no detectable changes in splicing patterns or in the alternative polyadenylation of the antisense transcripts (Figure 3—figure supplement 2J). In these experimental conditions, loss of FIO1 did not affect RNA-level expression or splicing of either CO or SVP (Figure 3—figure supplement 2K–L). In contrast, fio1-3 mutants exhibited defects in the splicing of pre-mRNA encoding MAF2, MAF3, FLM, and other genes that influence flowering time, as detailed further below.
Overall, we conclude that the predominant molecular phenotype resulting from loss of FIO1 function is a major disruption to patterns of splicing in a manner that is mostly independent of temperature. To understand the basis of these splicing differences, we analysed each of these splicing classes separately.
Alternative selection of 5’SSs in fio1-3 can be explained by U6 and U5 snRNA target sequences
We identified 2369 changes in 5’SS choice caused by loss of FIO1 function. Alternative 5’SSs were almost equally likely to be selected upstream or downstream, with many examples of 5’SS shifts of exactly −4 nt,+4 nt, and +5 nt (Figure 4—figure supplement 1A). We asked if differences in cis-element features could account for FIO1-dependent 5’SS selection. Comparing 5’SSs that exhibited reduced selection in fio1-3 with their corresponding alternative 5’SSs, which were more frequently selected in fio1-3, we found a significant difference in base composition at the −2 to +5 positions (G-test p<1 × 10−16). Specifically, we identified a //GURAG motif (R=A or G) at 5’SSs sensitive to the loss of FIO1 function (Figure 4A). This motif appeared more commonly in fio1-sensitive 5’SSs than in the full set of 5’SSs observed in the RNA-seq dataset (Figure 4—figure supplement 1B) and is consistent with recognition by the U6 snRNA ACAGA box (Kandels-Lewis and Séraphin, 1993; Lesser and Guthrie, 1993; Sawa and Abelson, 1992). In contrast, the alternative 5’SSs selected more frequently in fio1-3 were characterised by an AG//GU motif (Figure 4B). This motif is consistent with U5 snRNA loop 1 recognition of the upstream exon sequences (Galej et al., 2016; Newman and Norman, 1992). There was no significant difference between the positional base frequencies for 5’SSs that were preferred at lower and elevated temperatures in Col-0 (G-test p=0.28; Figure 4—figure supplement 1C–D), demonstrating that altered 5’SS sequence preference is specific to the fio1-3 mutant and not a general feature of alternative splicing.
The widespread changes in 5’SS selection are evident at individual loci. For example, at AtSAR1 (AT1G33410), loss of FIO1 function, at 20°C, results in a 52.2% reduction in the use of the normally near constitutive 5’SS in intron 21, UC//GUGAG, with a reciprocal increase in the use of a cryptic 5’SS UG//GUAUU 26 nt downstream (Figure 4C and D). At elevated temperatures, the selection of this cryptic 5’SS becomes increasingly preferred (Figure 4D). The change in usage of this 5’SS is confirmed in the fio1-1 allele using orthogonal nanopore DRS (Figure 4E). We could also identify cryptic alternative 5’SSs in the second intron of MAF2 (Figure 4—figure supplement 2A–B), and the third exon of another floral repressor, FLM, used in fio1-3 (Figure 4—figure supplement 2C–D). The overall effect of these (and other) splicing changes is to reduce the fraction of productive AtSAR1, MAF2, and FLM transcripts in fio1-3 (Figure 4—figure supplement 2E–G). Finally, at the gene encoding the METTL3-like writer complex component MTB (the orthologue of human METTL14), loss of FIO1 function resulted in a 40.4% increase in the use of a cryptic 5’SS in intron 4 that introduces a premature termination codon into the MTB open reading frame (Figure 4—figure supplement 2H–I). Consequently, loss of FIO1 may have indirect effects on m6A deposition by the METTL3-like writer complex.
Remarkably, the presence or absence of an A at the +4 position (A+4) correctly separates 78.6% of the 5’SSs exhibiting decreased and increased usage in fio1-3 (Figure 4F). In humans, where pairing of the 5’SS to U6 snRNA occurs before activation of the spliceosome, the 5’SS A+4 faces the m6A of the U6 snRNA ACm6AGA box in the B complex before docking of the 5’SS in the active site (Bertram et al., 2017). Of the 5’SSs with increased usage in fio1-3, 48.2% have U+4, which could make a Watson–Crick base pair with the corresponding unmethylated residue of U6 snRNA. In total, 61.8% of 5’SS changes in fio1-3 is associated with a switch from A+4 to B+4 (B=C, G or U), of which 58.8% are A→U+4 (Figure 4—figure supplement 3A). In comparison, only 3.8% of alternative 5’SS pairs are reciprocal B→A+4 switches, indicating that this shift is strongly directional. A further 8.0% of alternative 5’SS pairs are S→U+4 (S=C or G), suggesting that a Watson–Crick A–U base-pair is favoured when U6 snRNA is not m6-modified. Surprisingly, 22.3% of alternative 5’SS pairs have the same base at the +4 position (Figure 4—figure supplement 3A). However, of these 5’SS pairs, the majority are associated with a G→H+5 (H=A, C or U) and/or H→G–1 switch that weakens interactions with the U6 snRNA ACAGA box and strengthens U5 snRNA loop 1 interactions (Figure 4—figure supplement 3B).
We next examined how the effect size (absolute ΔPSI) of splicing changes correlated with the base at the +4 position of the 5’SS. The largest effect sizes were associated with A→U+4 (Figure 4G). In contrast, alternative 5’SS pairs with an A→A+4 shift had smaller effect sizes. We found that fio1-3-sensitive 5’SSs with 5’SS G+5 had larger effect sizes but that G+5 at the alternative 5’SS had less effect, suggesting that G+5 is only deleterious in fio1-3 when in combination with A+4 (Figure 4—figure supplement 3C). Finally, we found that fio1-3-sensitive 5’SSs with AG//GU motifs had smaller effect sizes, indicating that favourable interactions with U5 snRNA loop 1 can suppress the effect of unfavourable U6 snRNA interactions in fio1-3 (Figure 4H).
In summary, these analyses demonstrate widespread change in 5’SS selection in fio1 mutants. The 5’SSs with a strong match to the U6 snRNA ACAGA box and/or an A+4 are most sensitive to loss of FIO1. This sensitivity can be suppressed by a strong match in the 5’ exon to U5 snRNA loop 1. Alternative 5’SSs that are used more in fio1 have U+4, as well as 5’ exon sequence features that favour recognition by U5 snRNA loop 1.
FIO1-dependent changes in intron retention and exon skipping can be explained by U6 and U5 snRNA target sequences
In addition to altered 5’SS selection, we detected many instances of intron retention and exon skipping in fio1-3. We therefore asked if these FIO1-sensitive splicing phenotypes were associated with specific sequence motifs.
We identified 2576 introns with altered levels of retention in fio1-3, of which 55.8% had increased retention, and 44.2% had decreased retention (i.e. more splicing). Analysis of 5’SS sequences at introns with increased retention indicates that 76.9% had A+4, as well as a weaker interaction potential with U5 snRNA loop 1 (Figure 5A). For example, at WNK1 (AT3G04910), we detected increased retention of intron 5 (CG//GUGAG) in fio1-3 (Figure 5B–C) and fio1-1 (Figure 5D) using Illumina RNA-seq and nanopore DRS analysis, respectively. This indicates that these introns are normally recognised by strong U6 snRNA interactions and are less efficiently spliced in the absence of FIO1-dependent U6 snRNA m6A modification. The remaining introns with significantly altered retention are more efficiently spliced in fio1-3. These introns tend to have stronger matches to U5 snRNA loop 1 in the 5’ exon, and only 34.7% have A+4, whilst 44.3% have U+4 (Figure 5A). This suggests that loss of m6A from U6 snRNA actually increases the splicing efficiency of introns with U at the +4 position of the 5’SS.
A common form of alternative splicing involves the excision of an exon and flanking introns in a splicing event called exon skipping. These events involve two 5’SSs: those at the upstream intron and those at the downstream intron (i.e. at the 3’ end of the skipped exon; Figure 5E–F). We identified 442 exons with significantly different levels of skipping in fio1-3. Of these, 53.6% have increased levels of inclusion and 46.4% have increased levels of skipping.
In cases where exon skipping increased in fio1-3, more of the 5’SSs at downstream introns had A+4 compared to the 5’SSs at upstream introns. Conversely, more of the 5’SSs at upstream introns had strong U5 snRNA loop 1 interacting sequences (Figure 5E, Figure 5—figure supplement 1A). This is consistent with a relative weakening of the recognition of the 5’SS at the downstream intron in fio1-3. For example, at the gene encoding the floral repressor MAF3 (AT5G65060), there is an increase in skipping of exon 2 in fio1-3 (Figure 5—figure supplement 1B–C). The 5’SS at the 3’ end of this exon is UA//GUAAG, whereas the upstream 5’SS (AA//GUAAG) is a stronger match to U5 snRNA loop 1.
In cases where exon inclusion increased in fio1-3, we found that the majority of 5’SSs at upstream introns had A+4, whereas a plurality of 5’SSs at downstream introns had U+4 (Figure 5F, Figure 5—figure supplement 1A). For example, at PTB1 (AT3G01150), there is an increase in inclusion of exon 3 in fio1-3 (Figure 5—figure supplement 1D–E). The 5’SS at the 3’ end of this exon is AG//GUGUC, whereas the upstream 5’SS is UG//GUGAG. Although both the upstream and downstream 5’SSs are used when the exon is retained, it is possible that the relative strengthening of recognition of the downstream 5’SS in fio1-3 improves exon definition.
In summary, we can identify and characterise changes in intron retention and exon skipping sensitive to loss of FIO1 function. Although the outcome of these splicing events is different from alternative 5’SS selection, they are associated with the same changes in 5’SS sequence as those associated with alternative 5’SS selection.
Alternative 3’SS usage in fio1-3 can be explained by U6 snRNA-dependent interactions between 5’ and 3’SSs
We observed 1484 instances where 3’SS selection was altered in fio1-3. In some cases, these alternative 3’SSs were linked with altered 5’SS choice or intron retention, but this accounted for a relatively small proportion of events (Figure 6—figure supplement 1A). Alternative 3’SS selection was equally likely to switch in an upstream or downstream direction in fio1-3 (Figure 6A). There was a strong enrichment for very local switching of 3’SSs, with 37.9% of alternative 3’SSs occurring within 6 nt of the fio1-3-sensitive 3’SS, and 18.9% of examples occurring at exactly 3 nt upstream or downstream. These examples correspond to NAGNAG-like acceptors, which have previously been characterised in multiple species including human and Arabidopsis (Bradley et al., 2012; Hiller et al., 2004; Schindler et al., 2008).
To examine the 5’SSs associated with alternative 3’SSs used in fio1-3, we separated examples with increased upstream and downstream 3’SS usage and performed motif analysis. At introns with a relative increase in upstream 3’SS usage in fio1-3, we found that 80.5% of the corresponding 5’SSs were characterised by A+4 (Figure 6B). For example, at LHY (AT1G01060), we identified an AA//GUAAG 5’SS and a UAG\\CAG\\ alternative 3’SS pair in intron 7 (Figure 6D). In fio1-3, at 20°C, there is a shift in favour of the upstream UAG\\ 3’SS (Figure 6E). This 3’SS switch is supported by orthogonal nanopore DRS analysis of the fio1-1 allele (Figure 6F). Conversely, when we analysed the features of introns where a relative increase in downstream 3’SS usage was detected in fio1-3, we found that only 37.0% of these 5’SSs had A+4, whereas 48.8% had U+4 (Figure 6C). For example, at MAF2 intron 4, we identified an AG//GUAUU 5’SS and a UAG\\ACAG\\ alternative 3’SS pair (Figure 6—figure supplement 1B). In fio1-3, at 20°C, there is a shift in favour of the downstream CAG\\ 3’SS (Figure 6—figure supplement 1C).
3’SS choice involves scanning downstream from the branchpoint to the first available 3’SS motif (Smith et al., 1993). However, competition with downstream 3’SSs can occur within a short range, such as in NAGNAG acceptors. Downstream 3’SSs in fio1-sensitive alternative 3’SS pairs were more likely than upstream 3’SSs to have a cytosine as the −3 position (Figure 6B–C). However, this appears to reflect differences in the background rate of C−3 upstream and downstream (Figure 6—figure supplement 1D), rather than a change in 3’SSs motif competitiveness in fio1-3 (Bradley et al., 2012; Smith et al., 1993). To analyse whether the −3 position contributes to changes in alternative 3’SS usage in fio1, we separated fio1-sensitive 3’SS pairs where the base at the −3 position was the same at both the upstream and downstream 3’SS, from examples which had different bases at the −3 position. We found that among both sets of 3’SS pairs, 5’SS A+4 was still associated with increased upstream usage in fio1-3, whereas U+4 was associated with increased downstream usage (Figure 6G). This demonstrates that changes in 3’SS usage in fio1-3 are caused by a change in the competitiveness of distal 3’SSs, irrespective of 3’SS motif.
These findings link selection of the 3’SS and 5’SS. When U6 snRNA recognition of the 5’SS is favoured by either an m6A:A+4 interaction in WT Col-0 or A:U+4 interaction in fio1-3, this increases the usage of distal 3’SSs. When the U6 snRNA/5’SS interaction is less favoured, then upstream 3’SSs are more likely to be used. We conclude that the interactions of U6 snRNA with the 5’SS can influence usage of competing 3’SSs.
Reanalysis of RNA-seq data from independent fio1 alleles confirms disruption to splicing
In contrast to our findings, three recently published studies of fio1 mutants failed to identify an impact on splicing (Wang et al., 2022; Xu et al., 2022; Sun et al., 2022). We reanalysed these published data to understand the basis for this difference. Two studies performed Illumina RNA-seq analysis (Wang et al., 2022; Sun et al., 2022) but were not well powered with only two replicates per condition. The Illumina RNA-seq data from Wang et al. were of poor quality, with low mapping rates and one unusable replicate due to PCR-bias (Figure 6—figure supplement 2A). Nevertheless, we were able to detect 384 alternative 5’SS events by using our condition-specific transcriptome assembly and performing differential splicing analysis with SUPPA2. These events were associated with a shift from 5’SSs with A+4 to alternatives with U+4 (Figure 6—figure supplement 2B). We reanalysed the Illumina RNA-Seq data from Sun et al. (two replicates per condition and 20–26 million reads per replicate) to identify 923 alternative 5’SS events in fio1-1 and 947 alternative 5’SS events in a CRISPR-induced mutation (fio1-5; Figure 6—figure supplement 3A–B). There were 675 alternative 5’SS events common to both alleles, and 556 (82.4%) of these were identified in our own analysis of fio1-3 (Figure 6—figure supplement 3C). A third study analysed gene expression changes using nanopore DRS (Xu et al., 2022). Using these data, we could confirm splicing changes observed in our own data, for example at AtSAR1 intron 21, WNK1 intron 5, and MTB intron 4 (Figure 6—figure supplement 4A–C). These changes were also visible in nanopore DRS data generated by Sun et al. (Figure 6—figure supplement 4D–F). We conclude that a combination of under-powered datasets and methodological choices in data analysis explain previous failures to link the loss of U6 snRNA methylation in fio1 mutants to detectable changes in splicing.
Sequences targeted by U5 and U6 snRNAs are anticorrelated in the 5’SSs of developmentally complex eukaryotes
Our results indicate that a strong match to U5 snRNA loop 1 in the 3’ end of the upstream exon offsets the effects of 5’SS A+4 unfavourability in fio1-3. We, therefore, reasoned that if U5 and U6 cooperate in 5’SS selection, then strong U5 snRNA loop 1 interactions may globally compensate weaker U6 snRNA ACAGA box interactions and vice versa. To test this hypothesis, we used all annotated Arabidopsis 5’SS sequences to generate position-specific scoring matrices (PSSMs) for the −2 to −1 positions corresponding to the U5 snRNA loop 1 interacting region, and the +3 to +5 positions corresponding to the U6 snRNA ACAGA box interacting region. We then used these PSSMs to score the U5 and U6 snRNA interaction log-likelihood of each individual 5’SS. 5’SSs with AG//GU motifs have high-scoring U5 PSSM scores, and those with //GURAG motifs have high-scoring U6 PSSM scores. This analysis revealed that U5 and U6 snRNA PSSM scores are indeed negatively correlated in Arabidopsis (Spearman’s ρ = −0.36, p<1 × 10−16). We found that 5’SSs lacking a //GURAG motif had significantly higher U5 PSSM scores than //GURAG 5’SSs (Figure 7A and B), whereas BH//GU 5’SSs had significantly higher U6 PSSM scores than AG//GU 5’SSs (Figure 7C and D). We found similar anticorrelations in U5 and U6 PSSM scores in other species, including C. elegans (ρ = −0.38, p<1 × 10−16; Figure 7—figure supplement 1A), Drosophila melanogaster (ρ = −0.30, p<1 × 10−16; Figure 7—figure supplement 1B), Danio rerio (ρ = −0.36, p<1 × 10−16; Figure 7—figure supplement 1C), and Homo sapiens (ρ = −0.17, p<1 × 10−16; Figure 7—figure supplement 1D). These findings indicate that there are two major classes of 5’SS in distinct eukaryotes: //GURAG and AG//GU. These two classes occur with approximately equal frequency in the metazoan genomes that we analysed (Figure 7—figure supplement 1A–D).
Opposing U5 and U6 snRNA interaction potential is a feature of alternative 5’SS pairs in some eukaryote genomes
Given that many 5’SSs appear to have either strong U5 or U6 snRNA interacting sequences, we speculated that opposing U5 and U6 snRNA interaction strengths at pairs of alternative 5’SSs could contribute to alternative splicing. To test this hypothesis, we calculated the log odds ratio of U5 and U6 PSSM scores for all pairs of annotated alternative 5’SSs in Arabidopsis. We found that these ratios were negatively correlated (Spearman’s ρ = −0.25, p<1 × 10−16), indicating that alternative 5’SS pairs tend to have complementary strengths with respect to U5 and U6 snRNA recognition (Figure 7E). The strength of this complementary relationship varies in different organisms: in humans, we found a stronger negative correlation between relative U5 and U6 PSSM scores (ρ = −0.40, p<1 × 10−16; Figure 7F), whereas in C. elegans, the correlation was weaker (ρ = −0.15, p=8.61 × 10−7; Figure 7—figure supplement 2A), and in D. rerio, there was no correlation (ρ = −0.01, p=0.58; Figure 7—figure supplement 2B). We conclude that changing the relative favourability of 5’SSs with stronger U5 or U6 snRNA interacting sequences, as occurs in fio1-3, could be a mechanism contributing to alternative 5’SS choice.
Discussion
The splicing of mRNAs encoding regulators of flowering time is disrupted in fio1
We identified fio1 through a two-step mutant screen designed to reveal factors that control the splicing of MAF2 intron 3. Importantly, temperature-sensitive splicing of MAF2 intron 3 links ambient temperature to flowering time control. fio1 mutants have increased retention of MAF2 intron 3 compared to WT Col-0. However, splicing of MAF2 intron 3 remains responsive to temperature in fio1 mutants, suggesting that FIO1 is not the thermosensor in this process.
We found that disruption of FIO1 alters the splicing of mRNA encoding not only MAF2 but other floral repressors too, including FLM and MAF3. Detectable levels of sense and antisense RNAs at the locus encoding the floral repressor FLC were reduced in fio1, but no splicing changes were detected. FLM, MAF2, MAF3, and FLC function together with SVP in higher-order protein complexes to repress the expression of FT and SOC1 (Gu et al., 2013). We detected elevated transcript levels of FT and SOC1 in fio1 mutants, consistent with the idea that floral repressor activity had been compromised. Splicing changes were present in transcripts encoding circadian modulators such as LHY and WNK1 (Kumar et al., 2011; Mizoguchi et al., 2002), which may contribute to the lengthening of the circadian period observed in fio1 mutants (Kim et al., 2008). Finally, we found that aberrant splicing limits the functional expression of AtSAR1 in fio1 mutants. AtSAR1 encodes a nucleoporin that controls CO protein abundance (Dong et al., 2006; Li et al., 2020; Parry et al., 2006) and FLC expression (Jung et al., 2013). As a result, sar1 mutants are early flowering due to increased CO protein levels (Li et al., 2020) and reduced FLC expression (Jung et al., 2013; Li et al., 2020). Overall, these findings suggest that the early flowering phenotype of fio1 results from splicing changes that reduce the activity of floral repressors and increase the activity of factors that promote flowering, such as CO.
FIO1-dependent m6A modification of U6 snRNA determines splicing accuracy and efficiency
Our data indicate that the major effect of FIO1 occurs through m6A modification of U6 snRNA and the subsequent interaction of m6A-modified U6 snRNA with target 5’SSs. Recent suggestions that FIO1 does not impact splicing (Wang et al., 2022; Xu et al., 2022) can be explained by underpowered experiments and unsuitable methodological approaches. Our reanalysis of recently published RNA-seq and nanopore DRS data in independent fio1 mutant alleles confirms that splicing is indeed disrupted, consistent with the conserved role of FIO1 in methylating U6 snRNA (Wang et al., 2022; Xu et al., 2022; Sun et al., 2022).
We cannot rule out the possibility that FIO1 directly targets other RNA species not detected by our approaches. However, our nanopore DRS analyses did not reveal widespread FIO1-dependent m6A sites in poly(A)+mRNA. We found that 8.7% of FIP37-dependent m6A sites had slightly altered modification rates in fio1-1, suggesting that loss of FIO1 might indirectly affect m6A sites written by the METTL3-like complex. Consistent with this, splicing of MTB (METTL14) RNA is also defective in fio1 mutants, potentially compromising MTB function. Notably, although data from three recent studies (Wang et al., 2022; Xu et al., 2022; Sun et al., 2022) show almost no overlap in respect of RNAs putatively modified by FIO1-dependent addition of m6A (Sun et al., 2022), we find that all these datasets confirm the impact of FIO1 on splice site selection.
We can explain the vast majority of fio1-sensitive splicing events using simple rules governed by sequence features of 5’SSs. When U6 snRNA is m6Amodified, 5’SSs with A+4 are favoured. In the absence of U6 snRNA m6A, 5’SSs with U+4 and/or stronger interactions with U5 snRNA loop 1 are favoured. Cryo-electron microscopy (cryo-EM) structures of the human spliceosome suggest U6 snRNA m6A faces 5’SS+4 (Bertram et al., 2017). Remarkably, we can separate almost 80% of FIO1-dependent 5’SS choices by the identity of the base at the +4 position of 5’SSs alone, suggesting that many of the splicing changes we detect are direct. In contrast, no obvious difference in splice site sequences was associated with temperature-dependent alternative splicing detected in our study.
Our findings are similar to those recently reported for S. pombe, where disruption of the FIO1 orthologue MTL16 results in widespread intron retention (Ishigami et al., 2021). A more diverse set of splicing events is disrupted in Arabidopsis fio1 mutants. Nevertheless, in both species, most changes in splicing can be explained by the impact of m6A modification on U6 snRNA recognition of the 5’SS+4 position or the relative strength of U5 and U6 snRNA interactions at 5’SSs. Defective splicing of specific introns in S. pombe mtl16Δ strains can be experimentally rescued by expression of mutated U5 snRNAs designed to strengthen U5 loop 1 interactions with the upstream exon (Ishigami et al., 2021). Together, these findings reveal the impact of U6 snRNA m6A modification and indicate that U5 and U6 snRNAs have cooperative and compensatory roles in 5’SS selection dependent on 5’SS sequence (Figure 8A).
Two major classes of 5’SS linking cooperative and compensatory roles for U5 and U6 snRNA in 5’SS selection
Our analyses of annotated Arabidopsis 5’SSs identified anti-correlated biases in sequence composition at U5 and U6 snRNA interacting positions. We found the same to be true for 5’SSs in other eukaryotes, including humans. Compensatory patterns of base composition at human 5’SSs have been described (Burge and Karlin, 1997; Carmel et al., 2004; Sibley et al., 2016; Wong et al., 2018; Artemyeva-Isman and Porter, 2021). For example, a ‘seesaw linkage’ pattern was observed where −1 G permits any nucleotide at position +5, and conversely,+5 G permits any nucleotide at position −1. Our findings suggest that such compensatory base composition at 5’SSs can be accounted for by two major classes of 5’SS recognised mainly by interactions with either U5 snRNA loop 1 or the U6 snRNA ACAGA box. Furthermore, the cooperativity of U5 and U6 snRNA function in 5’SS selection may explain why having Gs at both −1 and +5 positions is highly preferential for efficient splicing (Wong et al., 2018). Compensatory interactions are likely to facilitate degeneracy in sequences that can be recognised as 5’SSs. Degeneracy may act as a buffer against deleterious mutations in splice site sequences, as well as lowering barriers to the evolution of new splicing structures. This is important because alternative splicing is clearly associated with developmental complexity (Bush et al., 2017; Nilsen and Graveley, 2010).
What is the role of U6 snRNA m6A modification in 5’ and 3’SS selection?
Our analysis of available cryo-EM structures of human spliceosomes (Bertram et al., 2017; Bertram et al., 2020; Fica et al., 2019; Zhang et al., 2018) does not provide clear evidence of a direct interaction between U6 snRNA m6A and a spliceosomal protein (Figure 8B; Figure 8—figure supplement 1). Instead, during the transfer of the 5’SS from U1 snRNA to U6 and U5 snRNAs in spliceosomal B complexes, U6 ACm6AGA m6A faces 5’SS A+4 in a trans Hoogsteen sugar edge interaction that could stabilise the 5’SS/U6 helix by capping (Figure 8B). This is important because 5’SS A+3 is not aligned to stack on A+4 at this stage. However, later in spliceosomal Bact and C complexes, 5’SS A+3 engages in a more stable interaction with the U6 snRNA G44 adjacent to m6A (ACm6AGA), stacking on top of the 5’SS/U6 helix and stabilising it (Figure 8C; Figure 8—figure supplement 1). Such stabilisation is important because degenerate 5’SS sequences mean that the helix formed with U6 snRNA ACAGA is relatively weak, short, and comprised mostly of non-canonical RNA-RNA interactions (Figure 8B).
Our global RNA-seq analyses reveal that m6A-modified U6 snRNA pairs less readily with 5’SSs containing U+4 than those containing A+4. Biophysical data from model RNAs provide a possible explanation for this finding: an m6A-U base pair can form in a duplex, but the methylamino group rotates from a syn geometry on the Watson–Crick face to a higher-energy anti-conformation, positioning the methyl group in the major groove (Roost et al., 2015). As a result, m6A has a destabilising effect on A−U basepairs in short RNA helices (Kierzek and Kierzek, 2003; Roost et al., 2015). Therefore, 5’SS U+4 may be selected less frequently by m6A-modified U6 snRNA because of destabilisation of an already weak helix. In contrast, the thermal stability of m6A:A is increased compared to A:A (Roost et al., 2015). Biophysical measurements also indicate m6A in unpaired positions base stacks more strongly than the unmodified base, adding substantial stabilisation to adjacent duplexes (Roost et al., 2015). In the absence of U6 snRNA m6A modification, our data reveal that 5’SS U+4 sites are preferred. This is consistent with a relatively strong interaction being important at the terminal position of the U6 snRNA ACAGA helix in the B complex to stabilise its formation.
The potential role of U6 snRNA m6A in 3’SS selection may be explained by U6 ACAGA interactions during the remodelling of the spliceosome for the second splicing reaction (Figure 8C; Figure 8—figure supplement 1). In humans, RNA rearrangements during the C to C* transition include stacking of the 3’SS G−1 onto U6 A45 (ACAGA), which remains paired to U+2 of the 5’SS (Fica, 2020; Fica et al., 2019; Wilkinson et al., 2017). This results in a continuous helix stack, involving the interaction between U6 A43 and 5’SS A+4, which forms the receptor onto which the 3’SS docks. Our data are consistent with a model in which a strong interaction between U6 and the 5’SS +4 position stabilises this receptor, enabling 3’SSs distal to the branchpoint to compete more efficiently with proximal 3’SSs that are favoured by scanning (Smith et al., 1993). A more stable receptor could allow the ATPase Prp22, which proofreads exon ligation (Mayas et al., 2006), to promote more efficient sampling and usage of distal 3’SSs. This has previously been observed in S. cerevisiae (Semlow et al., 2016), where the U6/5’SS helix is intrinsically stronger than in Arabidopsis and humans (Plaschka et al., 2019; Wilkinson et al., 2020).
A regulatory or adaptive role for U6 snRNA m6A modification?
Alternative splicing is the result of competition between multiple splicing choices. Our analysis of annotated alternative 5’SS pairs in Arabidopsis (and humans) revealed that when one of the alternative 5’SS pair exhibited strong U5 snRNA loop 1 recognition, the other was more likely to have stronger complementarity to U6 snRNA ACAGA box. This suggests that alternative splicing could be regulated by changing the relative favourability of U5 and U6 snRNA interactions.
It is an open question as to whether m6A modification of U6 snRNA is regulated. Either the activity of FIO1 might be controlled or demethylases might act directly upon U6 snRNA. There is a precedent for the control of splicing by regulation of U6 snRNA modifications because pseudouridylation of S. cerevisiae U6 snRNA, which controls entry into filamentous growth, can increase the splicing efficiency of suboptimal introns (Basak and Query, 2014). Human UsnRNAs are targeted by demethylases. For example, FTO targets cap adjacent methylation (m6Am) in UsnRNAs (except U6), and this activity may account for splicing differences detected in FTO knockout backgrounds (Mauer et al., 2019). Notably, Arabidopsis mutants defective in the function of the RNA demethylase ALKBH10b flower late (Duan et al., 2017).
Two genome-wide association studies have identified sequence differences at or near the FIO1 locus associated with natural variation in Arabidopsis ecotypes’ flowering time (Price et al., 2020; Sasaki et al., 2015). It will be important to test whether genetic variation alters the efficiency of U6 snRNA m6A modification by FIO1 in different ecotypes and what impact this has upon global splicing patterns, including for the regulators of flowering time that we characterise here. There are 13 genes encoding U6 snRNA in the Arabidopsis Col-0 genome (Wang and Brendel, 2004). Currently, we know little about the relative patterns of U6 snRNA gene expression or m6A modification status and how this may vary in different ecotypes. Indeed, this paucity of knowledge on U6 snRNA variants’ modification and expression applies to humans too. Our study on U6 snRNA m6A modification has important implications for understanding the mechanism of splicing and the evolution of alternative splicing. It will now be interesting to investigate the possibility that modulation of U6 snRNA m6A may be regulatory or adaptive.
Materials and methods
Plant material
Request a detailed protocolThe wild-type Col-0 accession, fio1-3 (SALK_084201), and maf2 (SALK_045623) were obtained from Nottingham Arabidopsis Stock Centre. The fio1-1 mutant was a gift from Prof Hong Gil Nam (Daegu Gyeongbuk Institute of Science and Technology), Republic of Korea. The fip37-4 mutant was a gift from Prof Rupert Fray (University of Nottingham), UK.
Plant growth conditions
Request a detailed protocolSeeds of wild-type Col-0, fip37-4, and fio1-1 used for nanopore DRS and m6A immunopurification were surface sterilised and sown on MS10 media plates supplemented with 2% agar, stratified at 4°C for 2 days, germinated in a controlled environment at 20°C under 16 hr light/8 hr dark conditions, and harvested 14 days after transfer to 20°C. Seeds of wild-type Col-0 and fio1-3 used for Illumina RNA-seq were surface sterilised and sown on 0.5 MS media plates supplemented with 2% agar, stratified at 4°C for two days, germinated in a controlled environment at 20°C under 16 hr light/8 hr dark conditions for 8 days. Seedlings were then transferred to either 28, 20, 12, or 4 °C for 4 hr under dark conditions before harvesting. Seeds used for MAF2 splicing analysis were sown on 0.5 MS media plates supplemented with 2% agar, stratified at 4°C for two days, germinated in a controlled environment at 20°C under 16 hr light/8 hr dark conditions, and harvested 14 days after transfer to 20°C.
Mutant screen
Request a detailed protocolThe two-step EMS mutant screen was conducted in a Col-0 line carrying a homozygous transgene in which the genomic MAF2 coding region was translationally fused to luciferase (gMAF2:LUC). A gateway cloning approach was used to introduce gMAF2:LUC into the Alligator vector pFP101, which features a recombination site downstream of a CaMV 35 S promoter and selection by Green Fluorescent Protein (GFP) fluorescence (Bensmihen et al., 2004). Agrobacterium tumefaciens strain GV3101 was used to transform Col-0 via the floral dip method (Clough and Bent, 1998). A homozygous line with strong LUC expression and a clear GFP fluorescence in the seed coat (Alligator selective marker) was identified and used for EMS mutagenesis.
Approximately, 20,000 M0 seeds of the gMAF2:LUC line were soaked overnight in 100 mL of phosphate buffer at 4°C. The buffer was replaced, with EMS added to a final concentration of 25 mM, and the seeds were incubated at room temperature with gentle agitation for 16 hr. The EMS was neutralised with 1 M NaOH, and the seeds were gently washed twice with 100 mM sodium thiosulphate for 15 min, followed by three 15 min washes in distilled water. The M1 seeds were air-dried on filter paper overnight, before planting on soil. The resulting plants were allowed to self-pollinate and grow to maturity. M2 seeds were collected from all plants. Approximately, 15,000 seeds were screened in two steps. In step 1, stratified seeds were planted on soil and grown in controlled environment chambers at 16°C, 16 hr light/8 hr dark, with light intensity of 200 µmol/m2/s. The first 100 plants to flower were selected for further screening. Splicing of MAF2 intron 3 was assessed by RT-qPCR to identify early flowering lines showing enhanced intron retention at 16°C. Although the parental line had the gMAF2:LUC reporter, luciferase activity was not monitored at any stage. Leaf material was frozen and homogenized using QIAGEN TissueLyser LT. RNA was extracted using a Nucleospin II RNA extraction kit (Machery–Nagel). Total RNA (1.5 µg per sample) was reverse-transcribed using the High Capacity cDNA Reverse Transcription Kit with Rnase Inhibitor (Applied Biosystems). Amplification of MAF2 was performed according to manufacturer’s guidelines for 35 PCR cycles using primers 3 and 4 from a previous study (Rosloski et al., 2013). These primers, which span MAF2 intron 3, generate two products, which correspond to MAF2 var2 (MAF2 intron 3 retained) and MAF2 var1 (MAF2 intron 3 excised), respectively. PCR amplified products were separated on a 2% w/v agarose gel to resolve the splice variants. M2 plants showing a significantly enhanced var2:var1 ratio compared to the parental line control were selected for further testing.
Bulk segregant analysis for mutant mapping
Request a detailed protocolThe EMS-129 line was backcrossed to the parental gMAF2:LUC reporter line to generate an F2 segregating population. Phenotypically early flowering and wild-type flowering time plants were identified from this F2 population and separated for bulk DNA sequencing. Genomic DNA was extracted using a Qiagen Plant DNA Maxi kit. Sequencing was carried out by the University of Leeds Next Generation Sequencing (NGS) facility using a HiSeq3000 sequencer with a 150 bp paired end library.
Sequenced reads were mapped to the TAIR10 reference sequence using bwa-mem2 (Vasimuddin et al., 2019). SNPs were called using bcftools and filtered for EMS-induced G>A transitions only (Danecek et al., 2021). Functional consequences of each SNP were predicted using SnpEff (Cingolani et al., 2012). Allele fraction plots were generated from vcf files using matplotlib. Significance testing was performed for each SNP using G-tests and smoothed using a tri-cube kernel with a window size of 2 megabases.
Nanopore DRS
Total RNA isolation
Request a detailed protocolTotal RNA was isolated using Rneasy Plant Mini kit (Qiagen) and treated with TURBO Dnase (ThermoFisher Scientific). The total RNA concentration was measured using a Qubit 1.0 Fluorometer and Qubit RNA BR Assay Kit (ThermoFisher Scientific). RNA quality and integrity were assessed using the NanoDrop 2000 spectrophotometer (ThermoFisher Scientific) and Agilent 2200 TapeStation System (Agilent).
Preparation of libraries for direct RNA sequencing of poly(A)+ mRNA using nanopores
Request a detailed protocolTotal RNA was isolated from the Col-0, fip37-4, and fio1-1 seedlings as detailed above. mRNA isolation and preparation of nanopore DRS libraries using the SQK-RNA002 nanopore DRS Kit (Oxford Nanopore Technologies) were performed as previously described (Parker et al., 2020). Libraries were loaded onto R9.4 SpotON Flow Cells (Oxford Nanopore Technologies) and sequenced using a minION device for a 48 hr runtime. Four biological replicates were performed for each genotype.
Nanopore DRS mapping
Request a detailed protocolNanopore DRS reads were basecalled using Guppy version 3.6.0 high accuracy RNA model. For mRNA modification analysis, reads were mapped to the Araport11 reference transcriptome (Cheng et al., 2017) using minimap2 version 2.17 (Li, 2018), with parameters -a -L –cs =short k14 –for-only –secondary =no. For other analyses, reads were mapped to the Arabidopsis TAIR10 genome (Arabidopsis Genome Initiative, 2000) using two pass alignment with minimap2 version 2.17 and 2passtools version 0.3 (Parker et al., 2021b). First pass minimap2 alignment was performed using the parameters -a -L –cs =short -x splice -G20000 –end-seed-pen 12 -uf. 2passtools score was then run with default parameters on each replicate to generate junctions, followed by 2passtools merge to combine them into a final set of guide junctions. Reads were remapped with minimap2 using the same parameters but with the addition of the guide junctions using –junc-bed and –juncbonus =10. Pipelines for processing of nanopore DRS data were built and executed using Snakemake version 6.15.3 (Köster and Rahmann, 2012). Public nanopore DRS data downloaded from ENA accessions PRJNA749003 (Xu et al., 2022) and PRJNA877932 (Sun et al., 2022) were processed in the same way as newly generated data.
Nanopore poly(A)+ mRNA modification analysis
Request a detailed protocolDifferential modification analysis was performed on Col-0, fip37-4, and fio1-1 data using the ‘n-sample’ GitHub branch of Yanocomp (Parker et al., 2021a). Kmer-level signal data were generated using f5c eventalign version 0.13.2 (Gamaarachchi et al., 2020; Loman et al., 2015) and Yanocomp prep. A three-way comparison between the genotypes was performed using Yanocomp gmmtest, with a minimum KS statistic of 0.25. A 5% false discovery rate threshold was used to identify transcriptomic sites with significant changes in modification rate. Motif analysis was performed using meme version 5.1.1 with the parameters -cons NNANN -minw 5 -maxw 5 -mod oops (Bailey et al., 2015). Differential poly(A) site usage was performed using d3pendr version 0.1 with default parameters and thresholded using a 5% false discovery rate and an effect size (measured using earth mover distance) of 25 (Parker et al., 2021c).
Illumina RNA sequencing
Preparation of libraries for Illumina RNA sequencing
Request a detailed protocolTotal RNA was isolated from Col-0 and fio1-3 seedlings using Nucleospin RNA kit (Macherey–Nagel, 740955) and treated with rDNase (Macherey–Nagel, CAS 9003-98-9) on columns according to manufacturers’ instructions. RNA concentration, quality, and integrity were assessed using the NanoDrop 1000 spectrophotometer (Labtech) and agarose gel electrophoresis. Poly(A)+RNA purification and Illumina RNA-seq library preparation were performed by Genewiz UK Ltd. Poly(A)+RNA was selected with NEBNext Poly(A) mRNA Magnetic Isolation Module. Preparation of the sequencing libraries was performed using the NEBNext Ultra II Directional RNA Library Prep Kit for Illumina (New England Biolabs). 150 bp paired-end sequencing was carried out using Illumina Novaseq 6000. Six biological replicates were performed for each genotype.
Illumina RNA sequencing data processing
Request a detailed protocolIllumina RNA-seq data were assessed for quality using FastQC version 0.11.9 and MultiQC version 1.8 (Andrews, 2017; Ewels et al., 2016). Reads were mapped to the TAIR10 genome using STAR version 2.7.3 a (Dobin et al., 2013) with a splice junction database generated from the Araport11 reference annotation (Cheng et al., 2017). We used Stringtie version 2.1.7 in mix mode (Kovaka et al., 2019) to generate condition-specific transcriptome assemblies from Illumina RNA-seq replicates and pooled Nanopore DRS data. All assemblies were merged with the Araport11 reference annotation using Stringtie merge to create a unified set of transcripts for quantification. Transcript open reading frames were annotated using Transuite version 0.2.2 (Entizne et al., 2020), which was also used to predict nonsense-mediated decay sensitivity. Pipelines were written and executed using Snakemake version 6.15.3 (Köster and Rahmann, 2012).
Transcripts were quantified using Salmon version 1.1.0 (Patro et al., 2017) with the TAIR10 genome assembly as decoys. SUPPA version 2.3 was used to generate event level PSIs from transcript level quantifications (Trincado et al., 2018). PSIs were loaded into Python 3.6.7 using pandas version 1.0, and generalised linear models (GLMs) were fitted per event using statsmodels version 0.11 (Harris et al., 2020; McKinney, 2010; Oliphant, 2007; Seabold and Perktold, 2010). GLMs were used to test the relationship of PSI with genotype, temperature, and genotype×temperature interaction. Calculated p values were adjusted for multiple testing using the Benjamini–Hochberg false discovery rate method. A false discovery rate of 5% was chosen to threshold events with significant changes in PSI which correlated with genotype, temperature, or genotype×temperature. Changes in motif composition were tested using G-tests of base frequencies at the −2 to +5 position of the splice site. To generate sequence logos, 5’SS or 3’SS from sets of alternative splicing events was filtered to remove duplicated positions, and probability logos were plotted using matplotlib version 3.3 and matplotlib_logo (Hunter, 2007; Parker, 2022b). Contingency tables of splice site classes at U5 and U6 interacting positions were generated using the difference of the −2 to −1 positions of the 5’SS from the consensus motif AG, and the difference of the +3 to +5 positions of the 5’SS from the consensus motif RAG, respectively. Classes were ordered by their log-likelihoods using PSSMs generated from all annotated 5’SSs (described below). Heatmaps of contingency tables were generated with seaborn version 0.11 (Waskom, 2021). Gene tracks using reads aligned to the TAIR10 reference genome were generated using pyBigWig version 0.3.17, pysam version 0.18, and matplotlib version 3.3 (Heger et al., 2014; Hunter, 2007; Ramírez et al., 2014).
Public Illumina RNA-seq datasets were downloaded from GSA accession CRA004052 (Wang et al., 2022) and ENA accession PRJNA877932 (Sun et al., 2022). Quality control was performed using FastQC and MultiQC as described above (Andrews, 2017; Ewels et al., 2016). Col-0 replicate 1 of dataset CRA004052 was discarded due to extreme GC-bias, likely caused by PCR overamplification and/or low RNA input. Transcript-level quantification and PSI estimation were performed as described above. PSIs were used to estimate splicing changes in fio1 mutant alleles using SUPPA2 diffSplice (Trincado et al., 2018).
Global splice site analyses
Request a detailed protocolTo measure the predicted strengths of U5 snRNA loop 1 and U6 snRNA ACAGA box interactions for individual sequences, all 5’SSs annotated in the Araport11 reference transcriptome (Cheng et al., 2017) were used to generate log transformed PSSMs. The U5 and U6 log-likelihood scores for individual sequences were then calculated using the −2 to −1 and +3 to +5 positions, respectively. Correlation of U5 and U6 PSSM scores was calculated using Spearman rank correlation coefficient. This analysis was repeated using 5’SSs from the H. sapiens GRCh38 (International Human Genome Sequencing Consortium, 2004), C. elegans Wbcel235 (The C. elegans Sequencing Consortium, 1998), D. melanogaster BDGP6 (dos Santos et al., 2015), and D. rerio GRCz11/danRer11 (Howe et al., 2013) genome assemblies and annotations, downloaded from Ensembl Genomes release 104 (Howe et al., 2021). To measure relative U5 and U6 interactions for pairs of alternative 5’SSs, SUPPA was used to identify pairs of alternative 5’SSs in each genome annotation (Trincado et al., 2018). These alternative 5’SS pairs were ordered by their genomic positions relative to the strand of the parent gene (i.e. upstream and downstream), and log-odds ratios were calculated using downstream 5’SS log-likelihoods as the denominator. Correlation of U5 and U6 log-odds ratios was calculated using Spearman rank correlation coefficient. Scatter plots with linear regression lines were plotted using seaborn (Waskom, 2021), and 95% CIs for regression lines were calculated using bootstrap sampling with replacement.
Immunopurification and detection of m6A modified RNA
Synaptic systems anti-m6A
Request a detailed protocolTotal RNA was purified from ~300 mg of frozen plant tissue using the miRVana miRNA isolation kit (Ambion) and treated with DNase I (New England Biolabs) according to the manufacturer’s instructions. The quantity and integrity of RNA were checked using a NanoDrop 2000 spectrophotometer (Thermo Fisher Scientific) and Agilent 2200 TapeStation System. Approximately 5 μg of RNA was suspended in 500 μL low-salt buffer (50 mM Tris-HCl pH 7.4, 150 mM NaCl, and 0.5% v/v NP-40), with 80 U RNAsin Plus RNase inhibitor, (Promega) and 10 μL (1 mg/mL) m6A-rabbit polyclonal purified antibody (#202 003 Synaptic Systems). The samples were mixed by rotation for 2 hr at 4°C. Protein A/G magnetic beads (Pierce, ThermoFisher) were washed twice with low-salt buffer, then added to each RNA/antibody sample and mixed with rotation at 4°C for 16 hr. The beads were washed twice with high-salt buffer (50 mM Tris-HCl pH 7.4, 0.5 M NaCl, and 1% v/v NP-40), twice with low-salt buffer, and twice with Polynucleotide kinase (PNK) wash buffer (20 mM Tris-HCl pH 7.4, 10 mM MgCl2, and 0.2% v/v Tween-20). Immunopurified RNA was eluted by digestion with proteinase K in 200 μL of proteinase K (PK) buffer (50 mM Tris-HCl pH 7.4, 150 mM NaCl, 1 mM EDTA, and 0.1% w/v SDS) at 37°C for 20 min with shaking at 300 rpm, followed by phenol-chloroform extraction and sodium acetate/ethanol precipitation.
Immunopurifed RNA (45 ng) was reverse transcribed (RT) with U6 and U2 snRNA reverse primers using SuperScript III Reverse Transcriptase (ThermoFisher), according to the manufacturer’s instructions. RT-qPCR was carried out using the SYBR Green I (Qiagen) mix with primers targeted to U2 and U6 snRNA (Supplementary file 2). The specificity of RT-qPCR was confirmed by sequencing of amplified products.
Millipore anti-m6A
Request a detailed protocolTotal RNA was purified from ~600 mg of frozen plant tissue using the miRVana miRNA isolation kit (Ambion) and treated with DNase I (New England Biolabs), according to the manufacturer’s instructions. The quantity and integrity of RNA were checked using a NanoDrop 1000 spectrophotometer (Thermo Fisher Scientific). After sodium acetate/ethanol precipitation, approximately 10 μg of RNA was resuspended in 20 μL nuclease free water, with 5% of the RNA being kept for the input sample.
Protein A/G magnetic beads (Millipore, 16–663) were washed twice with m6A wash buffer (10 mM TrisHCl pH7.4, 150 mM NaCl, and 0.1% v/v NP-40) and then resuspended in 50 μL m6A wash buffer. 5 μL (1 mg/mL) of Anti-N6-methyladenosine (m6A) antibody (Millipore, ABE572) or Rabbit Anti-IgG (Thermofisher, A16104) was coupled to the washed beads on a roller for 40 min at room temperature. After incubation, beads were washed three times in m6A wash buffer. m6A immunoprecipitation buffer was then added to the beads (10 mM TrisHcl pH7.4, 150 mM NaCl, 0.1% v/v NP-40, and 17.5 mM EDTA pH 8) with 80 U SUPERase•In RNase Inhibitor, and 10 μL of each RNA sample was added to both the anti-m6A and anti-IgG antibody/beads mixture. Samples were mixed by rotation at 4°C for 16 hr. The beads were washed five times with m6A wash buffer. Immunopurified and 5% input RNA was eluted by digestion with proteinase K in 150 μL of PK buffer (50 mM Tris-HCl pH 7.4, 150 mM NaCl, 0.1% v/v NP-40, and 0.1% w/v SDS) at 55°C for 30 min, followed by extraction with TRIzol LS (10296028, Thermofisher) and sodium acetate/ethanol precipitation.
Immunopurifed RNA (500 ng) was RT with random hexamers using ThermoFisher MultiscribeTM II Reverse Transcriptase (4311235, Thermo Fisher Scientific) according to the manufacturer’s instructions. RT-qPCR was carried out using the SYBR Green I (Qiagen) mix with primers targeted to U2 and U6 snRNA (Supplementary file 2). The specificity of RT-qPCR was confirmed by sequencing of amplified products.
m6A liquid chromatography with tandem mass spectrometry
Request a detailed protocolTotal RNA was isolated as described above. Poly(A)+RNA was purified using the Dynabeads mRNA purification kit (Thermo Fisher Scientific) according to the manufacturer’s instructions. The quality and quantity of mRNA were assessed using a NanoDrop 2000 spectrophotometer and Agilent 2200 TapeStation System. Samples for m6A LC-MS/MS were prepared as previously described (Huang et al., 2018) with modifications according to Parker et al., 2021c. LC-MS/MS was carried out by the FingerPrints Proteomics facility at the University of Dundee. m6A/A ratio quantification was performed in comparison with the curves obtained from pure adenosine (endogenous P1 receptor agonist, Abcam) and m6A (modified adenosine analog, Abcam) nucleoside standards. Statistical analysis was performed using a two-way t-test.
Code availability
Request a detailed protocolAll pipelines, scripts, and notebooks used to generate figures are available from GitHub at https://github.com/bartongroup/Simpson_Davies_Barton_U6_methylation and Zenodo at https://zenodo.org/record/6372644 (Parker, 2022a).
Data availability
Illumina sequencing data from the genetic screen that identified fio1-4 is available from ENA accession PRJEB51468. Col-0, fip37-4 and fio1-1 nanopore DRS data is available from ENA accession PRJEB51364. Col-0 and fio1-3 Illumina RNA-Seq data is available from ENA accession PRJEB51363.
-
European Nucleotide Archive (EMBL-EBI)ID PRJEB51468. Mutant screen of early flowering Arabidopsis EMS mutants with increased MAF2 intron retention.
-
European Nucleotide Archive (EMBL-EBI)ID PRJEB51364. Nanopore direct RNA sequencing of Col-0 and fio1-3 mutant Arabidopsis.
-
European Nucleotide Archive (EMBL-EBI)ID PRJEB51363. Illumina RNA sequencing of Col-0 and fio1-3 mutant Arabidopsis.
-
European Nucleotide Archive (EMBL-EBI)ID PRJNA749003. N6-Adenosine (m6A) Methylation profiling in Arabidopsis thaliana wild-type and fio1-2 by nanopore DRS sequencing.
-
NCBI Gene Expression OmnibusID GSE171928. Next Generation Sequencing of Wild Type (Col-0), fio1-1 and fio1-2 methylated RNA immunoprecipitation and transcriptome in Arabidopsis thaliana.
-
Genome Sequence Archive (CNCB)ID CRA004052. Identification and Function Research of Methyltransferase of N6 Methyladenine in Plants.
References
-
The genetic basis of flowering responses to seasonal cuesNature Reviews. Genetics 13:627–639.https://doi.org/10.1038/nrg3291
-
SoftwareFastQC: A quality control tool for high throughput sequence data, version v3Babraham Bioinformatics.
-
Mechanistic insights into M6a modification of U6 snRNA by human METTL16Nucleic Acids Research 48:5157–5168.https://doi.org/10.1093/nar/gkaa227
-
U5 snRNA interactions with exons ensure splicing precisionFrontiers in Genetics 12:676971.https://doi.org/10.3389/fgene.2021.676971
-
Prediction of complete gene structures in human genomic DNAJournal of Molecular Biology 268:78–94.https://doi.org/10.1006/jmbi.1997.0951
-
Alternative splicing and the evolution of phenotypic noveltyPhilosophical Transactions of the Royal Society of London. Series B, Biological Sciences 372:20150474.https://doi.org/10.1098/rstb.2015.0474
-
Contribution of major FLM isoforms to temperature-dependent flowering in Arabidopsis thalianaJournal of Experimental Botany 68:5117–5127.https://doi.org/10.1093/jxb/erx328
-
Correcting for differential transcript coverage reveals a strong relationship between alternative splicing and organism complexityMolecular Biology and Evolution 31:1402–1413.https://doi.org/10.1093/molbev/msu083
-
Araport11: a complete reannotation of the Arabidopsis thaliana reference genomeThe Plant Journal 89:789–804.https://doi.org/10.1111/tpj.13415
-
Twelve years of samtools and bcftoolsGigaScience 10:giab008.https://doi.org/10.1093/gigascience/giab008
-
Star: ultrafast universal RNA-seq alignerBioinformatics 29:15–21.https://doi.org/10.1093/bioinformatics/bts635
-
A putative Arabidopsis nucleoporin, atnup160, is critical for RNA export and required for plant tolerance to cold stressMolecular and Cellular Biology 26:9533–9543.https://doi.org/10.1128/MCB.01063-06
-
Cryo-Em snapshots of the human spliceosome reveal structural adaptions for splicing regulationCurrent Opinion in Structural Biology 65:139–148.https://doi.org/10.1016/j.sbi.2020.06.018
-
Localization of modified nucleotides in Schizosaccharomyces pombe spliceosomal small nuclear rnas: modified nucleotides are clustered in functionally important regionsRNA 2:909–918.
-
Matplotlib: a 2D graphics environmentComputing in Science & Engineering 9:90–95.https://doi.org/10.1109/MCSE.2007.55
-
Plant small nuclear rnas. II. U6 RNA and a 4.5SI-like RNA are present in plant nucleiNucleic Acids Research 15:543–560.https://doi.org/10.1093/nar/15.2.543
-
Snakemake -- a scalable bioinformatics workflow engineBioinformatics 28:2520–2522.https://doi.org/10.1093/bioinformatics/bts480
-
Rice WNK1 is regulated by abiotic stress and involved in internal circadian rhythmPlant Signaling & Behavior 6:316–320.https://doi.org/10.4161/psb.6.3.13063
-
The non-Watson-Crick base pairs and their associated isostericity matricesNucleic Acids Research 30:3497–3531.https://doi.org/10.1093/nar/gkf481
-
Minimap2: pairwise alignment for nucleotide sequencesBioinformatics 34:3094–3100.https://doi.org/10.1093/bioinformatics/bty191
-
Fto controls reversible m6am RNA methylation during snRNA biogenesisNature Chemical Biology 15:340–347.https://doi.org/10.1038/s41589-019-0231-8
-
Exon ligation is proofread by the DExD/H-box ATPase Prp22pNature Structural & Molecular Biology 13:482–490.https://doi.org/10.1038/nsmb1093
-
Data structures for statistical computing in pythonproceedings of the 9th python in science conferenceAustin, TX pp. 51–56.
-
Spliceosomal snRNA epitranscriptomicsFrontiers in Genetics 12:652129.https://doi.org/10.3389/fgene.2021.652129
-
The intronome of budding yeastsComptes Rendus Biologies 334:662–670.https://doi.org/10.1016/j.crvi.2011.05.015
-
Python for scientific computingComputing in Science & Engineering 9:10–20.https://doi.org/10.1109/MCSE.2007.58
-
Structural basis of nuclear pre-mrna splicing: lessons from yeastCold Spring Harbor Perspectives in Biology 11:a032391.https://doi.org/10.1101/cshperspect.a032391
-
In the presence of population structure: from genomics to candidate genes underlying local adaptationEcology and Evolution 10:1889–1904.https://doi.org/10.1002/ece3.6002
-
DeepTools: a flexible platform for exploring deep-sequencing dataNucleic Acids Research 42:W187–W191.https://doi.org/10.1093/nar/gku365
-
Structure and thermodynamics of N6-methyladenosine in RNA: a spring-loaded base modificationJournal of the American Chemical Society 137:2107–2115.https://doi.org/10.1021/ja513080v
-
ConferenceStatsmodels: Econometric and statistical modeling with pythonProceedings of the 9th Python in Science ConferencePresented at the Python in Science Conference. SciPy.https://doi.org/10.25080/majora-92bf1922-011
-
Lessons from non-canonical splicingNature Reviews. Genetics 17:407–421.https://doi.org/10.1038/nrg.2016.46
-
Scanning and competition between AGS are involved in 3’ splice site selection in mammalian intronsMolecular and Cellular Biology 13:4939–4952.https://doi.org/10.1128/mcb.13.8.4939-4952.1993
-
ConferenceEfficient Architecture-Aware Acceleration of BWA-MEM for Multicore Systems2019 IEEE International Parallel and Distributed Processing Symposium (IPDPS. pp. 314–324.https://doi.org/10.1109/IPDPS.2019.00041
-
Seaborn: statistical data visualizationJournal of Open Source Software 6:3021.https://doi.org/10.21105/joss.03021
-
Rna splicing by the spliceosomeAnnual Review of Biochemistry 89:359–388.https://doi.org/10.1146/annurev-biochem-091719-064225
Article and author information
Author details
Funding
Biotechnology and Biological Sciences Research Council (BB/M000338/1)
- Brendan H Davies
Biotechnology and Biological Sciences Research Council (BB/W007967/1)
- Brendan H Davies
Biotechnology and Biological Sciences Research Council (BB/T007222/1)
- Beth K Soanes
HORIZON EUROPE Marie Sklodowska-Curie Actions (799300)
- Katarzyna Knop
HORIZON EUROPE Marie Sklodowska-Curie Actions (896598)
- Nisha Joy
Wellcome Trust (220212/Z/20/Z)
- Sebastian M Fica
Global Challenges Research Fund (University of Dundee Global Challenges Research Fund)
- Geoffrey J Barton
- Gordon G Simpson
Biotechnology and Biological Sciences Research Council (BB/W007673/1)
- Gordon G Simpson
Biotechnology and Biological Sciences Research Council (BB/V010662/1)
- Geoffrey J Barton
- Gordon G Simpson
Biotechnology and Biological Sciences Research Council (BB/M010066/1)
- Geoffrey J Barton
- Gordon G Simpson
Biotechnology and Biological Sciences Research Council (BB/M004155/1)
- Geoffrey J Barton
- Gordon G Simpson
The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication. For the purpose of Open Access, the authors have applied a CC BY public copyright license to any Author Accepted Manuscript version arising from this submission.
Acknowledgements
We thank Mary McKay (School of Biology, University of Leeds) for excellent technical assistance with the EMS mutant screen. We are indebted to Oliver Manners for help in establishing m6A immunoprecipitation. We thank Abdelmadjid Atrih (Centre for Advanced Scientific Technologies) for the m6A LC-MS/MS analysis. We thank David Lilley, Tim Wilson and Adrian Whitehouse for helpful discussions. We thank Alper Akay, Martin Balcerowicz, James Lloyd, Anjil Srivastava and Carey Metheringham for comments on the manuscript. This work was supported by awards from the BBSRC (BB/V010662/1; BB/M010066/1; BB/M004155/1; BB/W007673/1) and the University of Dundee Global Challenges Research Fund to GGS and GJB; and awards from BBSRC (BB/M000338/1; BB/W007967/1) and ERA-CAPS FLOWPLAST to BD. BS is funded through the BBSRC DTP3 award, BB/T007222/1. KK and NJ were funded through the European Union Horizon 2020 research and innovation programme under Marie Skłodowska-Curie grant agreements 799300 and 896598, respectively. SMF is a Wellcome Trust and Royal Society Sir Henry Dale Fellow (grant number 220212/Z/20/Z). The FingerPrints Proteomics Facility of the University of Dundee is supported by a Wellcome Trust Technology Platform Award (097945/B/11/Z).
Copyright
© 2022, Parker et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.
Metrics
-
- 2,772
- views
-
- 353
- downloads
-
- 29
- citations
Views, downloads and citations are aggregated across all versions of this paper published by eLife.
Download links
Downloads (link to download the article as PDF)
Open citations (links to open the citations from this article in various online reference manager services)
Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)
Further reading
-
- Chromosomes and Gene Expression
- Developmental Biology
About 70% of human cleavage stage embryos show chromosomal mosaicism, falling to 20% in blastocysts. Chromosomally mosaic human blastocysts can implant and lead to healthy new-borns with normal karyotypes. Studies in mouse embryos and human gastruloids showed that aneuploid cells are eliminated from the epiblast by p53-mediated apoptosis while being tolerated in the trophectoderm. These observations suggest a selective loss of aneuploid cells from human embryos, but the underlying mechanisms are not yet fully understood. Here, we investigated the cellular consequences of aneuploidy in a total of 125 human blastocysts. RNA-sequencing of trophectoderm cells showed activated p53 pathway and apoptosis proportionate to the level of chromosomal imbalance. Immunostaining corroborated that aneuploidy triggers proteotoxic stress, autophagy, p53-signaling, and apoptosis independent from DNA damage. Total cell numbers were lower in aneuploid embryos, due to a decline both in trophectoderm and in epiblast/primitive endoderm cell numbers. While lower cell numbers in trophectoderm may be attributed to apoptosis, aneuploidy impaired the second lineage segregation, particularly primitive endoderm formation. This might be reinforced by retention of NANOG. Our findings might explain why fully aneuploid embryos fail to further develop and we hypothesize that the same mechanisms lead to the removal of aneuploid cells from mosaic embryos.
-
- Chromosomes and Gene Expression
- Developmental Biology
Transcription often occurs in bursts as gene promoters switch stochastically between active and inactive states. Enhancers can dictate transcriptional activity in animal development through the modulation of burst frequency, duration, or amplitude. Previous studies observed that different enhancers can achieve a wide range of transcriptional outputs through the same strategies of bursting control. For example, in Berrocal et al., 2020, we showed that despite responding to different transcription factors, all even-skipped enhancers increase transcription by upregulating burst frequency and amplitude while burst duration remains largely constant. These shared bursting strategies suggest that a unified molecular mechanism constraints how enhancers modulate transcriptional output. Alternatively, different enhancers could have converged on the same bursting control strategy because of natural selection favoring one of these particular strategies. To distinguish between these two scenarios, we compared transcriptional bursting between endogenous and ectopic gene expression patterns. Because enhancers act under different regulatory inputs in ectopic patterns, dissimilar bursting control strategies between endogenous and ectopic patterns would suggest that enhancers adapted their bursting strategies to their trans-regulatory environment. Here, we generated ectopic even-skipped transcription patterns in fruit fly embryos and discovered that bursting strategies remain consistent in endogenous and ectopic even-skipped expression. These results provide evidence for a unified molecular mechanism shaping even-skipped bursting strategies and serve as a starting point to uncover the realm of strategies employed by other enhancers.