qualimap对多组学数据比对结果进行QC

2018-07-01 本文已影响19人因地制宜的生信达人

英文文档本来就很清晰，但是需要一些时间来理解： http://qualimap.bioinfo.cipf.es/doc_html/analysis.html

测试数据也给的很全面： http://qualimap.bioinfo.cipf.es/doc_html/samples.html#counts-example-output

对外显子测序用qualimap

需要自己制作 mm10.exon.chr.bed 文件，制作方式如下：

wget ftp://ftp.ncbi.nlm.nih.gov/pub/CCDS//archive/21/CCDS.20161208.txt

cat CCDS.20161208.txt   |perl -alne '{/\[(.*?)\]/;next unless $1;$gene=$F[2];$exons=$1;$exons=~s/\s//g;$exons=~s/-/\t/g;print "$F[0]\t$_\t$gene" foreach split/,/,$exons;}'|sort -u |bedtools sort -i >mm10.exon

awk '{print $0"\t0\t+"}' mm10.exon >mm10.exon.bed

然后就可以运行 qualimap

## $1 for the config file:  bam_path.txt
## $2 and $3 for submit jobs.
exon_bed='/home/jianmingzeng/annotation/CCDS/mouse/mm10.exon.chr.bed'
qualimap='/home/jianmingzeng/biosoft/Qualimap/qualimap_v2.2.1/qualimap'
cat $1 |while read id 
do 
echo $id 
    if((i%$2==$3))
    then
    $qualimap bamqc --java-mem-size=20G -gff $exon_bed -bam $id 
    fi
    i=$((i+1))   
done

可以看外显子的测序情况。

对转录组数据进行QC

这里其实应该是首推RSeQC这个软件，可惜那是个python的，而且运行超慢，还具耗费内存。所以不得已转为

示例报告： http://kokonech.github.io/qualimap/kidney_rnaseqqc/qualimapReport.html

## $1 for the config file:  bam_path.txt
## $2 and $3 for submit jobs.
gtf='/home/jianmingzeng/reference/gtf/gencode/gencode.v25.annotation.gtf'
qualimap='/home/jianmingzeng/biosoft/Qualimap/qualimap_v2.2.1/qualimap'
cat $1 |while read id 
do 
file=$(basename $id )
sample=${file%%.*} 
echo $sample 
    if((i%$2==$3))
    then
    $qualimap rnaseq --java-mem-size=20G -gtf $gtf -bam $id -pe  -oc $sample 
    fi
    i=$((i+1))   
done

count qc

属于转录组数据质控的一部分，比如：6 samples in 2 conditions 的报告，这个时候的input数据是表达矩阵了：

综合比较多个bam文件

示例报告：

qualimap对多组学数据比对结果进行QC

对外显子测序用qualimap

对转录组数据进行QC

count qc

综合比较多个bam文件

猜你喜欢

热点阅读