标签: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/
标签:style http io 文件 for ar cti line sp
原文地址:http://www.cnblogs.com/lyhonk/p/3931769.html