请问这里错在哪里，为什么运行不了！FileNotFoundError: [Errno 2] No such file o...
import sysimport 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
self.POS = info
self.ID = info
self.REF = info
self.ALT = info
self.QUAL = info
self.FILTER = info
self.INFO = [{pair_lst: pair_lst if len(pair_lst)> 1 else ""} for pair_lst in .split(";")]]
self.FORMAT = info.split(":")
self.sample_num = len(info) -7
self.GT = []
for i in range(2):
GT_value = info.split(":")
GT_dict = {}
for g in range(len(GT_value)):
GT_dict] = GT_value
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()
gatk_result = "00007.vcf"
gatk = VCF(gatk_result)
gatk.header
record = gatk.record#这里record存储的是该Record类的地址
record.CHROM
record.line
record.ID #其他的属性同理
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.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.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'}]
for record in gatk:
# compare GATK tumor AF to 0.05
if float(record.GT['AF']) > 0.05:
print(record.line)
snv = "00007.vcf"
result = gatk.header
for record in gatk:
if record.FILTER == "PASS" and float(record.GT['AF']) > 0.05:
result.append(record.line)
# write out result
with open(snv, 'w+') as snvf:
for line in result:
print(line, file = snvf)
record = next(gatk)
print(record.line)
页:
[1]