NOseq: amplicon sequencing evaluation method for RNA m6A sites after chemical deamination

Abstract

Methods for the detection of m6A by RNA-Seq technologies are increasingly sought after. We here present NOseq, a method to detect m6A residues in defined amplicons by virtue of their resistance to chemical deamination, effected by nitrous acid. Partial deamination in NOseq affects all exocyclic amino groups present in nucleobases and thus also changes sequence information. The method uses a mapping algorithm specifically adapted to the sequence degeneration caused by deamination events. Thus, m6A sites with partial modification levels of ∼50% were detected in defined amplicons, and this threshold can be lowered to ∼10% by combination with m6A immunoprecipitation. NOseq faithfully detected known m6A sites in human rRNA, and the long non-coding RNA MALAT1, and positively validated several m6A candidate sites, drawn from miCLIP data with an m6A antibody, in the transcriptome of Drosophila melanogaster. Conceptually related to bisulfite sequencing, NOseq presents a novel amplicon-based sequencing approach for the validation of m6A sites in defined sequences.

Figure 1.

Deamination treatment and deamination sequencing (NOseq) output. (A) Reaction of adenosine, N⁶-methyladenosine (m⁶A) and guanosine under deamination treatment and their corresponding reaction products inosine, N⁶-methyl-N⁶-nitrosoadenosine (NOm⁶A) and xanthosine. (Supplement Figure S1 for cytidine and uridine). (B) Deamination treatment of a 53mer, containing m⁶A at position 33, with 72 different conditions (composed of different reaction temperatures, pH and reaction times) (Supplement Figure S4 for additional information). (C) Partial deamination sequencing output and required alignment adaptations (Example: native adenosine residues gave rise to either A or G signals, the latter deriving from the fraction having undergone deamination to I).

Figure 2.

General analysis pipeline of NOseq data (exemplary for deamination of the 53mer with a treatment of 50°C, pH 4.0 and 20 min). After deamination treatment the treated RNA is used for library preparation and sequenced (for comparison to two-step alignment - see Supplement Figure S17). (A) Reads are aligned to the reference with an asymmetrical substitution matrix (see Supplement Figure S5) to address partial deamination and the corresponding fraction of reads in % is plotted for each position in an alignment plot (A in green, G in orange, C in blue and T in red). (B) Adenosine sites are filtered from the alignment plot. (C) The fraction of reads in % is separated to generate adenosine and guanosine rate plots (separated adenosine and guanosine fraction of reads in % at A sites). (D) The A rate is then used for evaluation of the m⁶A site in a probability plot, which illustrates the probability of an A site being an m⁶A, calculated by the distance of the respective A rate to the average A rate, with a red line at 0.95 showing the detection threshold.

Figure 3.

NOseq data comparison (MeRIP enriched vs. non-MeRIP) to address its detection threshold for sub-stoichiometric m⁶A with an m⁶A content at position 33 (m⁶A/A) of the 53mer that was varied from 0 to 100%. (black dots MeRIP, grey triangles non-MeRIP). Data was averaged from duplicates (error bars show standard deviations). Red line shows detection threshold at 0.95.

Figure 4.

Library preparation scheme for amplicon sequencing of targeted m⁶A sites. The target region is amplified in a two-step process. First, reverse transcription of the deaminated RNA is performed with a target-specific RT primer. The latter contains the UMI and the binding sequence for the universal P7 primer including the barcode for Illumina sequencing. Afterwards, a target specific P5 primer and the P7 primer are used for PCR amplification.

Figure 5.

Adenosine rate (A fraction of reads in % at A sites) and probability plots of three biological amplicons, showing the progressed A-to-I conversion and corresponding m⁶A site evaluation (red line shows detection threshold at 0.95). (A) m⁶A at position 1832 in human 18S rRNA (HEK293T cells). (B) m⁶A at position 2577 in human lncRNA MALAT1 (HeLa cells). (C) m⁶A at position 349 in Drosophila melanogaster H mRNA after MeRIP enrichment. (Supplement Figures S14 & 15 for additional information, including non-treated references).