本文简单阐述miRNA分析中，如何对其数据进行比对
早先简单介绍如何对miRNA数据进行过滤miRNA分析--数据过滤（一）。就如何比对进行简单介绍。
在比对前，也可以将非miRNA去除，可比对到rRNA，tRNA等库中，将比对上的去掉即可。
比如：
bowtie rRNA.fa -v 0 -p 10 -q mirnaReads -S test.sam --un unmap.sam
# -v: 不错配
# --un: 没有比对上的reads
但是本次我并没有去掉。。
比对有2种方式，第一与miRbase库进行比对，根据gff定量；第二与参考基因组进行比对，还可以鉴定novel miRNA。
1. miRbase 比对
去miRbase 下载mature.fa 或者hairpin.fa文件即可。有些文章比对mature, 也有hairpin.fa 。还有将2者合并进行比较。本次只比较mature。
## 索引
bowtie-build mature.fa
## 比对
bowtie mature.fa -v 0 -p 10 -q miRNareads.fq -S test.sam
##定量
samtools view -F 4 test.sam |cut -f3 |sort |uniq -c | awk '{print $2"\t"$1}' >test.exp.xls
2. 参考基因组进行比较
本次选用mirDeep2
数据准备：
- miRNA reads
- genome
- 本物种的mature.fa （U替换为T）
- 其他物种的mature.fa
- 本物种前体序列
第一步 比对
mapper.pl miRNA_reads -e -j -l 18 -m -p genome -s test.collapsed.fa -t test.arf -v
# -e: 输入为fastq
# -c: 输入为fasta
# -j； 表示移除ATCGUNatcgun以外的字符
# -k: 表示去除3‘街头序列
# -l: 表示最小长度，默认18
# -m: 合併相同的reads
# -p: 表示所索引文件
# -s: 输出的fa
# -t: 输出后的比对文件
# -v: 输出进度报告
同样，mirdeep2也可以接受bam文件
## bam 转化为 sam
samtools view -h -o test.sam test.bam
# mapping
perl bwa_sam_converter.pl -i test.sam -c -o test_collapsed.fa -a test.arf
结果可以得到一个collapsed.fa文件，以及arf文件，其中arf包括序列名；长度；起始；终止；序列；map到哪条染色体；map到染色体的长度；map到染色体的起始终止位置；染色体上的序列；正负义链；错配数；如果有错配用M，否则就用m
第二步： 鉴定novel miRNA并定量
## 如果没有对应本物种的mature.fa 等，则全部用none
miRDeep2.pl test_collapsed.fa genome test.arf none none none 2>>report.log
## 如果有则
miRDeep2.pl test_collapsed.fa genome test.arf 本物种mature.fa other.mature.fa 本物种hairpin.fa -t 本物种 2>>report.log
结果可得到每一个miRNA的 reads count数量，用Deseq2分析即可。
网友评论