test_eurl_vtec_wgs_pt: scripts/ReMatCh/modules/rematch

author	cstrittmatter
date	Wed, 22 Jan 2020 08:41:44 -0500
parents
children	0cbed1c0a762

rev	line source
0 965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	1 import os.path
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	2 import multiprocessing
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	3 import utils
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	4 import functools
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	5 import sys
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	6 import pickle
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	7
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	8
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	9 def index_fasta_samtools(fasta, region_None, region_outfile_none, print_comand_True):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	10 command = ['samtools', 'faidx', fasta, '', '', '']
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	11 shell_true = False
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	12 if region_None is not None:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	13 command[3] = region_None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	14 if region_outfile_none is not None:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	15 command[4] = '>'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	16 command[5] = region_outfile_none
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	17 shell_true = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	18 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, shell_true, None, print_comand_True)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	19 return run_successfully, stdout
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	20
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	21
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	22 # Indexing reference file using Bowtie2
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	23 def indexSequenceBowtie2(referenceFile, threads):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	24 if os.path.isfile(str(referenceFile + '.1.bt2')):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	25 run_successfully = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	26 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	27 command = ['bowtie2-build', '--threads', str(threads), referenceFile, referenceFile]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	28 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, True)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	29 return run_successfully
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	30
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	31
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	32 # Mapping with Bowtie2
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	33 def mappingBowtie2(fastq_files, referenceFile, threads, outdir, conserved_True, numMapLoc, bowtieOPT):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	34 sam_file = os.path.join(outdir, str('alignment.sam'))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	35
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	36 # Index reference file
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	37 run_successfully = indexSequenceBowtie2(referenceFile, threads)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	38
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	39 if run_successfully:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	40 command = ['bowtie2', '-k', str(numMapLoc), '-q', '', '--threads', str(threads), '-x', referenceFile, '', '--no-unal', '', '-S', sam_file]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	41
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	42 if len(fastq_files) == 1:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	43 command[9] = '-U ' + fastq_files[0]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	44 elif len(fastq_files) == 2:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	45 command[9] = '-1 ' + fastq_files[0] + ' -2 ' + fastq_files[1]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	46 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	47 return False, None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	48
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	49 if conserved_True:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	50 command[4] = '--sensitive'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	51 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	52 command[4] = '--very-sensitive-local'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	53
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	54 if bowtieOPT is not None:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	55 command[11] = bowtieOPT
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	56
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	57 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, True)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	58
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	59 if not run_successfully:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	60 sam_file = None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	61
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	62 return run_successfully, sam_file
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	63
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	64
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	65 def split_cigar(cigar):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	66 cigars = ['M', 'I', 'D', 'N', 'S', 'H', 'P', '=', 'X']
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	67
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	68 splited_cigars = []
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	69 numbers = ''
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	70 for char in cigar:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	71 if char not in cigars:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	72 numbers += char
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	73 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	74 splited_cigars.append([int(numbers), char])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	75 numbers = ''
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	76
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	77 return splited_cigars
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	78
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	79
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	80 def recode_cigar_based_on_base_quality(cigar, bases_quality, softClip_baseQuality, mapping_position, direct_strand_true, softClip_cigarFlagRecode):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	81 cigar = split_cigar(cigar)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	82 soft_left = []
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	83 soft_right = []
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	84 cigar_flags_for_reads_length = ('M', 'I', 'S', '=', 'X')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	85 read_length_without_right_s = sum([cigar_part[0] for cigar_part in cigar if cigar_part[1] in cigar_flags_for_reads_length]) - (cigar[len(cigar) - 1][0] if cigar[len(cigar) - 1][1] == 'S' else 0)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	86 for x, base in enumerate(bases_quality):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	87 if ord(base) - 33 >= softClip_baseQuality:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	88 if x <= cigar[0][0] - 1:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	89 if cigar[0][1] == 'S':
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	90 soft_left.append(x)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	91 elif x > read_length_without_right_s - 1:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	92 if cigar[len(cigar) - 1][1] == 'S':
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	93 soft_right.append(x)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	94
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	95 left_changed = (False, 0)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	96 if len(soft_left) > 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	97 soft_left = min(soft_left) + 1
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	98 if soft_left == 1:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	99 cigar = [[cigar[0][0], softClip_cigarFlagRecode]] + cigar[1:]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	100 left_changed = (True, cigar[0][0])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	101 elif cigar[0][0] - soft_left > 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	102 cigar = [[soft_left, 'S']] + [[cigar[0][0] - soft_left, softClip_cigarFlagRecode]] + cigar[1:]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	103 left_changed = (True, cigar[0][0] - soft_left)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	104
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	105 right_changed = (False, 0)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	106 if len(soft_right) > 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	107 soft_right = max(soft_right) + 1
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	108 cigar = cigar[:-1]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	109 if soft_right - read_length_without_right_s > 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	110 cigar.append([soft_right - read_length_without_right_s, softClip_cigarFlagRecode])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	111 right_changed = (True, soft_right - read_length_without_right_s)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	112 if len(bases_quality) - soft_right > 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	113 cigar.append([len(bases_quality) - soft_right, 'S'])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	114
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	115 if left_changed[0]:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	116 # if direct_strand_true:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	117 mapping_position = mapping_position - left_changed[1]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	118 # if right_changed[0]:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	119 # if not direct_strand_true:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	120 # mapping_position = mapping_position + right_changed[1]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	121
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	122 return ''.join([str(cigar_part[0]) + cigar_part[1] for cigar_part in cigar]), str(mapping_position)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	123
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	124
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	125 def verify_is_forward_read(sam_flag_bit):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	126 # 64 = 1000000
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	127 forward_read = False
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	128 bit = format(sam_flag_bit, 'b').zfill(7)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	129 if bit[-7] == '1':
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	130 forward_read = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	131 return forward_read
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	132
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	133
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	134 def verify_mapped_direct_strand(sam_flag_bit):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	135 # 16 = 10000 -> mapped in reverse strand
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	136 direct_strand = False
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	137 bit = format(sam_flag_bit, 'b').zfill(5)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	138 if bit[-5] == '0':
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	139 direct_strand = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	140 return direct_strand
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	141
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	142
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	143 def verify_mapped_tip(reference_length, mapping_position, read_length, cigar):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	144 tip = False
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	145 if 'S' in cigar:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	146 cigar = split_cigar(cigar)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	147 if cigar[0][1] == 'S':
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	148 if mapping_position - cigar[0][0] < 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	149 tip = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	150 if cigar[len(cigar) - 1][1] == 'S':
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	151 if mapping_position + cigar[len(cigar) - 1][0] > reference_length:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	152 tip = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	153 return tip
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	154
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	155
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	156 def change_sam_flag_bit_mapped_reverse_strand_2_direct_strand(sam_flag_bit):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	157 bit = list(format(sam_flag_bit, 'b').zfill(5))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	158 bit[-5] = '0'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	159 return int(''.join(bit), 2)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	160
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	161
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	162 def change_sam_flag_bit_mate_reverse_strand_2_direct_strand(sam_flag_bit):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	163 bit = list(format(sam_flag_bit, 'b').zfill(6))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	164 bit[-6] = '0'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	165 return int(''.join(bit), 2)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	166
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	167
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	168 def move_read_mapped_reverse_strand_2_direct_strand(seq, bases_quality, sam_flag_bit, cigar):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	169 seq = utils.reverse_complement(seq)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	170 bases_quality = ''.join(reversed(list(bases_quality)))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	171 sam_flag_bit = change_sam_flag_bit_mapped_reverse_strand_2_direct_strand(sam_flag_bit)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	172 cigar = ''.join([str(cigar_part[0]) + cigar_part[1] for cigar_part in reversed(split_cigar(cigar))])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	173 return seq, bases_quality, str(sam_flag_bit), cigar
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	174
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	175
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	176 @utils.trace_unhandled_exceptions
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	177 def parallelized_recode_soft_clipping(line_collection, pickleFile, softClip_baseQuality, sequences_length, softClip_cigarFlagRecode):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	178 lines_sam = []
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	179 for line in line_collection:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	180 line = line.splitlines()[0]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	181 if len(line) > 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	182 if line.startswith('@'):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	183 lines_sam.append(line)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	184 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	185 line = line.split('\t')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	186 if not verify_mapped_tip(sequences_length[line[2]], int(line[3]), len(line[9]), line[5]):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	187 line[5], line[3] = recode_cigar_based_on_base_quality(line[5], line[10], softClip_baseQuality, int(line[3]), verify_mapped_direct_strand(int(line[1])), softClip_cigarFlagRecode)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	188 lines_sam.append('\t'.join(line))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	189 with open(pickleFile, 'wb') as writer:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	190 pickle.dump(lines_sam, writer)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	191
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	192
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	193 def recode_soft_clipping_from_sam(sam_file, outdir, threads, softClip_baseQuality, reference_dict, softClip_cigarFlagRecode):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	194 pickle_files = []
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	195 sequences_length = {}
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	196 for x, seq_info in reference_dict.items():
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	197 sequences_length[seq_info['header']] = seq_info['length']
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	198
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	199 with open(sam_file, 'rtU') as reader:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	200 pool = multiprocessing.Pool(processes=threads)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	201 line_collection = []
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	202 x = 0
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	203 for x, line in enumerate(reader):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	204 line_collection.append(line)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	205 if x % 10000 == 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	206 pickleFile = os.path.join(outdir, 'remove_soft_clipping.' + str(x) + '.pkl')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	207 pickle_files.append(pickleFile)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	208 pool.apply_async(parallelized_recode_soft_clipping, args=(line_collection, pickleFile, softClip_baseQuality, sequences_length, softClip_cigarFlagRecode,))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	209 line_collection = []
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	210 if len(line_collection) > 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	211 pickleFile = os.path.join(outdir, 'remove_soft_clipping.' + str(x) + '.pkl')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	212 pickle_files.append(pickleFile)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	213 pool.apply_async(parallelized_recode_soft_clipping, args=(line_collection, pickleFile, softClip_baseQuality, sequences_length, softClip_cigarFlagRecode,))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	214 line_collection = []
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	215 pool.close()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	216 pool.join()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	217
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	218 os.remove(sam_file)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	219
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	220 new_sam_file = os.path.join(outdir, 'alignment_with_soft_clipping_recoded.sam')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	221 with open(new_sam_file, 'wt') as writer:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	222 for pickleFile in pickle_files:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	223 if os.path.isfile(pickleFile):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	224 lines_sam = None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	225 with open(pickleFile, 'rb') as reader:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	226 lines_sam = pickle.load(reader)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	227 if lines_sam is not None:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	228 for line in lines_sam:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	229 writer.write(line + '\n')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	230 os.remove(pickleFile)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	231
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	232 return new_sam_file
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	233
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	234
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	235 # Sort alignment file
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	236 def sortAlignment(alignment_file, output_file, sortByName_True, threads):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	237 outFormat_string = os.path.splitext(output_file)[1][1:].lower()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	238 command = ['samtools', 'sort', '-o', output_file, '-O', outFormat_string, '', '-@', str(threads), alignment_file]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	239 if sortByName_True:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	240 command[6] = '-n'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	241 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, True)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	242 if not run_successfully:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	243 output_file = None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	244 return run_successfully, output_file
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	245
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	246
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	247 # Index alignment file
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	248 def indexAlignment(alignment_file):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	249 command = ['samtools', 'index', alignment_file]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	250 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, True)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	251 return run_successfully
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	252
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	253
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	254 def mapping_reads(fastq_files, reference_file, threads, outdir, conserved_True, numMapLoc, rematch_run, softClip_baseQuality, softClip_recodeRun, reference_dict, softClip_cigarFlagRecode, bowtieOPT):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	255 # Create a symbolic link to the reference_file
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	256 reference_link = os.path.join(outdir, os.path.basename(reference_file))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	257 os.symlink(reference_file, reference_link)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	258
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	259 bam_file = None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	260 # Mapping reads using Bowtie2
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	261 run_successfully, sam_file = mappingBowtie2(fastq_files, reference_link, threads, outdir, conserved_True, numMapLoc, bowtieOPT)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	262
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	263 if run_successfully:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	264 # Remove soft clipping
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	265 if rematch_run == softClip_recodeRun or softClip_recodeRun == 'both':
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	266 print 'Recoding soft clipped regions'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	267 sam_file = recode_soft_clipping_from_sam(sam_file, outdir, threads, softClip_baseQuality, reference_dict, softClip_cigarFlagRecode)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	268
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	269 # Convert sam to bam and sort bam
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	270 run_successfully, bam_file = sortAlignment(sam_file, str(os.path.splitext(sam_file)[0] + '.bam'), False, threads)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	271
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	272 if run_successfully:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	273 os.remove(sam_file)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	274 # Index bam
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	275 run_successfully = indexAlignment(bam_file)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	276
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	277 return run_successfully, bam_file, reference_link
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	278
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	279
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	280 def create_vcf(bam_file, sequence_to_analyse, outdir, counter, reference_file):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	281 gene_vcf = os.path.join(outdir, 'samtools_mpileup.sequence_' + str(counter) + '.vcf')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	282
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	283 command = ['samtools', 'mpileup', '--count-orphans', '--no-BAQ', '--min-BQ', '0', '--min-MQ', str(7), '--fasta-ref', reference_file, '--region', sequence_to_analyse, '--output', gene_vcf, '--VCF', '--uncompressed', '--output-tags', 'INFO/AD,AD,DP', bam_file]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	284
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	285 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, False)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	286 if not run_successfully:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	287 gene_vcf = None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	288 return run_successfully, gene_vcf
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	289
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	290
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	291 # Read vcf file
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	292 class Vcf():
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	293 def __init__(self, vcfFile):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	294 self.vcf = open(vcfFile, 'rtU')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	295 self.line_read = self.vcf.readline()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	296 while self.line_read.startswith('#'):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	297 self.line_read = self.vcf.readline()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	298 self.line = self.line_read
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	299
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	300 def readline(self):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	301 self.line_stored = self.line
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	302 self.line = self.vcf.readline()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	303 return self.line_stored
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	304
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	305 def close(self):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	306 self.vcf.close()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	307
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	308
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	309 def get_variants(gene_vcf):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	310 variants = {}
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	311
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	312 vfc_file = Vcf(gene_vcf)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	313 line = vfc_file.readline()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	314 while len(line) > 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	315 fields = line.splitlines()[0].split('\t')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	316 if len(fields) > 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	317 fields[1] = int(fields[1])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	318
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	319 info_field = {}
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	320 for i in fields[7].split(';'):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	321 i = i.split('=')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	322 if len(i) > 1:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	323 info_field[i[0]] = i[1]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	324 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	325 info_field[i[0]] = None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	326
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	327 format_field = {}
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	328 format_field_name = fields[8].split(':')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	329 format_data = fields[9].split(':')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	330
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	331 for i in range(0, len(format_data)):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	332 format_field[format_field_name[i]] = format_data[i].split(',')
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	333
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	334 fields_to_store = {'REF': fields[3], 'ALT': fields[4].split(','), 'info': info_field, 'format': format_field}
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	335 if fields[1] in variants:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	336 variants[fields[1]][len(variants[fields[1]])] = fields_to_store
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	337 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	338 variants[fields[1]] = {0: fields_to_store}
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	339
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	340 line = vfc_file.readline()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	341 vfc_file.close()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	342
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	343 return variants
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	344
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	345
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	346 def indel_entry(variant_position):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	347 entry_with_indel = []
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	348 entry_with_snp = None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	349 for i in variant_position:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	350 keys = variant_position[i]['info'].keys()
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	351 if 'INDEL' in keys:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	352 entry_with_indel.append(i)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	353 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	354 entry_with_snp = i
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	355
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	356 return entry_with_indel, entry_with_snp
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	357
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	358
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	359 def get_alt_noMatter(variant_position, indel_true):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	360 dp = sum(map(int, variant_position['format']['AD']))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	361 index_alleles_sorted_position = sorted(zip(map(int, variant_position['format']['AD']), range(0, len(variant_position['format']['AD']))), reverse=True)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	362 index_dominant_allele = None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	363 if not indel_true:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	364 ad_idv = index_alleles_sorted_position[0][0]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	365
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	366 if len([x for x in index_alleles_sorted_position if x[0] == ad_idv]) > 1:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	367 index_alleles_sorted_position = sorted([x for x in index_alleles_sorted_position if x[0] == ad_idv])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	368
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	369 index_dominant_allele = index_alleles_sorted_position[0][1]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	370 if index_dominant_allele == 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	371 alt = '.'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	372 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	373 alt = variant_position['ALT'][index_dominant_allele - 1]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	374
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	375 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	376 ad_idv = int(variant_position['info']['IDV'])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	377
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	378 if float(ad_idv) / float(dp) >= 0.5:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	379 if len([x for x in index_alleles_sorted_position if x[0] == index_alleles_sorted_position[0][0]]) > 1:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	380 index_alleles_sorted_position = sorted([x for x in index_alleles_sorted_position if x[0] == index_alleles_sorted_position[0][0]])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	381
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	382 index_dominant_allele = index_alleles_sorted_position[0][1]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	383 if index_dominant_allele == 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	384 alt = '.'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	385 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	386 alt = variant_position['ALT'][index_dominant_allele - 1]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	387 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	388 ad_idv = int(variant_position['format']['AD'][0])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	389 alt = '.'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	390
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	391 return alt, dp, ad_idv, index_dominant_allele
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	392
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	393
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	394 def count_number_diferences(ref, alt):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	395 number_diferences = 0
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	396
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	397 if len(ref) != len(alt):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	398 number_diferences += 1
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	399
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	400 for i in range(0, min(len(ref), len(alt))):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	401 if alt[i] != 'N' and ref[i] != alt[i]:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	402 number_diferences += 1
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	403
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	404 return number_diferences
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	405
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	406
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	407 def get_alt_correct(variant_position, alt_noMatter, dp, ad_idv, index_dominant_allele, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	408 alt = None
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	409 low_coverage = False
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	410 multiple_alleles = False
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	411
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	412 if dp >= minimum_depth_presence:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	413 if dp < minimum_depth_call:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	414 alt = 'N' * len(variant_position['REF'])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	415 low_coverage = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	416 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	417 if ad_idv < minimum_depth_call:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	418 alt = 'N' * len(variant_position['REF'])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	419 low_coverage = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	420 if float(ad_idv) / float(dp) < minimum_depth_frequency_dominant_allele:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	421 multiple_alleles = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	422 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	423 if float(ad_idv) / float(dp) < minimum_depth_frequency_dominant_allele:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	424 alt = 'N' * len(variant_position['REF'])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	425 if index_dominant_allele is not None:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	426 variants_coverage = [int(variant_position['format']['AD'][i]) for i in range(0, len(variant_position['ALT']) + 1) if i != index_dominant_allele]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	427 if sum(variants_coverage) > 0:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	428 if float(max(variants_coverage)) / float(sum(variants_coverage)) > 0.5:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	429 multiple_alleles = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	430 elif float(max(variants_coverage)) / float(sum(variants_coverage)) == 0.5 and len(variants_coverage) > 2:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	431 multiple_alleles = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	432 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	433 multiple_alleles = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	434 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	435 alt = alt_noMatter
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	436 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	437 low_coverage = True
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	438
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	439 return alt, low_coverage, multiple_alleles
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	440
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	441
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	442 def get_alt_alignment(ref, alt):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	443 if alt is None:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	444 alt = 'N' * len(ref)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	445 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	446 if len(ref) != len(alt):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	447 if len(alt) < len(ref):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	448 if alt == '.':
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	449 alt = ref
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	450 alt += 'N' * (len(ref) - len(alt))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	451 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	452 if alt[:len(ref)] == ref:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	453 alt = '.'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	454 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	455 alt = alt[:len(ref)]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	456
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	457 return alt
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	458
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	459
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	460 def get_indel_more_likely(variant_position, indels_entry):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	461 indel_coverage = {}
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	462 for i in indels_entry:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	463 indel_coverage[i] = int(variant_position['info']['IDV'])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	464 return indel_coverage.index(str(max(indel_coverage.values())))
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	465
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	466
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	467 def determine_variant(variant_position, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, indel_true):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	468 alt_noMatter, dp, ad_idv, index_dominant_allele = get_alt_noMatter(variant_position, indel_true)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	469
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	470 alt_correct, low_coverage, multiple_alleles = get_alt_correct(variant_position, alt_noMatter, dp, ad_idv, index_dominant_allele, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	471
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	472 alt_alignment = get_alt_alignment(variant_position['REF'], alt_correct)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	473
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	474 return variant_position['REF'], alt_correct, low_coverage, multiple_alleles, alt_noMatter, alt_alignment
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	475
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	476
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	477 def confirm_nucleotides_indel(ref, alt, variants, position_start_indel, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, alignment_true):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	478 alt = list(alt)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	479
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	480 for i in range(0, len(alt) - 1):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	481 if len(alt) < len(ref):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	482 new_position = position_start_indel + len(alt) - i - 1
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	483 alt_position = len(alt) - i - 1
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	484 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	485 if i + 1 > len(ref):
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	486 break
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	487 new_position = position_start_indel + 1 + i
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	488 alt_position = 1 + i
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	489
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	490 if alt[alt_position] != 'N':
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	491 if new_position not in variants:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	492 if alignment_true:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	493 alt[alt_position] = 'N'
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	494 else:
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	495 alt = alt[: alt_position]
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	496 break
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	497
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	498 entry_with_indel, entry_with_snp = indel_entry(variants[new_position])
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	499 new_ref, alt_correct, low_coverage, multiple_alleles, alt_noMatter, alt_alignment = determine_variant(variants[new_position][entry_with_snp], minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, False)
965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8 cstrittmatter parents: diff changeset	500 if alt_noMatter != '.' and alt[alt_position] != alt_noMatter:

0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1 import os.path

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

2 import multiprocessing

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

3 import utils

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

4 import functools

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

5 import sys

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

6 import pickle

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

7

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

8

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

9 def index_fasta_samtools(fasta, region_None, region_outfile_none, print_comand_True):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

10 command = ['samtools', 'faidx', fasta, '', '', '']

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

11 shell_true = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

12 if region_None is not None:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

13 command[3] = region_None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

14 if region_outfile_none is not None:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

15 command[4] = '>'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

16 command[5] = region_outfile_none

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

17 shell_true = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

18 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, shell_true, None, print_comand_True)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

19 return run_successfully, stdout

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

20

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

21

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

22 # Indexing reference file using Bowtie2

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

23 def indexSequenceBowtie2(referenceFile, threads):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

24 if os.path.isfile(str(referenceFile + '.1.bt2')):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

25 run_successfully = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

26 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

27 command = ['bowtie2-build', '--threads', str(threads), referenceFile, referenceFile]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

28 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, True)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

29 return run_successfully

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

30

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

31

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

32 # Mapping with Bowtie2

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

33 def mappingBowtie2(fastq_files, referenceFile, threads, outdir, conserved_True, numMapLoc, bowtieOPT):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

34 sam_file = os.path.join(outdir, str('alignment.sam'))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

35

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

36 # Index reference file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

37 run_successfully = indexSequenceBowtie2(referenceFile, threads)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

38

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

39 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

40 command = ['bowtie2', '-k', str(numMapLoc), '-q', '', '--threads', str(threads), '-x', referenceFile, '', '--no-unal', '', '-S', sam_file]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

41

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

42 if len(fastq_files) == 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

43 command[9] = '-U ' + fastq_files[0]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

44 elif len(fastq_files) == 2:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

45 command[9] = '-1 ' + fastq_files[0] + ' -2 ' + fastq_files[1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

46 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

47 return False, None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

48

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

49 if conserved_True:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

50 command[4] = '--sensitive'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

51 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

52 command[4] = '--very-sensitive-local'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

53

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

54 if bowtieOPT is not None:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

55 command[11] = bowtieOPT

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

56

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

57 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, True)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

58

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

59 if not run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

60 sam_file = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

61

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

62 return run_successfully, sam_file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

63

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

64

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

65 def split_cigar(cigar):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

66 cigars = ['M', 'I', 'D', 'N', 'S', 'H', 'P', '=', 'X']

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

67

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

68 splited_cigars = []

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

69 numbers = ''

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

70 for char in cigar:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

71 if char not in cigars:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

72 numbers += char

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

73 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

74 splited_cigars.append([int(numbers), char])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

75 numbers = ''

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

76

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

77 return splited_cigars

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

78

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

79

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

80 def recode_cigar_based_on_base_quality(cigar, bases_quality, softClip_baseQuality, mapping_position, direct_strand_true, softClip_cigarFlagRecode):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

81 cigar = split_cigar(cigar)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

82 soft_left = []

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

83 soft_right = []

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

84 cigar_flags_for_reads_length = ('M', 'I', 'S', '=', 'X')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

85 read_length_without_right_s = sum([cigar_part[0] for cigar_part in cigar if cigar_part[1] in cigar_flags_for_reads_length]) - (cigar[len(cigar) - 1][0] if cigar[len(cigar) - 1][1] == 'S' else 0)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

86 for x, base in enumerate(bases_quality):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

87 if ord(base) - 33 >= softClip_baseQuality:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

88 if x <= cigar[0][0] - 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

89 if cigar[0][1] == 'S':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

90 soft_left.append(x)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

91 elif x > read_length_without_right_s - 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

92 if cigar[len(cigar) - 1][1] == 'S':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

93 soft_right.append(x)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

94

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

95 left_changed = (False, 0)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

96 if len(soft_left) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

97 soft_left = min(soft_left) + 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

98 if soft_left == 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

99 cigar = [[cigar[0][0], softClip_cigarFlagRecode]] + cigar[1:]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

100 left_changed = (True, cigar[0][0])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

101 elif cigar[0][0] - soft_left > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

102 cigar = [[soft_left, 'S']] + [[cigar[0][0] - soft_left, softClip_cigarFlagRecode]] + cigar[1:]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

103 left_changed = (True, cigar[0][0] - soft_left)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

104

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

105 right_changed = (False, 0)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

106 if len(soft_right) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

107 soft_right = max(soft_right) + 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

108 cigar = cigar[:-1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

109 if soft_right - read_length_without_right_s > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

110 cigar.append([soft_right - read_length_without_right_s, softClip_cigarFlagRecode])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

111 right_changed = (True, soft_right - read_length_without_right_s)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

112 if len(bases_quality) - soft_right > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

113 cigar.append([len(bases_quality) - soft_right, 'S'])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

114

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

115 if left_changed[0]:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

116 # if direct_strand_true:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

117 mapping_position = mapping_position - left_changed[1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

118 # if right_changed[0]:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

119 # if not direct_strand_true:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

120 # mapping_position = mapping_position + right_changed[1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

121

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

122 return ''.join([str(cigar_part[0]) + cigar_part[1] for cigar_part in cigar]), str(mapping_position)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

123

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

124

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

125 def verify_is_forward_read(sam_flag_bit):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

126 # 64 = 1000000

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

127 forward_read = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

128 bit = format(sam_flag_bit, 'b').zfill(7)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

129 if bit[-7] == '1':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

130 forward_read = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

131 return forward_read

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

132

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

133

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

134 def verify_mapped_direct_strand(sam_flag_bit):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

135 # 16 = 10000 -> mapped in reverse strand

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

136 direct_strand = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

137 bit = format(sam_flag_bit, 'b').zfill(5)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

138 if bit[-5] == '0':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

139 direct_strand = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

140 return direct_strand

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

141

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

142

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

143 def verify_mapped_tip(reference_length, mapping_position, read_length, cigar):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

144 tip = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

145 if 'S' in cigar:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

146 cigar = split_cigar(cigar)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

147 if cigar[0][1] == 'S':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

148 if mapping_position - cigar[0][0] < 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

149 tip = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

150 if cigar[len(cigar) - 1][1] == 'S':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

151 if mapping_position + cigar[len(cigar) - 1][0] > reference_length:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

152 tip = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

153 return tip

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

154

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

155

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

156 def change_sam_flag_bit_mapped_reverse_strand_2_direct_strand(sam_flag_bit):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

157 bit = list(format(sam_flag_bit, 'b').zfill(5))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

158 bit[-5] = '0'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

159 return int(''.join(bit), 2)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

160

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

161

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

162 def change_sam_flag_bit_mate_reverse_strand_2_direct_strand(sam_flag_bit):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

163 bit = list(format(sam_flag_bit, 'b').zfill(6))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

164 bit[-6] = '0'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

165 return int(''.join(bit), 2)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

166

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

167

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

168 def move_read_mapped_reverse_strand_2_direct_strand(seq, bases_quality, sam_flag_bit, cigar):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

169 seq = utils.reverse_complement(seq)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

170 bases_quality = ''.join(reversed(list(bases_quality)))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

171 sam_flag_bit = change_sam_flag_bit_mapped_reverse_strand_2_direct_strand(sam_flag_bit)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

172 cigar = ''.join([str(cigar_part[0]) + cigar_part[1] for cigar_part in reversed(split_cigar(cigar))])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

173 return seq, bases_quality, str(sam_flag_bit), cigar

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

174

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

175

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

176 @utils.trace_unhandled_exceptions

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

177 def parallelized_recode_soft_clipping(line_collection, pickleFile, softClip_baseQuality, sequences_length, softClip_cigarFlagRecode):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

178 lines_sam = []

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

179 for line in line_collection:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

180 line = line.splitlines()[0]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

181 if len(line) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

182 if line.startswith('@'):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

183 lines_sam.append(line)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

184 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

185 line = line.split('\t')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

186 if not verify_mapped_tip(sequences_length[line[2]], int(line[3]), len(line[9]), line[5]):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

187 line[5], line[3] = recode_cigar_based_on_base_quality(line[5], line[10], softClip_baseQuality, int(line[3]), verify_mapped_direct_strand(int(line[1])), softClip_cigarFlagRecode)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

188 lines_sam.append('\t'.join(line))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

189 with open(pickleFile, 'wb') as writer:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

190 pickle.dump(lines_sam, writer)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

191

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

192

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

193 def recode_soft_clipping_from_sam(sam_file, outdir, threads, softClip_baseQuality, reference_dict, softClip_cigarFlagRecode):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

194 pickle_files = []

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

195 sequences_length = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

196 for x, seq_info in reference_dict.items():

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

197 sequences_length[seq_info['header']] = seq_info['length']

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

198

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

199 with open(sam_file, 'rtU') as reader:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

200 pool = multiprocessing.Pool(processes=threads)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

201 line_collection = []

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

202 x = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

203 for x, line in enumerate(reader):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

204 line_collection.append(line)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

205 if x % 10000 == 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

206 pickleFile = os.path.join(outdir, 'remove_soft_clipping.' + str(x) + '.pkl')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

207 pickle_files.append(pickleFile)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

208 pool.apply_async(parallelized_recode_soft_clipping, args=(line_collection, pickleFile, softClip_baseQuality, sequences_length, softClip_cigarFlagRecode,))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

209 line_collection = []

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

210 if len(line_collection) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

211 pickleFile = os.path.join(outdir, 'remove_soft_clipping.' + str(x) + '.pkl')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

212 pickle_files.append(pickleFile)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

213 pool.apply_async(parallelized_recode_soft_clipping, args=(line_collection, pickleFile, softClip_baseQuality, sequences_length, softClip_cigarFlagRecode,))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

214 line_collection = []

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

215 pool.close()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

216 pool.join()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

217

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

218 os.remove(sam_file)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

219

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

220 new_sam_file = os.path.join(outdir, 'alignment_with_soft_clipping_recoded.sam')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

221 with open(new_sam_file, 'wt') as writer:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

222 for pickleFile in pickle_files:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

223 if os.path.isfile(pickleFile):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

224 lines_sam = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

225 with open(pickleFile, 'rb') as reader:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

226 lines_sam = pickle.load(reader)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

227 if lines_sam is not None:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

228 for line in lines_sam:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

229 writer.write(line + '\n')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

230 os.remove(pickleFile)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

231

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

232 return new_sam_file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

233

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

234

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

235 # Sort alignment file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

236 def sortAlignment(alignment_file, output_file, sortByName_True, threads):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

237 outFormat_string = os.path.splitext(output_file)[1][1:].lower()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

238 command = ['samtools', 'sort', '-o', output_file, '-O', outFormat_string, '', '-@', str(threads), alignment_file]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

239 if sortByName_True:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

240 command[6] = '-n'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

241 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, True)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

242 if not run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

243 output_file = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

244 return run_successfully, output_file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

245

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

246

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

247 # Index alignment file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

248 def indexAlignment(alignment_file):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

249 command = ['samtools', 'index', alignment_file]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

250 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, True)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

251 return run_successfully

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

252

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

253

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

254 def mapping_reads(fastq_files, reference_file, threads, outdir, conserved_True, numMapLoc, rematch_run, softClip_baseQuality, softClip_recodeRun, reference_dict, softClip_cigarFlagRecode, bowtieOPT):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

255 # Create a symbolic link to the reference_file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

256 reference_link = os.path.join(outdir, os.path.basename(reference_file))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

257 os.symlink(reference_file, reference_link)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

258

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

259 bam_file = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

260 # Mapping reads using Bowtie2

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

261 run_successfully, sam_file = mappingBowtie2(fastq_files, reference_link, threads, outdir, conserved_True, numMapLoc, bowtieOPT)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

262

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

263 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

264 # Remove soft clipping

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

265 if rematch_run == softClip_recodeRun or softClip_recodeRun == 'both':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

266 print 'Recoding soft clipped regions'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

267 sam_file = recode_soft_clipping_from_sam(sam_file, outdir, threads, softClip_baseQuality, reference_dict, softClip_cigarFlagRecode)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

268

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

269 # Convert sam to bam and sort bam

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

270 run_successfully, bam_file = sortAlignment(sam_file, str(os.path.splitext(sam_file)[0] + '.bam'), False, threads)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

271

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

272 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

273 os.remove(sam_file)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

274 # Index bam

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

275 run_successfully = indexAlignment(bam_file)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

276

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

277 return run_successfully, bam_file, reference_link

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

278

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

279

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

280 def create_vcf(bam_file, sequence_to_analyse, outdir, counter, reference_file):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

281 gene_vcf = os.path.join(outdir, 'samtools_mpileup.sequence_' + str(counter) + '.vcf')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

282

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

283 command = ['samtools', 'mpileup', '--count-orphans', '--no-BAQ', '--min-BQ', '0', '--min-MQ', str(7), '--fasta-ref', reference_file, '--region', sequence_to_analyse, '--output', gene_vcf, '--VCF', '--uncompressed', '--output-tags', 'INFO/AD,AD,DP', bam_file]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

284

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

285 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

286 if not run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

287 gene_vcf = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

288 return run_successfully, gene_vcf

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

289

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

290

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

291 # Read vcf file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

292 class Vcf():

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

293 def __init__(self, vcfFile):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

294 self.vcf = open(vcfFile, 'rtU')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

295 self.line_read = self.vcf.readline()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

296 while self.line_read.startswith('#'):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

297 self.line_read = self.vcf.readline()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

298 self.line = self.line_read

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

299

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

300 def readline(self):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

301 self.line_stored = self.line

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

302 self.line = self.vcf.readline()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

303 return self.line_stored

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

304

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

305 def close(self):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

306 self.vcf.close()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

307

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

308

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

309 def get_variants(gene_vcf):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

310 variants = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

311

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

312 vfc_file = Vcf(gene_vcf)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

313 line = vfc_file.readline()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

314 while len(line) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

315 fields = line.splitlines()[0].split('\t')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

316 if len(fields) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

317 fields[1] = int(fields[1])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

318

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

319 info_field = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

320 for i in fields[7].split(';'):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

321 i = i.split('=')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

322 if len(i) > 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

323 info_field[i[0]] = i[1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

324 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

325 info_field[i[0]] = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

326

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

327 format_field = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

328 format_field_name = fields[8].split(':')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

329 format_data = fields[9].split(':')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

330

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

331 for i in range(0, len(format_data)):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

332 format_field[format_field_name[i]] = format_data[i].split(',')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

333

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

334 fields_to_store = {'REF': fields[3], 'ALT': fields[4].split(','), 'info': info_field, 'format': format_field}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

335 if fields[1] in variants:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

336 variants[fields[1]][len(variants[fields[1]])] = fields_to_store

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

337 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

338 variants[fields[1]] = {0: fields_to_store}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

339

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

340 line = vfc_file.readline()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

341 vfc_file.close()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

342

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

343 return variants

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

344

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

345

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

346 def indel_entry(variant_position):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

347 entry_with_indel = []

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

348 entry_with_snp = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

349 for i in variant_position:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

350 keys = variant_position[i]['info'].keys()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

351 if 'INDEL' in keys:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

352 entry_with_indel.append(i)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

353 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

354 entry_with_snp = i

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

355

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

356 return entry_with_indel, entry_with_snp

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

357

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

358

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

359 def get_alt_noMatter(variant_position, indel_true):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

360 dp = sum(map(int, variant_position['format']['AD']))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

361 index_alleles_sorted_position = sorted(zip(map(int, variant_position['format']['AD']), range(0, len(variant_position['format']['AD']))), reverse=True)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

362 index_dominant_allele = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

363 if not indel_true:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

364 ad_idv = index_alleles_sorted_position[0][0]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

365

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

366 if len([x for x in index_alleles_sorted_position if x[0] == ad_idv]) > 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

367 index_alleles_sorted_position = sorted([x for x in index_alleles_sorted_position if x[0] == ad_idv])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

368

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

369 index_dominant_allele = index_alleles_sorted_position[0][1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

370 if index_dominant_allele == 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

371 alt = '.'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

372 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

373 alt = variant_position['ALT'][index_dominant_allele - 1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

374

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

375 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

376 ad_idv = int(variant_position['info']['IDV'])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

377

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

378 if float(ad_idv) / float(dp) >= 0.5:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

379 if len([x for x in index_alleles_sorted_position if x[0] == index_alleles_sorted_position[0][0]]) > 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

380 index_alleles_sorted_position = sorted([x for x in index_alleles_sorted_position if x[0] == index_alleles_sorted_position[0][0]])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

381

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

382 index_dominant_allele = index_alleles_sorted_position[0][1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

383 if index_dominant_allele == 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

384 alt = '.'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

385 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

386 alt = variant_position['ALT'][index_dominant_allele - 1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

387 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

388 ad_idv = int(variant_position['format']['AD'][0])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

389 alt = '.'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

390

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

391 return alt, dp, ad_idv, index_dominant_allele

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

392

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

393

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

394 def count_number_diferences(ref, alt):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

395 number_diferences = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

396

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

397 if len(ref) != len(alt):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

398 number_diferences += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

399

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

400 for i in range(0, min(len(ref), len(alt))):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

401 if alt[i] != 'N' and ref[i] != alt[i]:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

402 number_diferences += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

403

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

404 return number_diferences

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

405

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

406

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

407 def get_alt_correct(variant_position, alt_noMatter, dp, ad_idv, index_dominant_allele, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

408 alt = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

409 low_coverage = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

410 multiple_alleles = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

411

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

412 if dp >= minimum_depth_presence:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

413 if dp < minimum_depth_call:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

414 alt = 'N' * len(variant_position['REF'])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

415 low_coverage = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

416 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

417 if ad_idv < minimum_depth_call:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

418 alt = 'N' * len(variant_position['REF'])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

419 low_coverage = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

420 if float(ad_idv) / float(dp) < minimum_depth_frequency_dominant_allele:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

421 multiple_alleles = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

422 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

423 if float(ad_idv) / float(dp) < minimum_depth_frequency_dominant_allele:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

424 alt = 'N' * len(variant_position['REF'])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

425 if index_dominant_allele is not None:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

426 variants_coverage = [int(variant_position['format']['AD'][i]) for i in range(0, len(variant_position['ALT']) + 1) if i != index_dominant_allele]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

427 if sum(variants_coverage) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

428 if float(max(variants_coverage)) / float(sum(variants_coverage)) > 0.5:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

429 multiple_alleles = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

430 elif float(max(variants_coverage)) / float(sum(variants_coverage)) == 0.5 and len(variants_coverage) > 2:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

431 multiple_alleles = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

432 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

433 multiple_alleles = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

434 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

435 alt = alt_noMatter

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

436 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

437 low_coverage = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

438

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

439 return alt, low_coverage, multiple_alleles

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

440

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

441

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

442 def get_alt_alignment(ref, alt):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

443 if alt is None:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

444 alt = 'N' * len(ref)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

445 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

446 if len(ref) != len(alt):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

447 if len(alt) < len(ref):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

448 if alt == '.':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

449 alt = ref

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

450 alt += 'N' * (len(ref) - len(alt))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

451 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

452 if alt[:len(ref)] == ref:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

453 alt = '.'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

454 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

455 alt = alt[:len(ref)]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

456

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

457 return alt

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

458

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

459

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

460 def get_indel_more_likely(variant_position, indels_entry):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

461 indel_coverage = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

462 for i in indels_entry:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

463 indel_coverage[i] = int(variant_position['info']['IDV'])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

464 return indel_coverage.index(str(max(indel_coverage.values())))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

465

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

466

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

467 def determine_variant(variant_position, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, indel_true):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

468 alt_noMatter, dp, ad_idv, index_dominant_allele = get_alt_noMatter(variant_position, indel_true)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

469

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

470 alt_correct, low_coverage, multiple_alleles = get_alt_correct(variant_position, alt_noMatter, dp, ad_idv, index_dominant_allele, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

471

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

472 alt_alignment = get_alt_alignment(variant_position['REF'], alt_correct)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

473

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

474 return variant_position['REF'], alt_correct, low_coverage, multiple_alleles, alt_noMatter, alt_alignment

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

475

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

476

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

477 def confirm_nucleotides_indel(ref, alt, variants, position_start_indel, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, alignment_true):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

478 alt = list(alt)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

479

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

480 for i in range(0, len(alt) - 1):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

481 if len(alt) < len(ref):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

482 new_position = position_start_indel + len(alt) - i - 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

483 alt_position = len(alt) - i - 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

484 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

485 if i + 1 > len(ref):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

486 break

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

487 new_position = position_start_indel + 1 + i

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

488 alt_position = 1 + i

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

489

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

490 if alt[alt_position] != 'N':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

491 if new_position not in variants:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

492 if alignment_true:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

493 alt[alt_position] = 'N'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

494 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

495 alt = alt[: alt_position]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

496 break

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

497

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

498 entry_with_indel, entry_with_snp = indel_entry(variants[new_position])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

499 new_ref, alt_correct, low_coverage, multiple_alleles, alt_noMatter, alt_alignment = determine_variant(variants[new_position][entry_with_snp], minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

500 if alt_noMatter != '.' and alt[alt_position] != alt_noMatter:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

501 alt[alt_position] = alt_noMatter

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

502

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

503 return ''.join(alt)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

504

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

505

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

506 def snp_indel(variants, position, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

507 entry_with_indel, entry_with_snp = indel_entry(variants[position])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

508

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

509 if len(entry_with_indel) == 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

510 ref, alt_correct, low_coverage, multiple_alleles, alt_noMatter, alt_alignment = determine_variant(variants[position][entry_with_snp], minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

511 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

512 ref_snp, alt_correct_snp, low_coverage_snp, multiple_alleles_snp, alt_noMatter_snp, alt_alignment_snp = determine_variant(variants[position][entry_with_snp], minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

513

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

514 indel_more_likely = entry_with_indel[0]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

515 if len(entry_with_indel) > 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

516 indel_more_likely = get_indel_more_likely(variants[position], entry_with_indel)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

517

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

518 ref, alt_correct, low_coverage, multiple_alleles, alt_noMatter, alt_alignment = determine_variant(variants[position][indel_more_likely], minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, True)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

519

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

520 if alt_noMatter == '.':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

521 ref, alt_correct, low_coverage, multiple_alleles, alt_noMatter, alt_alignment = ref_snp, alt_correct_snp, low_coverage_snp, multiple_alleles_snp, alt_noMatter_snp, alt_alignment_snp

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

522 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

523 if alt_correct is None and alt_correct_snp is not None:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

524 alt_correct = alt_correct_snp

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

525 elif alt_correct is not None and alt_correct_snp is not None:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

526 if alt_correct_snp != '.' and alt_correct[0] != alt_correct_snp:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

527 alt_correct = alt_correct_snp + alt_correct[1:] if len(alt_correct) > 1 else alt_correct_snp

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

528 if alt_noMatter_snp != '.' and alt_noMatter[0] != alt_noMatter_snp:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

529 alt_noMatter = alt_noMatter_snp + alt_noMatter[1:] if len(alt_noMatter) > 1 else alt_noMatter_snp

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

530 if alt_alignment_snp != '.' and alt_alignment[0] != alt_alignment_snp:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

531 alt_alignment = alt_alignment_snp + alt_alignment[1:] if len(alt_alignment) > 1 else alt_alignment_snp

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

532

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

533 # if alt_noMatter != '.':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

534 # alt_noMatter = confirm_nucleotides_indel(ref, alt_noMatter, variants, position, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

535 # if alt_correct is not None and alt_correct != '.':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

536 # alt_correct = confirm_nucleotides_indel(ref, alt_correct, variants, position, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

537 # if alt_alignment != '.':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

538 # alt_alignment = confirm_nucleotides_indel(ref, alt_alignment, variants, position, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, True)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

539

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

540 return ref, alt_correct, low_coverage, multiple_alleles, alt_noMatter, alt_alignment

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

541

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

542

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

543 def get_true_variants(variants, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, sequence):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

544 variants_correct = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

545 variants_noMatter = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

546 variants_alignment = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

547

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

548 correct_absent_positions = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

549 correct_last_absent_position = ''

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

550

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

551 noMatter_absent_positions = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

552 noMatter_last_absent_position = ''

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

553

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

554 multiple_alleles_found = []

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

555

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

556 counter = 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

557 while counter <= len(sequence):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

558 if counter in variants:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

559 noMatter_last_absent_position = ''

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

560

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

561 ref, alt_correct, low_coverage, multiple_alleles, alt_noMatter, alt_alignment = snp_indel(variants, counter, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

562

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

563 if alt_alignment != '.':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

564 variants_alignment[counter] = {'REF': ref, 'ALT': alt_alignment}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

565

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

566 if alt_noMatter != '.':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

567 variants_noMatter[counter] = {'REF': ref, 'ALT': alt_noMatter}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

568

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

569 if alt_correct is None:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

570 if counter - len(correct_last_absent_position) in correct_absent_positions:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

571 correct_absent_positions[counter - len(correct_last_absent_position)]['REF'] += ref

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

572 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

573 correct_absent_positions[counter] = {'REF': ref, 'ALT': ''}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

574 correct_last_absent_position += ref

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

575 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

576 if alt_correct != '.':

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

577 if len(alt_correct) < len(ref):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

578 if len(alt_correct) == 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

579 correct_absent_positions[counter + 1] = {'REF': ref[1:], 'ALT': ''}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

580 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

581 correct_absent_positions[counter + 1] = {'REF': ref[1:], 'ALT': alt_correct[1:]}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

582

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

583 correct_last_absent_position = ref[1:]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

584 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

585 variants_correct[counter] = {'REF': ref, 'ALT': alt_correct}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

586 correct_last_absent_position = ''

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

587 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

588 correct_last_absent_position = ''

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

589

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

590 if multiple_alleles:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

591 multiple_alleles_found.append(counter)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

592

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

593 counter += len(ref)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

594 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

595 variants_alignment[counter] = {'REF': sequence[counter - 1], 'ALT': 'N'}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

596

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

597 if counter - len(correct_last_absent_position) in correct_absent_positions:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

598 correct_absent_positions[counter - len(correct_last_absent_position)]['REF'] += sequence[counter - 1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

599 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

600 correct_absent_positions[counter] = {'REF': sequence[counter - 1], 'ALT': ''}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

601 correct_last_absent_position += sequence[counter - 1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

602

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

603 if counter - len(noMatter_last_absent_position) in noMatter_absent_positions:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

604 noMatter_absent_positions[counter - len(noMatter_last_absent_position)]['REF'] += sequence[counter - 1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

605 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

606 noMatter_absent_positions[counter] = {'REF': sequence[counter - 1], 'ALT': ''}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

607 noMatter_last_absent_position += sequence[counter - 1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

608

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

609 counter += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

610

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

611 for position in correct_absent_positions:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

612 if position == 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

613 variants_correct[position] = {'REF': correct_absent_positions[position]['REF'], 'ALT': 'N'}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

614 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

615 if position - 1 not in variants_correct:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

616 variants_correct[position - 1] = {'REF': sequence[position - 2] + correct_absent_positions[position]['REF'], 'ALT': sequence[position - 2] + correct_absent_positions[position]['ALT']}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

617 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

618 variants_correct[position - 1] = {'REF': variants_correct[position - 1]['REF'] + correct_absent_positions[position]['REF'][len(variants_correct[position - 1]['REF']) - 1:], 'ALT': variants_correct[position - 1]['ALT'] + correct_absent_positions[position]['ALT'][len(variants_correct[position - 1]['ALT']) - 1 if len(variants_correct[position - 1]['ALT']) > 0 else 0:]}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

619

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

620 for position in noMatter_absent_positions:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

621 if position == 1:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

622 variants_noMatter[position] = {'REF': noMatter_absent_positions[position]['REF'], 'ALT': 'N'}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

623 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

624 if position - 1 not in variants_noMatter:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

625 variants_noMatter[position - 1] = {'REF': sequence[position - 2] + noMatter_absent_positions[position]['REF'], 'ALT': sequence[position - 2] + noMatter_absent_positions[position]['ALT']}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

626 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

627 variants_noMatter[position - 1] = {'REF': variants_noMatter[position - 1]['REF'] + noMatter_absent_positions[position]['REF'][len(variants_noMatter[position - 1]['REF']) - 1:], 'ALT': variants_noMatter[position - 1]['ALT'] + noMatter_absent_positions[position]['ALT'][len(variants_noMatter[position - 1]['ALT']) - 1 if len(variants_noMatter[position - 1]['ALT']) > 0 else 0:]}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

628

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

629 return variants_correct, variants_noMatter, variants_alignment, multiple_alleles_found

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

630

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

631

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

632 def clean_variant_in_extra_seq_left(variant_dict, position, length_extra_seq, multiple_alleles_found, number_multi_alleles):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

633 number_diferences = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

634

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

635 if position + len(variant_dict[position]['REF']) - 1 > length_extra_seq:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

636 if multiple_alleles_found is not None and position in multiple_alleles_found:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

637 number_multi_alleles += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

638

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

639 temp_variant = variant_dict[position]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

640 del variant_dict[position]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

641 variant_dict[length_extra_seq] = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

642 variant_dict[length_extra_seq]['REF'] = temp_variant['REF'][length_extra_seq - position:]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

643 variant_dict[length_extra_seq]['ALT'] = temp_variant['ALT'][length_extra_seq - position:] if len(temp_variant['ALT']) > length_extra_seq - position else temp_variant['REF'][length_extra_seq - position]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

644 number_diferences = count_number_diferences(variant_dict[length_extra_seq]['REF'], variant_dict[length_extra_seq]['ALT'])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

645 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

646 del variant_dict[position]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

647

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

648 return variant_dict, number_multi_alleles, number_diferences

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

649

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

650

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

651 def clean_variant_in_extra_seq_rigth(variant_dict, position, sequence_length, length_extra_seq):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

652 if position + len(variant_dict[position]['REF']) - 1 > sequence_length - length_extra_seq:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

653 variant_dict[position]['REF'] = variant_dict[position]['REF'][: - (position - (sequence_length - length_extra_seq)) + 1]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

654 variant_dict[position]['ALT'] = variant_dict[position]['ALT'][: - (position - (sequence_length - length_extra_seq)) + 1] if len(variant_dict[position]['ALT']) >= - (position - (sequence_length - length_extra_seq)) + 1 else variant_dict[position]['ALT']

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

655

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

656 number_diferences = count_number_diferences(variant_dict[position]['REF'], variant_dict[position]['ALT'])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

657

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

658 return variant_dict, number_diferences

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

659

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

660

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

661 def cleanning_variants_extra_seq(variants_correct, variants_noMatter, variants_alignment, multiple_alleles_found, length_extra_seq, sequence_length):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

662 number_multi_alleles = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

663 number_diferences = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

664

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

665 counter = 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

666 while counter <= sequence_length:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

667 if counter <= length_extra_seq:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

668 if counter in variants_correct:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

669 variants_correct, number_multi_alleles, number_diferences = clean_variant_in_extra_seq_left(variants_correct, counter, length_extra_seq, multiple_alleles_found, number_multi_alleles)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

670 if counter in variants_noMatter:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

671 variants_noMatter, ignore, ignore = clean_variant_in_extra_seq_left(variants_noMatter, counter, length_extra_seq, None, None)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

672 if counter in variants_alignment:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

673 variants_alignment, ignore, ignore = clean_variant_in_extra_seq_left(variants_alignment, counter, length_extra_seq, None, None)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

674 elif counter > length_extra_seq and counter <= sequence_length - length_extra_seq:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

675 if counter in variants_correct:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

676 if counter in multiple_alleles_found:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

677 number_multi_alleles += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

678 variants_correct, number_diferences_found = clean_variant_in_extra_seq_rigth(variants_correct, counter, sequence_length, length_extra_seq)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

679 number_diferences += number_diferences_found

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

680 if counter in variants_noMatter:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

681 variants_noMatter, ignore = clean_variant_in_extra_seq_rigth(variants_noMatter, counter, sequence_length, length_extra_seq)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

682 if counter in variants_alignment:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

683 variants_alignment, ignore = clean_variant_in_extra_seq_rigth(variants_alignment, counter, sequence_length, length_extra_seq)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

684 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

685 if counter in variants_correct:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

686 del variants_correct[counter]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

687 if counter in variants_noMatter:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

688 del variants_noMatter[counter]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

689 if counter in variants_alignment:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

690 del variants_alignment[counter]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

691

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

692 counter += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

693

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

694 return variants_correct, variants_noMatter, variants_alignment, number_multi_alleles, number_diferences

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

695

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

696

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

697 def get_coverage(gene_coverage):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

698 coverage = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

699

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

700 with open(gene_coverage, 'rtU') as reader:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

701 for line in reader:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

702 line = line.splitlines()[0]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

703 if len(line) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

704 line = line.split('\t')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

705 coverage[int(line[1])] = int(line[2])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

706

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

707 return coverage

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

708

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

709

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

710 def get_coverage_report(coverage, sequence_length, minimum_depth_presence, minimum_depth_call, length_extra_seq):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

711 if len(coverage) == 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

712 return sequence_length - 2 * length_extra_seq, 100.0, 0.0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

713

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

714 count_absent = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

715 count_lowCoverage = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

716 sum_coverage = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

717

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

718 counter = 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

719 while counter <= sequence_length:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

720 if counter > length_extra_seq and counter <= sequence_length - length_extra_seq:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

721 if coverage[counter] < minimum_depth_presence:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

722 count_absent += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

723 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

724 if coverage[counter] < minimum_depth_call:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

725 count_lowCoverage += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

726 sum_coverage += coverage[counter]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

727 counter += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

728

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

729 mean_coverage = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

730 percentage_lowCoverage = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

731 if sequence_length - 2 * length_extra_seq - count_absent > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

732 mean_coverage = float(sum_coverage) / float(sequence_length - 2 * length_extra_seq - count_absent)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

733 percentage_lowCoverage = float(count_lowCoverage) / float(sequence_length - 2 * length_extra_seq - count_absent) * 100

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

734

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

735 return count_absent, percentage_lowCoverage, mean_coverage

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

736

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

737

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

738 # Get genome coverage data

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

739 def compute_genome_coverage_data(alignment_file, sequence_to_analyse, outdir, counter):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

740 genome_coverage_data_file = os.path.join(outdir, 'samtools_depth.sequence_' + str(counter) + '.tab')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

741 command = ['samtools', 'depth', '-a', '-q', '0', '-r', sequence_to_analyse, alignment_file, '>', genome_coverage_data_file]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

742 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, True, None, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

743 return run_successfully, genome_coverage_data_file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

744

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

745

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

746 def write_variants_vcf(variants, outdir, sequence_to_analyse, sufix):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

747 vcf_file = os.path.join(outdir, str(sequence_to_analyse + '.' + sufix + '.vcf'))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

748 with open(vcf_file, 'wt') as writer:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

749 writer.write('##fileformat=VCFv4.2' + '\n')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

750 writer.write('#' + '\t'.join(['SEQUENCE', 'POSITION', 'ID_unused', 'REFERENCE_sequence', 'ALTERNATIVE_sequence', 'QUALITY_unused', 'FILTER_unused', 'INFO_unused', 'FORMAT_unused']) + '\n')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

751 for i in sorted(variants.keys()):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

752 writer.write('\t'.join([sequence_to_analyse, str(i), '.', variants[i]['REF'], variants[i]['ALT'], '.', '.', '.', '.']) + '\n')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

753

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

754 compressed_vcf_file = vcf_file + '.gz'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

755 command = ['bcftools', 'convert', '-o', compressed_vcf_file, '-O', 'z', vcf_file]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

756 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

757 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

758 command = ['bcftools', 'index', compressed_vcf_file]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

759 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

760

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

761 if not run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

762 compressed_vcf_file = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

763

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

764 return run_successfully, compressed_vcf_file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

765

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

766

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

767 def parse_fasta_inMemory(fasta_memory):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

768 fasta_memory = fasta_memory.splitlines()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

769 sequence_dict = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

770 for line in fasta_memory:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

771 if len(line) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

772 if line.startswith('>'):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

773 sequence_dict = {'header': line[1:], 'sequence': ''}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

774 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

775 sequence_dict['sequence'] += line

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

776

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

777 return sequence_dict

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

778

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

779

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

780 def compute_consensus_sequence(reference_file, sequence_to_analyse, compressed_vcf_file, outdir, sufix):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

781 sequence_dict = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

782

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

783 gene_fasta = os.path.join(outdir, str(sequence_to_analyse + '.fasta'))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

784

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

785 run_successfully, stdout = index_fasta_samtools(reference_file, sequence_to_analyse, gene_fasta, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

786 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

787 command = ['bcftools', 'norm', '-c', 's', '-f', gene_fasta, '-Ov', compressed_vcf_file]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

788 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

789 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

790 command = ['bcftools', 'consensus', '-f', gene_fasta, compressed_vcf_file, '-H', '1']

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

791 run_successfully, stdout, stderr = utils.runCommandPopenCommunicate(command, False, None, False)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

792 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

793 sequence_dict = parse_fasta_inMemory(stdout)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

794

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

795 return run_successfully, sequence_dict

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

796

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

797

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

798 def create_sample_consensus_sequence(outdir, sequence_to_analyse, reference_file, variants, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, sequence, length_extra_seq):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

799 variants_correct, variants_noMatter, variants_alignment, multiple_alleles_found = get_true_variants(variants, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, sequence)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

800

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

801 variants_correct, variants_noMatter, variants_alignment, number_multi_alleles, number_diferences = cleanning_variants_extra_seq(variants_correct, variants_noMatter, variants_alignment, multiple_alleles_found, length_extra_seq, len(sequence))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

802

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

803 run_successfully = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

804 consensus = {'correct': {}, 'noMatter': {}, 'alignment': {}}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

805 for variant_type in ['variants_correct', 'variants_noMatter', 'variants_alignment']:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

806 run_successfully, compressed_vcf_file = write_variants_vcf(eval(variant_type), outdir, sequence_to_analyse, variant_type.split('_', 1)[1])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

807 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

808 run_successfully, sequence_dict = compute_consensus_sequence(reference_file, sequence_to_analyse, compressed_vcf_file, outdir, variant_type.split('_', 1)[1])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

809 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

810 consensus[variant_type.split('_', 1)[1]] = {'header': sequence_dict['header'], 'sequence': sequence_dict['sequence'][length_extra_seq:len(sequence_dict['sequence']) - length_extra_seq]}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

811

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

812 return run_successfully, number_multi_alleles, consensus, number_diferences

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

813

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

814

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

815 @utils.trace_unhandled_exceptions

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

816 def analyse_sequence_data(bam_file, sequence_information, outdir, counter, reference_file, length_extra_seq, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, threads):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

817 count_absent = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

818 percentage_lowCoverage = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

819 meanCoverage = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

820 number_diferences = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

821

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

822 # Create vcf file (for multiple alleles check)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

823 run_successfully, gene_vcf = create_vcf(bam_file, sequence_information['header'], outdir, counter, reference_file)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

824 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

825 # Create coverage tab file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

826 run_successfully, gene_coverage = compute_genome_coverage_data(bam_file, sequence_information['header'], outdir, counter)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

827

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

828 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

829 variants = get_variants(gene_vcf)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

830

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

831 coverage = get_coverage(gene_coverage)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

832

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

833 run_successfully, number_multi_alleles, consensus_sequence, number_diferences = create_sample_consensus_sequence(outdir, sequence_information['header'], reference_file, variants, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, sequence_information['sequence'], length_extra_seq)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

834

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

835 count_absent, percentage_lowCoverage, meanCoverage = get_coverage_report(coverage, sequence_information['length'], minimum_depth_presence, minimum_depth_call, length_extra_seq)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

836

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

837 utils.saveVariableToPickle([run_successfully, counter, number_multi_alleles, count_absent, percentage_lowCoverage, meanCoverage, consensus_sequence, number_diferences], outdir, str('coverage_info.' + str(counter)))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

838

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

839

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

840 def clean_header(header):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

841 problematic_characters = ["|", " ", ",", ".", "(", ")", "'", "/", ":"]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

842 new_header = str(header)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

843 if any(x in header for x in problematic_characters):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

844 for x in problematic_characters:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

845 new_header = new_header.replace(x, '_')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

846 return header, new_header

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

847

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

848

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

849 def get_sequence_information(fasta_file, length_extra_seq):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

850 sequence_dict = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

851 headers = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

852 headers_changed = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

853

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

854 with open(fasta_file, 'rtU') as reader:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

855 blank_line_found = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

856 sequence_counter = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

857 temp_sequence_dict = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

858 for line in reader:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

859 line = line.splitlines()[0]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

860 if len(line) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

861 if not blank_line_found:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

862 if line.startswith('>'):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

863 if len(temp_sequence_dict) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

864 if temp_sequence_dict.values()[0]['length'] - 2 * length_extra_seq > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

865 sequence_dict[temp_sequence_dict.keys()[0]] = temp_sequence_dict.values()[0]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

866 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

867 print '{header} sequence ignored due to length <= 0'.format(header=temp_sequence_dict.values()[0]['header'])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

868 del headers[temp_sequence_dict.values()[0]['header']]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

869 temp_sequence_dict = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

870

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

871 original_header, new_header = clean_header(line[1:])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

872 if new_header in headers:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

873 sys.exit('Found duplicated sequence headers: {original_header}'.format(original_header=original_header))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

874

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

875 sequence_counter += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

876 temp_sequence_dict[sequence_counter] = {'header': new_header, 'sequence': '', 'length': 0}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

877 headers[new_header] = str(original_header)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

878 if new_header != original_header:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

879 headers_changed = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

880 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

881 temp_sequence_dict[sequence_counter]['sequence'] += line

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

882 temp_sequence_dict[sequence_counter]['length'] += len(line)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

883 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

884 sys.exit('It was found a blank line between the fasta file above line ' + line)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

885 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

886 blank_line_found = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

887

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

888 if len(temp_sequence_dict) > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

889 if temp_sequence_dict.values()[0]['length'] - 2 * length_extra_seq > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

890 sequence_dict[temp_sequence_dict.keys()[0]] = temp_sequence_dict.values()[0]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

891 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

892 print '{header} sequence ignored due to length <= 0'.format(header=temp_sequence_dict.values()[0]['header'])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

893 del headers[temp_sequence_dict.values()[0]['header']]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

894

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

895 return sequence_dict, headers, headers_changed

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

896

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

897

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

898 def determine_threads_2_use(number_sequences, threads):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

899 if number_sequences >= threads:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

900 return 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

901 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

902 return threads / number_sequences

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

903

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

904

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

905 def sequence_data(sample, reference_file, bam_file, outdir, threads, length_extra_seq, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, debug_mode_true, notWriteConsensus):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

906 sequence_data_outdir = os.path.join(outdir, 'sequence_data', '')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

907 utils.removeDirectory(sequence_data_outdir)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

908 os.mkdir(sequence_data_outdir)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

909

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

910 sequences, headers, headers_changed = get_sequence_information(reference_file, length_extra_seq)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

911

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

912 threads_2_use = determine_threads_2_use(len(sequences), threads)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

913

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

914 pool = multiprocessing.Pool(processes=threads)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

915 for sequence_counter in sequences:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

916 sequence_dir = os.path.join(sequence_data_outdir, str(sequence_counter), '')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

917 utils.removeDirectory(sequence_dir)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

918 os.makedirs(sequence_dir)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

919 pool.apply_async(analyse_sequence_data, args=(bam_file, sequences[sequence_counter], sequence_dir, sequence_counter, reference_file, length_extra_seq, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, threads_2_use,))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

920 pool.close()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

921 pool.join()

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

922

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

923 run_successfully, sample_data, consensus_files, consensus_sequences = gather_data_together(sample, sequence_data_outdir, sequences, outdir.rsplit('/', 2)[0], debug_mode_true, length_extra_seq, notWriteConsensus)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

924

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

925 return run_successfully, sample_data, consensus_files, consensus_sequences

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

926

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

927

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

928 def chunkstring(string, length):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

929 return (string[0 + i:length + i] for i in range(0, len(string), length))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

930

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

931

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

932 def write_consensus(outdir, sample, consensus_sequence):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

933 consensus_files = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

934 for consensus_type in ['correct', 'noMatter', 'alignment']:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

935 consensus_files[consensus_type] = os.path.join(outdir, str(sample + '.' + consensus_type + '.fasta'))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

936 with open(consensus_files[consensus_type], 'at') as writer:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

937 writer.write('>' + consensus_sequence[consensus_type]['header'] + '\n')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

938 fasta_sequence_lines = chunkstring(consensus_sequence[consensus_type]['sequence'], 80)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

939 for line in fasta_sequence_lines:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

940 writer.write(line + '\n')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

941 return consensus_files

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

942

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

943

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

944 def gather_data_together(sample, data_directory, sequences_information, outdir, debug_mode_true, length_extra_seq, notWriteConsensus):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

945 run_successfully = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

946 counter = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

947 sample_data = {}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

948

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

949 consensus_files = None

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

950 consensus_sequences_together = {'correct': {}, 'noMatter': {}, 'alignment': {}}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

951

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

952 write_consensus_first_time = True

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

953

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

954 genes_directories = [d for d in os.listdir(data_directory) if not d.startswith('.') and os.path.isdir(os.path.join(data_directory, d, ''))]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

955 for gene_dir in genes_directories:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

956 gene_dir_path = os.path.join(data_directory, gene_dir, '')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

957

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

958 files = [f for f in os.listdir(gene_dir_path) if not f.startswith('.') and os.path.isfile(os.path.join(gene_dir_path, f))]

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

959 for file_found in files:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

960 if file_found.startswith('coverage_info.') and file_found.endswith('.pkl'):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

961 file_path = os.path.join(gene_dir_path, file_found)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

962

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

963 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

964 run_successfully, sequence_counter, multiple_alleles_found, count_absent, percentage_lowCoverage, meanCoverage, consensus_sequence, number_diferences = utils.extractVariableFromPickle(file_path)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

965

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

966 if not notWriteConsensus:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

967 for consensus_type in consensus_sequence:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

968 consensus_sequences_together[consensus_type][sequence_counter] = {'header': consensus_sequence[consensus_type]['header'], 'sequence': consensus_sequence[consensus_type]['sequence']}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

969

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

970 if write_consensus_first_time:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

971 for consensus_type in ['correct', 'noMatter', 'alignment']:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

972 file_to_remove = os.path.join(outdir, str(sample + '.' + consensus_type + '.fasta'))

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

973 if os.path.isfile(file_to_remove):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

974 os.remove(file_to_remove)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

975 write_consensus_first_time = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

976 consensus_files = write_consensus(outdir, sample, consensus_sequence)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

977

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

978 gene_identity = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

979 if sequences_information[sequence_counter]['length'] - 2 * length_extra_seq - count_absent > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

980 gene_identity = 100 - (float(number_diferences) / (sequences_information[sequence_counter]['length'] - 2 * length_extra_seq - count_absent)) * 100

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

981

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

982 sample_data[sequence_counter] = {'header': sequences_information[sequence_counter]['header'], 'gene_coverage': 100 - (float(count_absent) / (sequences_information[sequence_counter]['length'] - 2 * length_extra_seq)) * 100, 'gene_low_coverage': percentage_lowCoverage, 'gene_number_positions_multiple_alleles': multiple_alleles_found, 'gene_mean_read_coverage': meanCoverage, 'gene_identity': gene_identity}

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

983 counter += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

984

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

985 if not debug_mode_true:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

986 utils.removeDirectory(gene_dir_path)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

987

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

988 if counter != len(sequences_information):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

989 run_successfully = False

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

990

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

991 return run_successfully, sample_data, consensus_files, consensus_sequences_together

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

992

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

993

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

994 rematch_timer = functools.partial(utils.timer, name='ReMatCh module')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

995

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

996

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

997 @rematch_timer

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

998 def runRematchModule(sample, fastq_files, reference_file, threads, outdir, length_extra_seq, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, minimum_gene_coverage, conserved_True, debug_mode_true, numMapLoc, minimum_gene_identity, rematch_run, softClip_baseQuality, softClip_recodeRun, reference_dict, softClip_cigarFlagRecode, bowtieOPT, gene_list_reference, notWriteConsensus):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

999 rematch_folder = os.path.join(outdir, 'rematch_module', '')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1000 utils.removeDirectory(rematch_folder)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1001 os.mkdir(rematch_folder)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1002

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1003 # Map reads

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1004 run_successfully, bam_file, reference_file = mapping_reads(fastq_files, reference_file, threads, rematch_folder, conserved_True, numMapLoc, rematch_run, softClip_baseQuality, softClip_recodeRun, reference_dict, softClip_cigarFlagRecode, bowtieOPT)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1005

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1006 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1007 # Index reference file

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1008 run_successfully, stdout = index_fasta_samtools(reference_file, None, None, True)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1009 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1010 print 'Analysing alignment data'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1011 run_successfully, sample_data, consensus_files, consensus_sequences = sequence_data(sample, reference_file, bam_file, rematch_folder, threads, length_extra_seq, minimum_depth_presence, minimum_depth_call, minimum_depth_frequency_dominant_allele, debug_mode_true, notWriteConsensus)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1012

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1013 if run_successfully:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1014 print 'Writing report file'

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1015 number_absent_genes = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1016 number_genes_multiple_alleles = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1017 mean_sample_coverage = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1018 with open(os.path.join(outdir, 'rematchModule_report.txt'), 'wt') as writer:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1019 writer.write('\t'.join(['#gene', 'percentage_gene_coverage', 'gene_mean_read_coverage', 'percentage_gene_low_coverage', 'number_positions_multiple_alleles', 'percentage_gene_identity']) + '\n')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1020 for i in range(1, len(sample_data) + 1):

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1021 writer.write('\t'.join([gene_list_reference[sample_data[i]['header']], str(round(sample_data[i]['gene_coverage'], 2)), str(round(sample_data[i]['gene_mean_read_coverage'], 2)), str(round(sample_data[i]['gene_low_coverage'], 2)), str(sample_data[i]['gene_number_positions_multiple_alleles']), str(round(sample_data[i]['gene_identity'], 2))]) + '\n')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1022

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1023 if sample_data[i]['gene_coverage'] < minimum_gene_coverage or sample_data[i]['gene_identity'] < minimum_gene_identity:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1024 number_absent_genes += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1025 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1026 mean_sample_coverage += sample_data[i]['gene_mean_read_coverage']

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1027 if sample_data[i]['gene_number_positions_multiple_alleles'] > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1028 number_genes_multiple_alleles += 1

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1029

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1030 if len(sample_data) - number_absent_genes > 0:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1031 mean_sample_coverage = float(mean_sample_coverage) / float(len(sample_data) - number_absent_genes)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1032 else:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1033 mean_sample_coverage = 0

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1034

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1035 writer.write('\n'.join(['#general', '>number_absent_genes', str(number_absent_genes), '>number_genes_multiple_alleles', str(number_genes_multiple_alleles), '>mean_sample_coverage', str(round(mean_sample_coverage, 2))]) + '\n')

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1036

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1037 print '\n'.join([str('number_absent_genes: ' + str(number_absent_genes)), str('number_genes_multiple_alleles: ' + str(number_genes_multiple_alleles)), str('mean_sample_coverage: ' + str(round(mean_sample_coverage, 2)))])

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1038

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1039 if not debug_mode_true:

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1040 utils.removeDirectory(rematch_folder)

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1041

965517909457 planemo upload commit 15239f1674081ab51ab8dd75a9a40cf1bfaa93e8

cstrittmatter

parents:

diff changeset

1042 return run_successfully, sample_data if 'sample_data' in locals() else None, {'number_absent_genes': number_absent_genes if 'number_absent_genes' in locals() else None, 'number_genes_multiple_alleles': number_genes_multiple_alleles if 'number_genes_multiple_alleles' in locals() else None, 'mean_sample_coverage': round(mean_sample_coverage, 2) if 'mean_sample_coverage' in locals() else None}, consensus_files if 'consensus_files' in locals() else None, consensus_sequences if 'consensus_sequences' in locals() else None

Mercurial > repos > cstrittmatter > test_eurl_vtec_wgs_pt

annotate scripts/ReMatCh/modules/rematch_module.py @ 0:965517909457 draft