Given an experiment consisting of an input (baseline RNA) and IP (pulldown to find RNAs bound to certain protein of interest)... Is a DE analysis performed over the RNA-seq data from the samples (lets say with EdgeR or DESEQ2) suitable to reveal the preferentially bound RNAs? What other software tools would you recommend?
3 Answers
You need to be careful of terminology. To me, a RIP-seq experiment involves a pull down, followed by a RNAseq library prep. Thus, the whole transcript is captured, not just the binding site of the protein (as in CLIP-Seq, HITS-CLIP, PAR-CLIP, iCLIP or eCLIP). Thus "peak-callers" whether they be designed for calling protein-DNA or protein-RNA binding sites are not suitable as the signal will not be punctuate or peaky. Of the methods mentioned by @Devon Ryan, only RIPSeeker seems setup to deal with this sort of data.
In some of the papers that talk about calling RNA-protein interactions they lump RIP-seq in with the CLIP techniques (see piranha and ASPeak papers referenced by @Devon Ryan). They appear to be talking about a technique where the RNA is fragmented before it is pulled down. Thus you would capture peaky binding sites (actually, the protocol would be remarkably similar to CLIP).
You should work out which of these you have.
As for using DESeq/EdgeR etc, in the RIPSeeker paper they say:
Furthermore, programs for de novo transcript assembly followed by differential expression (DE) analysis, such as the Cufflinks/Cuffdiff suite (15,16), and for DE on a set of known transcripts, such as DESeq (17), may appear applicable to RIP-seq analysis. Unlike peak-calling strategy, however, the transcript-based methods assume the full transcriptome being sequenced at a fairly deep coverage (as usually the case in RNA-seq) and thus may be sensitive to background noise typical to the IP-based protocols, which is due to both the non-specific RNA interactions with a protein of interest and the non-specific RNA input from the pull-down of the (mutant) control (Supplementary Figures S3 and S4).
I don't know to what extent RIPSeeker performs better than the naive DESeq approach as its not included as one of the comparators in the paper.
-
$\begingroup$ How important is peakiness though? Macs2 seems to deal well with broad peaks. Couldn’t it be used here? $\endgroup$ Aug 25, 2017 at 12:49
-
$\begingroup$ There are number of reasons not to use a protein-DNA peak caller, many detailed in the RIPSeeker paper. Many of these are to do with splicing: fragment sizes estimated from mapped RIPseq will not reflect genuine fragment sizes, using either extended reads or fragment midpoints will likely leave you in the middle of an intron. MACS2 broadpeak works by merging narrow peaks connected by sub-significant sequence. Its not clear this is a good description of the data you'd get in RIP. The RIPSkeeker paper compares to Macs, and looks like it out-performs by quite a way to me. $\endgroup$ Aug 25, 2017 at 13:04
You can try doing standard differential expression, but I worry that the between-sample normalization will work poorly. Personally, I would do peak calling instead, followed by diffBind. You have a few tools to choose from when it comes to this. In the past, I've rolled my own methods for this using MACS2 and genomic alignments (I then converted those to bedGraph files where entries are transcripts, so peak calling is then in transcript-coordinates).
It's probably easier, though, to use premade software. There are a number of packages out there for RIP-seq. These include piranha, RIPSeeker and ASPeak. Particularly if you're just starting out, you'd be well served with using one of those tools.
I've got the RIP-seq analysis outputs. We used an isotype control antibody VS antibody against an RBP to pull down corresponding RNAs!
I did the following pipeline: Hisat2 (mapping the reads after QC) -> htseq-count (getting read counts) -> DESeq2 (getting enriched RNAs in test VS isotype control) -> volcano plot shows up-regulated and down-regulated genes in the test sample.
Logically, I should have seen just enriched RNAs (up-regulated genes), but surprisingly there are some RNAs down-regulated in the test sample (or up-regulated in the isotype control sample)!!!!! How can we explain and interpret this???? it seems some genes are pull down in isotype control sample as well. Are these background? Is my computational pipeline right?
BTW, it's worth to mention that all enriched RNA (called up-regulated genes in DESeq2 pipeline) are significantly relevant and meaningful to my experiment. However, I am not sure how much confident I can be on this output? And how should I explain down-regulated RNAs? because normally we just expect enriched RNAs in RIP-seq experiment.
Any thoughts?