Python教程

python读取bam文件

本文主要是介绍python读取bam文件,对大家解决编程问题具有一定的参考价值,需要的程序猿们随着小编来一起学习吧!

1、读取fasta文件:

(1)方法1:Bio库

from Bio import SeqIO

# 读取包含单个序列 Fasta 格式文件
fa_seq = SeqIO.read("res/sequence1.fasta", "fasta")

seq = str(fa_seq.seq)

# 一个多序列文件中的所有序列
seqs = [fa.seq for fa in SeqIO.parse("res/multi.fasta", "fasta")]

(2)方法2:pysam库

fa = pysam.FastaFile(‘filename’)

# 打印所有的references
fa.references

# 打印所有的references的长度
fa.lengths

# 获取指定references的[a, b)的序列ACTGN
seq = fa.fetch(reference=’chr21’, start=a, end=b) # 范围是前闭后开

# 提取chr21的整条序列
seq = fa.fetch(reference=’chr21’)

2、将序列转化为字符串:str(fa.seq)

# 单个核苷酸计数
print("G Counts: ", fa.count("G"))

# 获取反向序列
print("reverse: ", fa[::-1])

# 获取反向互补序列
print("Reverse complement: ", fa.complement())

3、pysam读取bam文件时,报错:

File "pysam/libcalignmentfile.pyx", line 742, in pysam.libcalignmentfile.AlignmentFile.__cinit__

File "pysam/libcalignmentfile.pyx", line 952, in pysam.libcalignmentfile.AlignmentFile._open

File "pysam/libchtslib.pyx", line 365, in pysam.libchtslib.HTSFile.check_truncation

OSError: no BGZF EOF marker; file may be truncated

可以这样读取文件,即加入ignore_truncation=True:

bam = pysam.AlignmentFile(‘**.bam’, "rb", ignore_truncation=True)

4、CRAM格式的文件:具有高压缩的特性,并且比BAM有更高的压缩率,以后可能多数会使用CRAM格式的文件

5、比对数据(BAM/CRAM/SAM),变异数据:(VCF/BCF)

6、获得每一条read的情况

for read in bam:
    read.reference_name # 比对参考序列染色体名称;

    read.pos # read比对的位置

    read.mapq # read的比对质量值

    read.query_qualities # read sequence base qualities

    read.query_sequence # read sequence bases

    read.reference_length # 在reference genome上read比对的长度

7、读取cram文件

cf = pysam.AlignmentFile(‘*.cram’, ‘rc’)

8、读取sam文件

samfile = pysam.AlignmentFile(‘**.sam’, ‘r’)

9、获取bam文件中某一区域

for r in pysam.AlignmentFile(‘*.bam’, ‘rb’).fetch(‘chr21’, 300, 310):
    pass  # 这样做的前提是:*.bam文件必须有索引

这篇关于python读取bam文件的文章就介绍到这儿,希望我们推荐的文章对大家有所帮助,也希望大家多多支持为之网!