干货|CircRNA 预测软件之find_circ使用流程(二)
在前面的文章中我们讨论了目前研究火热的环状RNA预测软件,通过对同一个样本采用RNAse R+/-两种建库方式,我们比较了5种常用环状RNA预测软件精确度和敏感度。在今天的推文里,我将重点选取find_circ进行详细介绍。
Find_circ工具是最早利用高通量测序数据预测环状RNA的开山鼻祖。该工具是Memczak等人2013年在权威Nature杂志上发表题为“Circular RNAs are a large class of animal RNAs with regulatory potency”的文章时首次发布的,从而掀起了环状RNA的研究热潮。
基于高通量测序数据的环状RNA预测的关键步骤是寻找不能连续比对到基因组或者转录组上的junction read。想要完成这项工作,第一步就是将RNA reads mapping到基因组上,然后去寻找mapping不上的序列。Find_circ将这些mapping不上的reads各取两头20bp(保证可以唯一比对到基因组上),再次mapping到基因组上。接下来,通过短序列比对来判断GU/AG剪切位点,从而推测出潜在的环状RNA序列。图1给出环状RNA的预测过程:
图1 环状RNA预测过程
下面给大家介绍find_circ的工作流程和命令行参数
Find_circ需要运行在装有python 2.7的64位系统上,同时需要安装numpy和pysam这两个python模块。其运行需要借助bowtie2和samtools来完成基因组mapping的过程。
bowtie2 -p 16 --very-sensitive --score-min=C,-15,0 --mm -x /path/to/bowtie2_index -q -1 mate1.fq -2 mate2.fq 2 > map.bowtie2.log | samtools view -hbuS - | samtools sort – output.bam
###bowtie2参数介绍###
-p 使用多线程; --very-sensitive 允许多重比对,报告出最好的一个; --score-min=C,-15,0 设置比对分数函数;--mm 设置I/O模式。
###samtools view参数介绍###
-h 文件包含header line; -b 输出bam格式; -u 输出非压缩的bam格式 –S 忽略版本兼容
samtools view -hf 4 output.bam | samtools view -Sb - > unmapped.bam
python unmapped2anchors.py unmapped.bam | gzip > anchors.qfa.gz
bowtie2 -p 16 --reorder --mm -M20 --score-min=C,-15,0 -q -x /path/to/bowtie2_index -U anchors.qfa.gz | python find_circ.py -G /path/to/chomosomes.fa -p prefix -s find_circ.sites.log > find_circ.sites.bed 2 > find_circ.sites.reads
grep CIRCULAR find_circ.sites.bed | grep -v chrM | gawk '$5>=2' | grep UNAMBIGUOUS_BP | grep ANCHOR_UNIQUE | $path/maxlength.py 100000 > finc_circ.candidates.bed
###利用grep工具筛选出符合以下条件的circRNA: (1)线粒体染色体除外;(2)至少有2个junction read;(3)最大的不超过100kb
Find_circ输出BED文件格式,由18列组成,包含了预测到的环状RNA的各类指标。
关于CricRNA分析的精彩内容,我们将继续分享。