码迷,mamicode.com
首页 > 其他好文 > 详细

How to use segemehl

时间:2014-08-23 21:37:01      阅读:396      评论:0      收藏:0      [点我收藏+]

标签:style   http   io   文件   for   ar   cti   line   sp   

setp 1. 建立索引

./segemehl.x -x chr.ctidx -y chr.gaidx -d chr.fa -F [1 or 2]

setp 2. 生成sam文件

./segemehl.x -i chr.ctidx -j chr.gaidx -d chr.fa -q myreads.fa -o mymap.sam -F [1 or 2]

setp 3. 利用callMajorMethyl.pl脚本进一步分析

samtools view -hS target.sam | awk ‘/^@/ || /XB:Z:F..CT/‘ > tmp.target.ct.sam

samtools view -bS tmp.target.ct.sam > tmp.target.ct.bam

samtools sort tmp.target.ct.bam tmp.target.ct.sorted

samtools mpileup -Bf chr.fa tmp.target.ct.sorted.bam > tmp.target.test

perl callMajorMethyl.pl -s 1 -i tmp.target.test -o tmp.target.test.txt

========================================================================

在setp3中需循环执行samtools mpileup,需用shell脚本执行

#!/bin/bash
if [ $UID -ne 0 ];
then
echo "Please run as root"
exit 1
fi
##定义.fa文件的路径
fa_path=$1
##定义.sam文件的路径
sam_path=$2
##在sam文件中匹配所需要条目
samtools view -hS ${sam_path} | awk ‘/^@/ || /XB:Z:F..CT/‘ > tmp.ct.sam
##将sam文件转换成bam文件
samtools view -bS tmp.ct.sam > tmp.ct.bam
##将bam文件sort处理
samtools sort tmp.ct.bam tmp.ct.sorted
##对fa文件夹中不同染色体进行比对
for fa in ${fa_path}*.fa;do
name=${fa##*/}
samtools mpileup -Bf ${fa} tmp.ct.sorted.bam > tmp
perl /home/dklv/segemehl/callMajorMethyl.pl -s 1 -i tmp -o {$name}.txt
echo "${name} has been processed "
done
##删除临时文件
rm -rf tmp*
rm -rf ${fa_path}*.fai
echo ""
echo ""
echo ""
echo "=========================================="
echo "All .fa file have been processed"
echo ""

 

==================================================

附录1 segemehl相关手册

http://pan.baidu.com/s/1nt5bJXf

附录2 segemehl项目主页

http://www.bioinf.uni-leipzig.de/Software/segemehl/

How to use segemehl

标签:style   http   io   文件   for   ar   cti   line   sp   

原文地址:http://www.cnblogs.com/lyhonk/p/3931769.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!