抱歉,您的浏览器无法访问本站
本页面需要浏览器支持(启用)JavaScript
了解详情 >

输入是两个字符串,输出是对齐后的两个字符串。

简介

使用过pysamsamtools的小伙伴肯定了解 pileup的操作,如果把BAM文件看作表格的话,那么通常我们是按行去解析它的record,进而获得一些信息,例如比对到哪条染色体,比对的开始位置和结束位置等. 另一种情况下,我们想要按照列去循环解析,得到这个列上的具体信息,典型的就是这个列上比对序列的碱基是什么?比对序列的位置是什么?以及是Match or Mismatch or indel 等。那么,该操作就需要引入pileup操作了。

首先在shell中测试如下命令

1
2
#!/bin/sh
time gzip -d -c risearch_chr1:143971112-143971134:+:FAM72C.out.gz > risearch_chr1:143971112-143971134:+:FAM72C.out

我们在运行bwa mem比对的时候,由于某些不明的原因会造成程序中断,例如内存超了,IO错误,计算节点崩溃等,然而BAM是否完整很难察觉,最终导致后续流程无法运行。这里,我们通过一段简短的代码来检查BAM文件的完整性,代码如下:

如题,官方已经提供了一个R的版本createGCcontentFile.R ,但是根据代码就能看出这个版本非常占内存了,首先要把基因组整个序列都load入内存中去,每次计算出的矫正数据也是储存dataframe中。为了降低内存占用,也为了提高计算速度,我写了一个julia版本的。代码如下:

Julia读取BAM的库

想要计算Insert size,需要提供一个基因组比对后的文件,sam也好,bam也罢。那么,使用julia语言计算该值的第一步便是了解如何读取和解析BAM文件格式。

julia本身是一门很快速的语言,但是现代计算机往往具有多核心多线程设计,因此,充分发挥硬件,能进一步提高效率

拆分原理

  • 软件的逻辑是首先获取barcode列表。然后采用多线程分别在fastq文件中并行提取对应barcode的reads。

  • WGS的下机数据经常出现在fastq2里。所以程序会从fastq中自动查找是否存在对应barcode。

  • 程序可以自动检测barcode始于开始还是末尾,计算hanming距离,运行1bp的mismatch。

众所周知,计算相关性非常的简单,因为R 语言中有函数cor.test(),该函数可以计算多种方法的相关性检验,返回相关性,Pvalue等检验值,但是这个函数在Julia中并不存在,让Julia作为一门科学计算语言显得并不完美。

前言

自己写的好几种算法企图实现bedtools的功能,虽然julia性能足够好,但都难以在效率上达到bedtools的性能,于是最后只能借助轮子了。