python通讯录文件读取_python通用读取vcf文件的类(可以直接复制粘贴使用)
前言
处理vcf文件的时候,需要多种切割,正则匹配,如果要自己写其实会比较麻烦,并且每次还得根据vcf文件格式或者需要读取的值不同要修改相应的代码。因此很多人会选择一些python的vcf的库,但是首先你得安装这个库, 并且有一些库它固定了能够读的内容,如果你的vcf的信息不在它固定的里面,就读不出来。比如最近我想读一个样本的AF,但是它放在最后样本的GT那列,不在INFO那一列,有一些库竟然无能为力。
因此我写了这个通用的读vcf的类,直接复制粘贴这部分代码就可以方便的用这个类进行vcf文件的读取,过滤,写出等操作。
使用说明
首先复制类的代码,后面就可以直接用了
import sys
import os
import subprocess
class Record(object):
'''
One line information in vcf file
'''
def __init__(self, line):
info = line.split("\t")
self.line = line
self.CHROM = info[0]
self.POS = info[1]
self.ID = info[2]
self.REF = info[3]
self.ALT = info[4]
self.QUAL = info[5]
self.FILTER = info[6]
self.INFO = [{pair_lst[0]: pair_lst[1] if len(pair_lst)> 1 else ""} for pair_lst in [pair.split("=") for pair in info[7].split(";")]]
self.FORMAT = info[8].split(":")
self.sample_num = len(info) -7
self.GT = []
for i in range(2):
GT_value = info[8 + i +1].split(":")
GT_dict = {}
for g in range(len(GT_value)):
GT_dict[self.FORMAT[g]] = GT_value[g]
self.GT.append(GT_dict)
class VCF(object):
'''
VCF class, read VCF, write VCF, get VCF information
'''
def __init__(self, uncompress_vcf):
self.header = []
self.reader = open(uncompress_vcf, 'r')
self.line = self.reader.readline().strip()
while self.line.startswith('#'):
self.header.append(self.line)
self.line = self.reader.readline().strip()
self.record = Record(self.line)
def __iter__(self):
return self
def __next__(self):
self.line = self.reader.readline().strip()
if self.line != "":
self.record = Record(self.line)
return self.record
else:
self.reader.close()
raise StopIteration()
def reader_close(self):
self.reader.close()
主要有两个类,一个是VCF类,存储的是vcf的信息,及对vcf文件的操作,一个是Record类,它包括vcf某一行存储的全部信息
读入vcf文件
gatk_result = "realignment.vcf"
gatk = VCF(gatk_result)
查看vcf的header
gatk.header
查看vcf当前行中储存的信息,一开始是首行。它以Record这个类保存的。注意VCF类是个迭代器类,可以用next和for循环来读入每一行的信息
record = gatk.record #这里record存储的是该Record类的地址
查看该record的属性,包括line(行的内容,方便写出某行), CHROM, POS, ID,REF,ALT, QUAL, FILTER, INFO(字典的形式存储), FORMAT, sample_num(多少个样本),GT(样本的基因型信息,这里在vcf一般是在后面用样本名表示的列)
record.CHROM
record.line
record.ID #其他的属性同理
INFO的读取
这是vcf中INFO的原始表示
CONTQ=28;DP=38;ECNT=1;GERMQ=76;MBQ=20,37;MFRL=171,229;MMQ=60,60;MPOS=26;NALOD=1.16;NLOD=3.91;
POPAF=6.00;RCNTS=0,0;ROQ=14;SEQQ=1;STRANDQ=11;TLOD=4.56
它在record中的存储形式
record.INFO
[{'CONTQ': '28'}, {'DP': '38'}, {'ECNT': '1'}, {'GERMQ': '76'}, {'MBQ': '20,37'}, {'MFRL': '171,229'}, {'MMQ': '60,60'}, {'MPOS': '26'}, {'NALOD': '1.16'}, {'NLOD': '3.91'}, {'POPAF': '6.00'}, {'RCNTS': '0,0'}, {'ROQ': '14'}, {'SEQQ': '1'}, {'STRANDQ': '11'}, {'TLOD': '4.56'}]
GT的读取
这是GT在vcf的存储形式,FORMAT对应着GT的值
GT:AD:AF:DP:F1R2:F2R1:OBAM:OBAMRC:OBF:OBP:OBQ:OBQRC:SB 0/1:21,2:0.120:23:7,1:13,1:false:true:0.500:0.078:100.00:41.80:12,9,1,1 0/0:13,0:0.065:13:7,0:6,0:false:false:.:.:.:.:10,3,0,0
分别是FORMAT, tumor样本GT, normal样本GT对应的值
这是在record中的存储形式
record.GT
[{'GT': '0/1', 'OBQRC': '41.80', 'SB': '12,9,1,1', 'DP': '23', 'OBF': '0.500', 'OBAM': 'false', 'OBP': '0.078', 'AD': '21,2', 'F2R1': '13,1', 'F1R2': '7,1', 'AF': '0.120', 'OBQ': '100.00', 'OBAMRC': 'true'}, {'GT': '0/0', 'OBQRC': '.', 'SB': '10,3,0,0', 'DP': '13', 'OBF': '.', 'OBAM': 'false', 'OBP': '.', 'AD': '13,0', 'F2R1': '6,0', 'F1R2': '7,0', 'AF': '0.065', 'OBQ': '.', 'OBAMRC': 'false'}]
第一个字典就是tumor的GT,第二个字典就是normal的GT,当然,根据你的样本数量会有多个字典,这里可以按索引取出比如要取出第一个样本的,只需要record.GT[0]就行
把tumor AF大于0.5的line打印出来
for record in gatk:
# compare GATK tumor AF to 0.05
if float(record.GT[0]['AF']) > 0.05:
print(record.line)
把FILTER为PASS的并且tumor AF>0.05写入列表并写出最后的VCF文件
snv = "filter.vcf"
result = gatk.header
for record in gatk:
if record.FILTER == "PASS" and float(record.GT[0]['AF']) > 0.05:
result.append(record.line)
# write out result
with open(snv, 'w+') as snvf:
for line in result:
print(line, file = snvf)
查看gatk的下一个record, 因为VCF类是可迭代的,因此除了for也支持next
record = next(gatk)
print(record.line)
python通讯录文件读取_python通用读取vcf文件的类(可以直接复制粘贴使用)相关推荐
- python open读取_Python:open的文件读取操作,utf-8,UnicodeDecodeError
简要目录: open函数 将文件设置为utf-8编码格式 UnicodeDecodeError f.read() 和 f.read(size) f.readline() 和 f.readlines() ...
- python播放本地视频_python opencv 读取本地视频文件 修改ffmpeg的方法
Python + opencv 读取视频的三种情况: 情况一:通过摄像头采集视频 情况二:通过本地视频文件获取视频 情况三:通过摄像头录制视频,再读取录制的视频 摄像头采集.本地视频文件的读取.写视频 ...
- python处理pdf实例_Python实现读取PDF文件案例
最近有一个本地客户需求是读取PDF文件,然后做自动化处理.这其实是一种典型的RPA自动化需求,简单而言就是模拟人工来操作文件,网页,客户端系统等,只要操作规则定义清楚,就可以实施这种RPA应用,而如果 ...
- python pandas读取文件内容_python pandas读取大型文件
数据处理:pandas处理大型csv文件,使用pandas分块处理大文件 1. 读取限定列 一个 csv文件中,有很多列,而我们只关心其中的某些列是,如果把每行数据都读取出来,在提取信息,显然会增加I ...
- python怎么读数据_Python如何读取数据
一.从文本中读取数据 Python可以读取任何格式的文本数据,使用Python读取文本数据一般分为三个步骤.分别是:定义数据文件.获取文件对象.读取文件内容.下面让小编来分别解释它们是如何实现的. 定 ...
- python读取excel一行一行的读取_python逐行读取excel
python读取excel文件,将每一行都保存为一个列表.比如:文件是a,b,d f,e,g u,w,o 想要变成list1=[a,b,c] list2=[f,e,g] lpython读写excel文 ...
- python导入图片数据_Python中读取图片的6种方式
Python进行图片处理,第一步就是读取图片,这里给大家整理了6种图片的读取方式,并将读取的图片装换成numpy.ndarray()格式.首先需要准备一张照片,假如你有女朋友的话,可以用女朋友的,没有 ...
- python创建空文本文件_Python干货:「文件处理整合」 创建、打开、写入和删除...
什么是Python文件处理? Python文件处理就是如何在python中创建.打开.读取.写入和删除文件的示例.从文件(EXCELSheet.doc文件.文本文件.cv等)或Add中读取数据.是一个 ...
- python 写txt 换行_python中写入txt文件需要换行,以及\r 和\n
在Python中,用open()函数打开一个txt文件,写入一行数据之后需要一个换行 如果直接用 f.write('\n') 只会在后面打印一个字符串'\n',而不是换行' 需要用 f.write(' ...
- python创建文件对象_python基础教程:文件读写
在Linux系统中,一切都是文件.但我们通常说的文件是保存在磁盘上的图片.文档.数据.程序等等.而在程序的IO操作中,很多时候就是从磁盘读写文件.本节我们讲解Python中的文件对象如何操作文件. 创 ...
最新文章
- html景图片怎么设置百分比,8个风光照片拍摄技巧
- java多线程之API初探(一)
- MAC安装mysql8.0.11以及修改root密码
- thinkphp mysql 日志_thinkphp5 日志
- 图像信噪比计算公式_CT 科研设计之图像质量对比研究三步曲
- pytorch自带网络_PyTorch机器学习笔记(1)整好环境
- 前端学习(686):for循环
- python 串口助手 简书_python用pyserial读取串口问题解决
- 第 11 章 日志管理 - 089 - 初探 ELK
- mysql with as 用法_Python之图解with语句
- 男子刷机多次拒绝来电 维修师傅的做法被网友怒赞
- linux socket 104 错误,linux socket连接中 ERRNO错误
- IE如何兼容placeholder属性
- 【python游戏编程之旅】第八篇---pygame游戏开发常用数据结构
- 人人开源代码生成器启报错
- 弱电工程标书制作,从入门到精通
- html全屏banner轮播代码,jquery自动全屏轮播图banner代码
- 打印机生产食物?3D打印技术推出人造肉口感佳
- 下载频道2013上半年超人气精华资源汇总---全都是免积分下载
- 英特尔第十代处理器为什么不支持win7_10代cpu能不能装win7?10代cpu装win7全面分析(支持十代cpu)...