【直播】我的基因组61:scalpel软件找indel
那么现在正式的开始第61讲:
其实这次的call variation的软件,不仅仅是找到SNV,也顺便找到了indel,只是可能不太准确。一般业界的公认标准是 GATK的best practice,不过那个我已经做了,现在来一点新的,我正好看到了这个scalpel软件。
当然,为什么使用它,完全是随心所欲,也可以选择Pindel等其它软件。我在这里只是为了秀一个软件的用法,生信工程师该如何持续学习。
Scalpel is available here: http://scalpel.sourceforge.net/
文章:
http://www.nature.com/nmeth/journal/v11/n10/full/nmeth.3069.html
软件说明书写的也比较详细:http://scalpel.sourceforge.net/manual.html
他提供了3种情况的找INDELs变异,我目前需要的就是对我的全基因组测序数据来找,所以用single模式。
为了节省对计算资源的消耗,作者建议我单独对每条染色体分别处理。
软件安装是:
## Download and install Scalpel
cd ~/biosoft
mkdir Scalpel && cd Scalpel
wget [url]https://downloads.sourceforge.net/project/scalpel/scalpel-0.5.3.tar.gz[/url]
tar zxvf scalpel-0.5.3.tar.gz
cd scalpel-0.5.3
make
~/biosoft/Scalpel/scalpel-0.5.3/scalpel-discovery --help
~/biosoft/Scalpel/scalpel-0.5.3/scalpel-export --help
它需要自己指定--bed参数来选择染色体运行,而且不是给一个chr1就可以了,需要指定染色体及其起始终止坐标:single region in format chr:start-end (example: 1:31656613-31656883),所以就比较考验shell编程技巧啦!
制作 ~/reference/genome/hg19/hg19.chr.bed 这个文件,我就不多说了,前面我们已经讲过了!
chr10 1 135534747
chr11 1 135006516
chr12 1 133851895
chr13 1 115169878
chr14 1 107349540
chr15 1 102531392
chr16 1 90354753
chr17 1 81195210
chr18 1 78077248
chr19 1 59128983
chr1 1 249250621
chr20 1 63025520
chr21 1 48129895
chr22 1 51304566
chr2 1 243199373
chr3 1 198022430
chr4 1 191154276
chr5 1 180915260
chr6 1 171115067
chr7 1 159138663
chr8 1 146364022
chr9 1 141213431
区分染色体分别运行scalpel软件代码如下:
cat ~/reference/genome/hg19/hg19.chr.bed |while read id
do
arr=($id)
# arr=($a) will split the $a to $arr , ${arr[0]} ${arr[1]} ~~~, but ${arr[@]} is the whole array .
# OLD_IFS="$IFS"
# IFS=","
# arr=($a)
# IFS="$OLD_IFS"
#arr=($a)用于将字符串$a分割到数组$arr ${arr[0]} ${arr[1]} ... 分别存储分割后的数组第1 2 ... 项 ,${arr[@]}存储整个数组。
#变量$IFS存储着分隔符,这里我们将其设为逗号 "," OLD_IFS用于备份默认的分隔符,使用完后将之恢复默认。
echo ${arr[0]}:${arr[1]}-${arr[2]}
date
start=`date +%s`
~/biosoft/Scalpel/scalpel-0.5.3/scalpel-discovery --single
--bam ~/data/project/myGenome/fastq/bamFiles/jmzeng.filter.rmdup.bam
--ref ~/reference/genome/hg19/hg19.fa
--bed ${arr[0]}:${arr[1]}-${arr[2]}
--window 600 --numprocs 5 --dir ${arr[0]}
end=`date +%s`
runtime=$((end-start))
echo "Runtime for ${arr[0]}:${arr[1]}-${arr[2]} was $runtime"
done
最后得到的是每一条染色体一个vcf文件记录着INDEL情况,暂时我还没进行下一步处理。
这里我其实主要是想讲如何用shell进行并行,查看原文可以看到我们的题目及视频讲解,关于这个软件的并行使用!
顺便预告一下,我在wegene测得的芯片数据已经完成了全流程,下载是wegene专题。
还有,我们生信菜鸟团热心群友指出了我前面用常染色体做祖源分析的不足之处,希望我可以继续用Y染色体和线粒体DNA来做下去,给了我几个网址,我估计要学习两个月左右才能完全搞明白,毕竟是孤家寡人兼职学习,有点累,有兴趣的可以学习下面的内容,跟我交流,我的email是jmzeng1314@163.com
https://isogg.org/tree/index.html
http://www.ybrowse.org/gbrowse2/gff/
http://www.phylotree.org/tree/index.htm
https://dna.jameslick.com/mthap/
文:Jimmy
图文编辑:吃瓜群众
- iOS内存管理:从MRC到ARC实践
- MySQL错误修复:Table xx is marked as crashed and last (automatic?) repair failed
- PHP跨站脚本攻击(XSS)漏洞修复方法(一)
- Windows下获取网络连线实际名称,加强IP类设置脚本的兼容性
- Android APP 快速 Pad 化实现
- PHP彩蛋还是漏洞?expose_php彩蛋的触发和屏蔽方法
- 深入源码探索 ReactNative 通信机制
- PHP跨站脚本攻击(XSS)漏洞修复思路(二)
- WordPress发布文章自动同步到新浪微博(带特色图片)
- go http 服务器编程(1)
- Linux系统内存监控、性能诊断工具vmstat命令详解
- go http 服务器编程(2)
- 利用placeholder属性来添加输入框默认文字提示,提高用户体验
- Linux系统监控、诊断工具之top命令详解
- JavaScript 教程
- JavaScript 编辑工具
- JavaScript 与HTML
- JavaScript 与Java
- JavaScript 数据结构
- JavaScript 基本数据类型
- JavaScript 特殊数据类型
- JavaScript 运算符
- JavaScript typeof 运算符
- JavaScript 表达式
- JavaScript 类型转换
- JavaScript 基本语法
- JavaScript 注释
- Javascript 基本处理流程
- Javascript 选择结构
- Javascript if 语句
- Javascript if 语句的嵌套
- Javascript switch 语句
- Javascript 循环结构
- Javascript 循环结构实例
- Javascript 跳转语句
- Javascript 控制语句总结
- Javascript 函数介绍
- Javascript 函数的定义
- Javascript 函数调用
- Javascript 几种特殊的函数
- JavaScript 内置函数简介
- Javascript eval() 函数
- Javascript isFinite() 函数
- Javascript isNaN() 函数
- parseInt() 与 parseFloat()
- escape() 与 unescape()
- Javascript 字符串介绍
- Javascript length属性
- javascript 字符串函数
- Javascript 日期对象简介
- Javascript 日期对象用途
- Date 对象属性和方法
- Javascript 数组是什么
- Javascript 创建数组
- Javascript 数组赋值与取值
- Javascript 数组属性和方法
- mysqldump命令详解 Part 3- 备份全库
- mysqldump命令详解 Part 4-备份单表
- mysqldump命令详解 Part 5-按条件备份表数据
- mysqldump命令详解 Part 6-导出事件,函数和存储过程
- StampedLock知识点梳理
- ReentrantReadWriteLock知识点梳理
- 每天都在用String,你真的了解吗?
- [Go] Golang练习项目-快速排序的GO语言实现
- MySQL学习笔记【基础篇】
- 设计模式~缺省适配模式
- 设计模式~合成模式
- mknod | Linux 后门系列
- Nmap NSE 库分析 >>> base64
- 视频远程通话会议EasyRTC通过SSH部署,关闭SSH后进程停止运行如何解决?
- MapReduce之ReduceJoin案例