创建python2.7的虚拟环境：
# 创建lumpy-sv使用的python2.7虚拟环境
conda create -n Py27 python=2.7
conda activate Py27
-----------------------------------------------------------------------------------补知识点专用分割线-------------------------------------------------------------------------------------------
-------------------------------------------------------------------conda activate 与 source activate的关系---------------------------------------------------------------------------------
上网查找了一下,发现是版本更新后出现的两条指令,现在使用conda activate激活环境更方便
conda activate：环境激活的逻辑和机制已经被重新设计。在conda4.4中，condaactivate和condadeactivate现在是激活和停用环境的首选命令。您会发现它们比以前conda版本中的source activate和source deactivate命令要快得多。condaactivate命令还有以下优点：（1）在所有操作系统、shell和平台上都是通用的；（2）与来自其他包（如pythonvirtualenv的activate脚本）的脚本没有路径冲突。
顺便提一句,我在使用中出现了如下问题:
报错:CommandNotFoundError: Your shell has not been properly configured to use 'conda deactivate'.
原因是我上一次使用了source activate激活环境后,直接关闭了终端窗口,没有使用source deactivate 退出,因此现在不能使用conda activate 激活我想用的环境,因此只能先执行source deactivate或者conda deactivate(最好两个指令都依次试一下,很玄学),接着就能从新使用source activate命令啦
————————————————
版权声明：本文为CSDN博主「他的北方」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/weixin_44098238/article/details/111946924
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
安装lumpy-sv和samblaster
conda install -c bioconda lumpy-sv
conda install -c bioconda samblaster
参考文献：https://www.cnblogs.com/YlnChen/p/12667431.html
使用方法
LUMPY是一款基于概率框架检测结构变异(structure variants)的软件, 它根据read-pair, split-read, read-depth和其他先验知识寻找基因组上可能的结构变异。
软件在编译的时候会先安装HTSLIB，根据Makefile, 需要预先安装好curl和zlib. 此外还推荐安装Python的Pysam和Numpy，Samtools(0.1.18+)，SAMBLASTER(0.1.19+), sambamba。安装完成之后可以用测试数据进行测试, 数据下载地址为http://layerlab.org/lumpy/data.tar.gz, 要存放在/lumpy-sv/data下
lumpy基于paired-end reads比对后得到的三类信息推断SV，局部异常的测序深度，不一致(discordant)的联配和断裂的联配(split-read alignment)。局部异常的测序深度比较容易理解，平均30X测序的地方，如果深度大于100X，意味着存在着拷贝数变异，如果深度程度非常低，可能意味着这里存在大片段缺失。不一致的联配和断裂的联配能够提供的信息更多，如果基因组一个区域齐刷刷的截断(如下图)，就意味着这个区域可能存在插入/缺失。当然也有其他可能，当两个read在不同链或者不同染色体时，可能是易位或倒置。
read截断
因此，运行lumpy需要预先整理出不一致的短读以及断裂的联配, 如下是lumpy提供的数据预处理方式
# Align the data
bwa mem -R "@RG\tID:id\tSM:sample\tLB:lib" human_g1k_v37.fasta sample.1.fq sample.2.fq \ | samblaster --excludeDups --addMateTags --maxSplitCount 2 --minNonOverlap 20 \ | samtools view -S -b - \ > sample.bam
# Extract the discordant paired-end alignments.
samtools view -b -F 1294 sample.bam > sample.discordants.unsorted.bam
# Extract the split-read alignments
samtools view -h sample.bam \ | scripts/extractSplitReads_BwaMem -i stdin \ | samtools view -Sb - \ > sample.splitters.unsorted.bam (亲测无效，只生成head，里面没有任何东西，具体啥原因没搞明白，下了软件的测试数据，好像也只有一个head，分析测试数据命令为samtools view -h pe.pos_sorted.bam|/data/kongwl/miniconda3/envs/Py27/bin/extractSplitReads_BwaMem -i stdin|samtools view -Sb - | samtools sort - > test1.bam &，此处提取和排序一块做了，难道本来就只有head？？）
# Sort both alignments
samtools sort sample.discordants.unsorted.bam sample.discordantssamtools sort sample.splitters.unsorted.bam sample.splitters
-----------------------------------------------------------------------------------补知识点专用分割线-------------------------------------------------------------------------------------------
-------------------------------------------------------------------samblaster的使用，一步实现去重复，提取split reads和discordants read pairs--------------------------------------------------------------------------------
嫌picard去重复太慢？快来试试这两个高效快速的工具吧（http://glb-biotech.com/zixun1.php?id=350&lm=12）
今天给大家推荐两款sam文件处理小工具samblaster和sambamba，它们具有排序、比对信息查看等常用功能之外，最棒的是可以用来代替picard去除重复序列，在筛选标准不变的前提下速度能提升30倍以上，赶紧来试试吧~
这两款软件比picard快30倍
SAMBAMBA
http://lomereiter.github.io/sambamba/
功能介绍
sambamba主要有filter，merge,slice和duplicate等七个功能来处理sam/bam文件。
一、安装 (支持mac OS/linux 64位)
git clone --recursive https://github.com/lomereiter/sambamba.git
cd sambamba
make
也可以在download里面手动下载，上传到系统解压安装sambamba_v0.6.7_linux.tar.bz2
二、使用方法
1.排序
sambamba sort OPTIONS <input.bam>
主要参数：
-o, --out 设置输出文件的名字（默认 .sorted.bam）
-n, --sort-by-name 按reads id排序（默认按照在参考基因组上的位置排序）
2.建立索引
sambamba index [-p|--show-progress] [-n|--threads=NTHREADS] <input.bam> [<output.bai>]
示例：
$ sambamba index example.bam
#显示处理过程
$ sambamba index --show-progress example.bam /tmp/example.bam.bai
3.提取文件的信息
sambamba view OPTIONS <input.bam | input.sam> [region1 [...]]
主要参数：
-S 输入文件为sam（默认为bam）
-F, --filter=FILTER 过滤提取bam
-f, --format=FORMAT 指定输出文件格式（默认是sam, 还支持bam, json, or msgpack ）
-h, --with-header 保留header
示例：
#显示参考基因组序列基本信息
$ sambamba view --reference-info ex1_header.bam
[{"name":"chr1","length":1575},{"name":"chr2","length":1584}]
#计算3号染色体上质量值大于5且序列长大于80bp的reads个数
$ sambamba view -c -F "ref_id == 3 and mapping_quality >= 50 and sequence_length >= 80" ex1_header.bam
3124
4.合并多个bam文件
sambamba merge OPTIONS <output.bam> <input1.sorted.bam> <input2.sorted.bam> [...] #automatically like in Picard merging tool
主要参数：
-t, --nthreads=NTHREADS 设置线程数
-H, --header 合并sam中的header
-l, --compression-level 按0 to 9设置文件压缩的程度
5.查看reads flag的比对结果
sambamba flagstat OPTIONS <input.bam>
显示以下信息：
First line contains numbers of QC-passed and QC-failed reads. Then come pairs of numbers, the former for QC-passed reads, the latter for QC-failed ones:
· duplicates
· mapped reads (plus percentage relative to the numbers from the first line)
· reads with 'is_paired' flag set
· paired reads which are first mates
· paired reads which are second mates
· paired reads with 'proper_pair' flag set (plus percentage relative to the numbers of QC-passed/failed reads with 'is_paired' flag set)
· paired reads where both mates are mapped
· paired reads where read itself is unmapped but mate is mapped
· paired reads where mate is mapped to a different chromosome
· the same as previous but mapping quality is not less than
6.查重复序列
sambamba markdup OPTIONS <input.bam> <output.bam>
主要参数：
-r, --remove-duplicates 去掉重复序列，否则仅mark重复序列
-t, --nthreads=NTHREADS
-l, --compression-level=N
--tmpdir=TMPDIR 设置中间文件生成目录(默认为/tmp)
此外，还可以提取sam文件的某一段，sambamba slice OPTIONS<input.bam> region
SAMBLASTER
https://github.com/GregoryFaust/samblaster
https://doi.org/10.1093/bioinformatics/btu314
一、安装 (支持linux/mac OS Version 10.7以上)
git clone git://github.com/GregoryFaust/samblaster.git
cd samblaster
make
cp samblaster /usr/local/bin/.
二、使用方法
主要参数：
-i --input 输入sam文件（必须包含header且按reads id排序）
-o --output 输出sam文件
-d --discordantFile 输出discordant read pairs
-s --splitterFile 输出split reads
-u --unmappedFile 输出unmapped/clipped reads
其他参数:
-a --acceptDupMarks 不去重
-e --excludeDups 去掉discordant, splitter, and/or unmapped等重复（具体定义详见samblaster主页）
-r --removeDups 去掉重复(-e --excludeDups类似)
--addMateTags 添加MC and MQ tags
-M 与bwa mem -M 类似
命令组合
bowtie2 -x $index -U $id | samblaster -e -d $sample.disc.sam -s $sample.split.sam | samtools view -Sb - > $sample.clean.bam
示例：
#自动输出discordant read pairs和split read alignments:
bwa mem <idxbase> samp.r1.fq samp.r2.fq | samblaster -e -d samp.disc.sam -s samp.split.sam | samtools view -Sb - > samp.out.bam
#从bam文件中提取 split reads和discordants read pairs
samtools view -h samp.bam | samblaster -a -e -d samp.disc.sam -s samp.split.sam -o /dev/null
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
让人感兴趣的是-F 1294用来提取不一致的联配，用samtools flags 1294可以发现1294表示"PROPER_PAIR,UNMAP,MUNMAP,SECONDARY,DUP"，带上-F意味着以上这些标记在我们筛选的联配记录中都不会出现，也就意味着筛选的记录要符合下面要求
不能是PROPER_PAIR: 就是比对工具认为都正确比对到基因组上，在同一条染色体，在同一条链的情况，常见的就是83,147和99,163
不能是UNMAP和MUNMAP，也就是配对的短读至少有一个能够比对到参考基因组上
也不能是SECONDARY， 也就是他必须是主要联配
光学重复，DUP, 就更加不能要了
于是，经过上一步，那就得到了包含所有数据的sample.bam，不一致的联配sample.discordants.bam 和断裂联配sample.splitters.bam, 使用作者封装好的调用函数进行结构变异检测。
lumpyexpress \ -B sample.bam \ -S sample.splitters.bam \ -D sample.discordants.bam \ -o sample.vcf
在得到的结构变异基础上，作者推荐是用SVTyper进行基因型确定。
提高准确度的方法：根据先验剔除已知低复杂区域和高覆盖的区域，见参考资料的lumpy-sv教程。
----------------------------------------------------引用别人的声明--------------------------------------------------------- --------
---------------------------------------------------仅供自己学习记录使用，如有侵权，联系我删除-------------------------
作者：xuzhougeng
链接：https://www.jianshu.com/p/1844870f3b25
来源：简书
著作权归作者所有。商业转载请联系作者获得授权，非商业转载请注明出处。
网友评论