*Filter variants per region (in this example, print out only variants mapped to chr1 and chr2)
初步设想在Bioinfo板块中分享一些常见的生信分析软件的使用,原则就是有现成的轮子就不去自己造了。
bcftools 是samtools 的开发者提供的一款专门操作VCF文件的工具,它可以处理VCF格式,也可以处理VCF对应的二进制文件。
bcftools也可以进行SNP calling。在之前的版本中,通常都是和samtools的mpileup命令结合使用, 命令如下
当然了,如何提问,就需要一点点背景知识啦, 比如知道什么是变异位点,什么是过滤,然后就可以很简单的两个提问即可:
VCF 是生物信息分析中非常重要的一种格式。主要用来描述基因组突变的信息,无论是检测出来的 SNP,indel,cnv,还是 SV,都可以存储格式都为 vcf 格式。从比对生成的 bam 文件中,将潜在变异信息筛选出来,就是 vcf 格式。vcf 是一种列表格式,里面包含很多的内容。需要掌握每一列的信息,并能使用相对应的软件对 vcf 进行处理。处理 VCF 格式软件主要包括 bcftools,vcftools,gatk,python pyvcf,plink 等。
https://www.nature.com/articles/s41586-023-06062-z
export BCFTOOLS_PLUGINS=/bi/software/bcftools-1.16/plugins;
随便找一个SNP-calling的综述就可以找到一大堆软件的评价,我这里强烈推荐A survey of tools for variant analysis of next-generation genome sequencing data : http://bib.oxfordjournals.org/content/15/2/256.short 如果你在找变异这一个知识点还很模糊,甚至可以花点时间来翻译一下,同时欢迎大家交流自己的学习心得(http://www.biotrainee.com/threa
在上一次直播中,我们说到了一个不符合我们的认知的问题。就是我的全基因组测序数据里找到的SNV的纯合杂合比例失衡,这着实让我非常纠结。在朋友圈大量求助中,肿瘤所的朋友非常热心的帮我检查了她手头的几百个外
本篇主要介绍annotate, concat, merge, isec, stats这五个命令。
在基因组分析中,处理流程从上游测序数据到下游突变分析,中间的关键就是call突变。
本篇主要介绍index, view, query, sort, reheader这五个命令。
写在前面:『思考问题的熊』专栏上次更新还要追溯到4月19号的 Variant 分析阶段小结1-基础碎碎念,过去接近一个月的时间里我分别经历了两次长途出差和电脑无法连接特定网络的持续尴尬。特定网络正是所有以 qq.com 结尾的网站,当然包括微信公众平台,所以文章都编辑不了。身体和心理的双重袭击让我只能选择围笑:)
PSMC 模型使用单个个体的完整二倍体序列中的信息来推断种群规模变化的历史。它最初于 2011 年发布,现已成为基因组学领域非常流行的工具。在本教程中,我们将逐步完成为 PSMC 生成必要的输入数据的步骤,并在发布的猛犸象数据上运行它。
Sentieon 致力于解决生物信息数据分析中的速度与准确度瓶颈,通过算法的深度优化和企业级的软件工程,大幅度提升NGS数据处理的效率、准确度和可靠性。
组装得到基因组序列之后,进一步的工作就是探究基因结构。Augusust是一款预测真核生物基因结构的软件,官网如下:
但是,我们经常对某些细胞系进行有针对性的设计变异,比如BAF155的R1064K呀,H3F3A的K27呀,那我我们拿到高通量测序数据的时候,就肯定希望可以快速的看看这个基因是否被突变成功了。假设,我们已经得到了所有样本的sort好的bam文件,想看看自己设计的基因突变是否成功了,可以有针对性的只call 某个基因的突变!
我们以这三年疫情的微生物SARS– CoV-2为例,文章:《Genomic Diversity of Severe Acute Respiratory Syndrome– Coronavirus 2 in Patients With Coronavirus Disease 2019》,它就列出来了微生物的DNA测序找变异位点的流程,主要是4个软件,步骤如下所示:
这个snakemake workflow 主要包括:mapping, sort >> index >> call variants
非常多已经造好的轮子可以完成,包括bcftools,vcftools,还有大名鼎鼎的GATK,随便举例如下:
原地址 几点说明 1.非简单翻译,所有代码均可运行,为了辅助理解,基本每步代码都有结果,需要比较的进行了整合 2.原文中的软件都下载最新版本 3.原文中有少量代码是错误的,这里进行了修正 4.对于需要的一些知识背景,在这里进行了注释或链接到他人博客
数据使用GWAS-Cookbook中的GWASdat1中的数据,将数据变为vcf格式。
在PCA(Principal Component Analysis)分析中,常用的工具有EIGENSOFT工具的smartpca,GCTA工具的PCA模块和R包中做PCA分析的princomp函数或glPCA功能。EIGENSOFT工具只支持linux系统,从安装到使用都很复杂。GCTA工具支持不同平台(wins/linux/mac),常用于群体遗传相关分析。在群体遗传中,R包从读取vcf文件、PCA分析到可视化,对内存要求较高。
Snakemake是一款流行的生物信息学工作流管理系统,由Johannes Köster及其团队开发。它旨在降低复杂数据分析的复杂性,使生物信息学工作流的创建和执行变得更加容易和可重复。Snakemake的设计灵感来自于Makefile,但它是专门为生物信息学和数据密集型科学工作流设计的,使用Python语言进行工作流的定义,这使得它在生物信息学社区中特别受欢迎。
找到了一份种群基因组学数据分析的教程,原文用的数据是2015年发表在science上的一篇论文Genomic islands of speciation separate cichlid ecomorphs in an East African crater lake。这份教程利用这篇文章的数据分析了部分内容。
菜鸟团公众号肯定讲过annovar的使用了。比如Nickier的vcf文件的注释及ANNOVAR的使用。
csq命令可以分析SNP位点在基因组上的位置,同时还会预测基因突变对编码蛋白的影响。
如果当前目录下只有vcf格式文件,会遇到报错Failed to open .vcf.gz: could not load index,可以参考 https://www.cnblogs.com/chenwenyan/p/11945445.html
今天的我们,还是继续探究那一个困扰我这么久的问题。为什么我作为堂堂正正的男性,明明X,Y染色体都只有一条,可是却测到了那么多的杂合突变的问题。 在之前,我们在QC阶段详细的探究了X,Y染色体的覆盖度和测序深度,其中X的平均测序深度才16x,而Y却高达60x,我们完全有理由怀疑测序深度对SNV的准确性影响甚大!而且Y染色体总共长度才60M,就有一半是N碱基,有效长度就30M不到,却找到了近3万个SNV,这有着很明显的问题,太密集了~ 所以从测序深度和位点间距来看SNV分布情况是非常有必要的! 对于我的X染色体
背景知识 男性只有一条X染色体和一条Y染色体,所以,理论上它们上面的SNV都应该是纯合的! X,Y除了同源区域外,其它地方差异很大。所以在女性样本里面即使是混入了极低量的男性样本,也很容易检测出来。同理,男性样本里面混入了女性样本,会给男性带来大量的X染色体的杂合SNV,也很容易检测出来。 我的测序结果 我对前面步骤call到的vcf格式的变异位点文件进行了X,Y染色体的简单统计,代码如下: cat jmzeng.freebayes.vcf |grep -w 'chrY'|grep -v "^#" |cu
将测序数据与参考基因组进行比对之后,得到排序并建立索引的 bam 之后,就可以进行 SV检测了。检测 SV 的工具众多,这里面我们推荐使用 sniffles 与 cuteSV 两款软件来进行处理。
基因组结构变异(structure variant, SV)是基因组变异的重要组成部分,大片段插入(Insertion, INS)、缺失(Deletion, DEL)、倒位(Inversion, INV)、易位(Translocation)、重复(Duplication, DUP)等类型的变异。第三代基因组测序因其读长较长,可轻松跨越重复区域和基因组复杂区域,能够更全面的检测基因组的SV。结构变异往往会对基因结构和表达产生更大的影响,在遗传病和肿瘤的发生发展中扮演了重要角色,因此发现和正确注释结构变异对于疾病的诊断有着至关重要的意义。
mpileup是samtools的一个命令,用来生存bcf文件,然后再用bcftools进行SNP和Indel的分析。另外,bcftools是samtools的附带软件。
https://genome.cshlp.org/content/early/2022/09/15/gr.276839.122
未找到原文所用数据,本文使用GATK4.0和全基因组数据分析实践(上)文章中的大肠杆菌基因组作为参考序列,使用wgsim软件模拟生成双端150bp测序数据
构建生信分析流程是生物信息学从业人员必备的技能之一,对该项能力的评估常常是各大公司招录人员的参考项目之一。
前面我们扯到bam文件的各种操作,vcf文件的各种操作,基础知识不牢固的同学可能已经云里雾里了。这次我们来讲一个简单的。就是拿到了fastq的测序数据,如何把全基因组分析给跑一遍。(不谈细节!) 首先就是fastq文件比对到参考基因组变成sam文件: head -40 read1.fq >tmp/read1.fq head -40 read2.fq >tmp/read2.fq ~/biosoft/bwa/bwa-0.7.15/bwa mem -t 20 -M ~/reference/index/bwa/
下载gatk4,用迅雷下载比较快。GATK4下载地址 或者直接wget下载(我的速度慢)
https://www.biorxiv.org/content/10.1101/2023.06.27.545624v1
A1_R1.fa.gz,A1_R2.fa.gz: 二代重测序clean_reads(个体B)【注:必须两个个体A与B】
最初开发 ANNOVAR 时,几乎所有 call 突变的软件都有自己的一套输出格式(SamTools,SOAPSNP,SOLiD BioScope,Illumina CASAVA,CG ASM-var,CG ASM-masterVAR 等),因此 ANNOVAR 就决定采用一种最简单的格式(仅包含 chr, start, end, ref, alt 以及 optional fields)作为输入。现将其称为 avinput 文件。我们也在 ANNOVAR 软件包中提供了 convert2annovar.pl 程序,以方便进行格式转换。
原文出处:https://www.danielecook.com/using-gnu-parallel-for-bioinformatics/
问:有许多朋友后台留言问我,为什么找变异的步骤不用GATK的best practice呢?而是选取bcftools,freebayes这种小众软件。 答:我这里统一回复一下,不同软件找到的variation的差别我前面已经说了,它们小众,并不代表不堪大用,对一个真正的variation来说,不管是什么软件,都是可以找到的,对一个模棱两可的variation,我一定会去去IGV用肉眼查看的,毕竟这可事关我的健康呀,不能马虎的!不同软件就是在sensitivity和specificity之间找平衡,而我早期并不
为了其它相关软件的顺利运行,我们根据教程来设置默认的安装目录及变量环境:Ensembl's VEP , If you don't have VEP installed, then follow this gist.
本文档描述了如何利用Sentieon®基因组学工具的分片能力将DNAseq®流程分布到多台服务器上;将其他流程(如TNseq®)进行分布遵循相同原则,因为所有Sentieon®基因组学工具都具有相同的内置分布式处理能力。这种分布的目标是为了减少流程的总运行时间,以更快地生成结果;然而,这种分布也会带来一些额外的开销,使计算成本增加。 利用分布能力,流程的每个阶段被分成小任务;每个任务处理基因组的一部分,并可以在不同的服务器上并行运行。每个任务生成一个部分结果,需要按顺序合并为最终的单一输出;这种合并需要仔细进行,以确保考虑到边界并生成与没有分片运行的流程相同的结果。 分布的执行框架不在本文档的范围内,用户需要在保持正确的数据依赖关系的同时,分发数据/文件并启动正确的进程。
可重复的生信分析一直是未来的趋势。如果实现可重复的生信分析,关键在于分析软件版本的控制,一致的环境设置还有良好的分析流程的记录。Conda可以说是版本控制和生信工具安装的一大神器。相信大家对它了解肯定不少,但是又该怎么样利用它,进行可重复的分析呢?今天继续讲第二部分 Conda的介绍。
本文描述了使用Sentieon® DNAscope进行PacBio® HiFi数据胚系突变检测。PacBio® HiFi技术产⽣质量值超过Q20的高质量长读段,平均长度在10-25kb之间。准确的长读段可以对短读段和高噪音长读段方法无法检测的基因组重复区域进行精准的变异检测。
领取专属 10元无门槛券
手把手带您无忧上云