生信媛公众号编辑、生信必修课之软件安装课程作者
在「Bionano系列」光学图谱混合组装应该怎么做?这篇文章中,我展示了下面这张图。 和之前的图不同的是,我加了几个箭头,这些箭头所指向的区域的特征就是,这些区域并未被Bionano所覆盖。如果不去思考这些区域到底是什么,直接进行混合组装,那么这其实对最后结果的不负责任。
官方并没有一个很详细的文档描述Bionano的从头组装流程的具体过程,所以我只能根据自己实际项目进行介绍: AutoNoise + SplitBNX: 这一步会将bnx和参考的cmap文件进行比对,估算出噪声系数,然后把bnx进行拆分便与后续比对 Pairwse: 这一步进行molecules...
从这部分开始,就开始涉及一些软件的操作和数据分析,因此在进入正文之前,我们需要准备好环境。 环境准备 第一步:从 https://bionanogenomics.com/library/datasets/下载人类测试数据集,以及对应的NA12878人类基因组。
本文的sort命令是GNU版本(8.22), 和BSD的sort不同 sort是我最常用Linux命令之一,它的功能就是排序,一般后面还会和uniq搭配,对数据进行去重。 下面的操作假设你有一个文件,叫做chr.
如果只有一个样本,或者样本量不大的情况下,我会选择一次性投递所有的任务。但是如果有100个以上的样本,那我就得谨慎考虑。 用 snakemake 很好解决这个问题,它会按照你给定的任务数和CPU数,确定每次投递多少任务。
本文只是按照自己的需求翻译了HaploMerger2提供的手册部分内容。HaploMerger2的帮助文档写的非常好,一定要花点时间去读啊! HaploMerger2的分析流程如下 重建单倍体组装中的等位基因关系 检测并纠正二倍体组装中的错连(mis-join) 重建2个单倍型组装 进一步对单倍...
都靠时间的积累
果子老师做了一个非常详细的新手入门R语言的安装策略,叫做新手第1课,无敌无脑的R语言环境配置教程。基本上,你只要照着他的说的做,一字一句的阅读他的文档里的内容(注意,一定要一字一句),基本上R语言就能顺利用起来了。
基因组某些区域可能有着比较高的杂合度,这会导致基因组该区域的两个单倍型被分别组装成primary contig, 而不是一个为primary contig, 另一个是associated haplotig. 如果下游分析主要关注于单倍型,这就会导致一些问题。
软件安装 官方提供了编译好的jar包,方便使用 wget https://github.com/broadinstitute/pilon/releases/download/v1.22/pilon-1.22.jar java -Xmx16G -jar pilon-1.22.jar 如果要顺利运行程序,要求JAVA > 1.7, 以及根据基因组大小而定的内存,一般而言是1M大小的基因对应1GB的内存。
之前写过一篇文章Fastq-dump: 一个神奇的软件, 详细介绍了fastq-dump的用法。 虽然fastq-dump参数很多,而且一直被吐槽参数说明写的太差,但是如果真的要用起来其实也就是一行代码 fastq-dump --gzip --split-3 --defline-qual '+' --defline-seq '@$ac-$si/$ri' SRRXXXXX| SRRXXXX.sra # 加上--gzip后需要时间进行文件压缩 当然除了参数问题,还有一个让人诟病的地方就是他只能单个线程,所以速度特别的慢。
目前尚没有能力直接去阅读htslib的源代码,看到bioawk的代码稍微简单点,因此准备先从这里下手,bioawk的项目地址为https://github.com/lh3/bioawk。
问题描述 序列 指的是一组对象的集合,其中允许重复。序列分为有限序列和无限序列两种类型,我们通常用 表示序列中的第n个对象。 递归其实就是当前的序列依赖于之前的序列。
这里的新版指的是PacBio公司在2018年9月发布pb-assembly, 而这篇文章是在2018年9月30日发的。 今年早些时候在参加三代培训时,听说PacBio会在今年对Falcon进行一些改变。
准备参考序列 注意:这一步必须在后续步骤之前运行。 通常,我们需要准备一个物种的基因组fasta文件,当然RNA和protein都是没有问题。通过prepare-refseqs.pl格式化生成的track,这为后续所有文件提供一个坐标,一直放大后参考序列的碱基也会显示出来。
JBrowse is a fast, scalable genome browser built completely with JavaScript and HTML5. It can run on your desktop, or be embedded in your website. 如果你想要使用JBrowse,一定要有管理员权限,否则建议使用IGV。
但JBrowse要展示数据达到一定量级之后,如何方便管理这些track就成了一个问题,JBrowse支持三种track选择器: JBrowse/View/TrackList/Simple:已经作古的选择方法 JBrowse/View/TrackLis...
DNA核苷酸计数 问题描述: 给定一行核苷酸序列,长度最长为1000 nt, 返回其中'A', 'T', 'C', 'G'出现的次数 C代码如下: #include #include #include #include "dbg.
题图来源: https://pexels.com 只要开始学习,就会出现疑问。即便一个作者认为自己讲的多么具体,但是由于“知识的诅咒”存在,必然有一些点是他认为理所当然,但是读者/听众却从未听过的现象。
Peak Calling 关于MACS2的使用方法, 我写了如何使用MACS进行peak calling详细地介绍了它的参数,在用MACS2之前尽量去阅读下。
样本间相关性评估 上一步得到各个样本的BAM文件之后,就可以在全基因组范围上看看这几个样本之间是否有差异。也就是先将基因组分成N个区间,然后用统计每个区间上比对上的read数。
文件重命名 我们需要对下载的SRRXXXXX文件进行重命名,毕竟有意义的命名才能方便后续展示。那么,应该如何做呢? 首先,你需要将GSE97072页面的中Samples这部分的内容复制到一个文本文件中(我将其命名为sample_name.txt),分为两列,第一列是GSM编号,第二列是样本的命名。
提高自己分析能力的一个好的方法就是重复别人文章里的分析策略,所以这里会尝试对第一篇介绍R-ChIP技术文章"R-ChIP Using Inactive RNase H Reveals Dynamic Coupling of R-loops with T...
作图的目的是希望在图里面发现问题或者解释问题,当然更本质一点就是你想解决什么问题? 前几天做了一个PCA的图,图是画出来了,但是问题有很多,比如说主成分是是啥意思,图里面的箭头有什么含义?为了不做无意义的重复,所以写一篇文章尝试做一个解释。
在《鸟哥的Linux私房菜》第三版的329-330,鸟哥介绍了数据流的重定向,里面提到一点,如何你想把正确与错误数据都写到一个文件中,这个时候用下面第一行代码是错误的,而使用第二行和第三行代码才是正确的 find /home -name .
训练 ab initio 基因预测工具(以SNAP为例) 对于一个新的物种而言,你大概率是没有一个高质量的基因模型去进行基因预测。但是我们可以利用EST序列(少部分物种估计有)、二代测序数据、同源物种蛋白序列,先直接用Maker做基因注释,尽管得到的模型可能不是特别的完美,但可以作为输入反复迭代运行Maker,从而提高最终的表现。
我用的一个服务器上装了一个集群管理工具(SGE, Sun Grid Engine), 用于从登陆节点上向计算节点进行任务投递。一开始,不太会用,但是经过一段时间的摸索学习后,终于能顺手的用起来了。
基因组组装完成后,或者是完成了草图,就不可避免遇到一个问题,需要对基因组序列进行注释。注释之前首先得构建基因模型,有三种策略: 从头注释(de novo prediction):通过已有的概率模型来预测基因结构,在预测剪切位点和UTR区准确性较低 同...
maker 在基因组注释上,MAKER算是一个很强大的分析流程。能够识别重复序列,将EST和蛋白序列比对到基因组,进行从头预测,并在最后整合这三个结果保证结果的可靠性。
最近在服务器上运行matplotlib相关的脚本时遇到了"Invalid DISPLAY variable"报错,从报错中就可以知道这是因为没有显示设备导致的报错。
许久之前,师弟问了我一个问题,为什么shell中添加环境变量的写法是下面这种方式 PATH=~/.lib:$PATH; export PATH 而下面这种会报错呢? $PATH=~/.lib:$PATH; export PATH 当时我的回答是,"shell就是这样子规定的呀"。
最近用简书的markdown写作时,偶尔发现简书又支持几个markdown扩展语法。 第一个是上标, 实现方法是用两个^包围需要上标的字符,如hoptop^TM^, 效果为hoptopTM 第二个是下标, 实现方法是用两个~包围需要下标的字符,如hoptop~TM~, 效果就是hoptopTM 这两个语法在简书官方markdown推文-献给写作者的 Markdown 新手指南-中并没有提到, 算时最近写作的一些小惊喜吧。
每当我写个脚本涉及文件输入时,一般写法都是下面这种 import sys file_in = sys.argv[1] for line in open(file_in, 'r'): commands 这个操作有一个缺点就是,如果我想从管道里面传入数据给Python的话,就会报错,因为原代码要求文件而不是标准输入。
R语言的数据结构原来可以这样理解 这是R数据科学的读书笔记之一,《R数据科学》是一本教你如何用R语言进行数据分析的书。即便我使用R语言快2年多了,但是读这本书还是受益颇多。
Mikado是基于Python3写的基因组结构注释工具,它主要做的是从多个转录组组装工具得到的转录本里挑选出最好的结果作为基因组的结构注释。此外,它还会基于同源蛋白比对结果对转录本打分。
Bioconductor开发的物种注释包系列集合了一个物种不同来源的注释信息,能够根据基因ID对其进行多种来源的注释,比如说基因的别名,基因的功能等。 我之前也写过一篇文章用Bioconductor对基因组注释介绍如何使用AnnotationHub下载注释数据库, 使用select(), mapIds等函数进行注释操作。
R语言中的管道操作 这是R数据科学的读书笔记之一,《R数据科学》是一本教你如何用R语言进行数据分析的书。即便我使用R语言快2年多了,但是读这本书还是受益颇多。
redundans的目标是辅助杂合基因组的组装,输入文件可以是组装的contig,测序文库以及额外的参考基因组,最后用于搭建出scaffold级别的纯合基因组组装结果。
GMAP最早用于讲EST/cDNA序列比对到参考基因组上,可以用于基因组结构注释。后来高通量测序时代,又开发了GSNAP支持高通量数据比对,这篇文章主要介绍GMAP,毕竟高通量转录组数据比对大家更喜欢用STAR, HISTA2等软件。
这篇文章承接R语言中GCC编译的问题,这篇文章主要解决我在Linux系统上安装"expm"出现的问题。 出现的问题 这个问题非常的有趣,因为我在两台服务器分别安装同一个包,其中一台没有任何问题,而另一台则失败,尽管操作系统都是“CentOS Linux release 7.4.1708 (Core)”。
DRAGEN 这是一篇工具介绍贴,考虑这个工具是要钱的,那些动不动就说别人忘了初心的用户肯定认为我写的是软文,所以这些人就不要继续往下看了。
基因组组装完成之后,就需要对最后的质量进行评估。我们希望得到的contig文件中,每个contig都能足够的长,能够有一个完整的基因结构,归纳一下就是3C原则: 连续性(Contiguity): 得到的contig要足够的长 正确性(Correctness): 组装的contig错误率要低 完整性(Completeness):尽可能包含整个原始序列 但是这三条原则其实是相互矛盾的,连续性越高,就意味着要处理更多的模糊节点,会导致整体错误率上升,为了保证完全的正确,那么就会导致contig非常的零碎。
不知不觉间,我给别人上了好多次课程,有些是线上课程,有些是线下课程,有些是收费的,有些是免费的,如果有什么感悟要说,那就是“学习从来都不是一件容易的事。 思考一个问题,你是如何学会数学的乘法和除法的?为什么我要为这个问题呢,这要从这次回家休假,我接到了一个非常高难度的任务说起。
我再次来到了北京,距离上次,时间已过去3年 三年前的那个暑假,机缘巧合,我去参加北京植物所开的光合夏令营。当时来回只能报销火车的硬座票,所以我硬生生从发车的下午坐到了隔天的上午,坐着实在是不舒服,一宿未眠。
我在公众号发过很多编程语言的学习笔记,但是一直没有发Perl语言的编程教程。我大学的时候,学过一段时间的Perl语言,所以和Perl也有点缘分。这次去北京参加培训时发现他们教的Perl,所以接着机会发一波我现场的学习记录。