bugs

README.md

@@ -27,8 +27,6 @@ June 2021
 5. [License](#license)
 
 
-
-
 ## <a name="install"></a>Installation
 *detettore* is written in Python 3 and available on PyPI. To avoid conflicts
 with dependencies, it is best to install it in a virtual environment:
@@ -153,9 +151,6 @@ parallel -j10 < run_detettore.cmds 2> err.log > stdout.log
 
 
 
-
-
-
 ## <a name="version"></a>**New in version 2**
 - all output in vcf format
 - output invariant sites

detettore.py

@@ -195,7 +195,6 @@ def main():
     elif "taps" in parameters.modus and parameters.paired_end:
         combined_vcf = taps_vcf
 
-
     date = time.strftime("%d/%m/%Y")
 
     metainfo = [
@@ -270,7 +269,6 @@ def main():
 
     log.info(args)
 
-
     # Clean up
     if not args.keep:
         shutil.rmtree(args.outname + '_tmp')

scripts/strumenti.py

@@ -335,7 +335,6 @@ def output_vcf(self, parameters):
             POS = site[1]
             REF = parameters.ref_contigs[CHROM][POS - 1]
 
-
             discordant = site[2]
             split = site[3]
 
@@ -558,10 +557,10 @@ def output_vcf(self, parameters):
                 INFO = 'MEINFO=%s;SVTYPE=INS;HOMSEQ=%s;HOMLEN=%i;DPADJ=%i;DR=%i;SR=%i;AL=%i;IMPRECISE;CIPOS=%i,%i;AP=%f;BPIQR=%i' \
                     % (MEINFO, TSD, TSDLEN, DPADJ, DR, SR, AL, CI_lower, CI_upper, round(PROP, 2), IQR)
 
-
             FORMAT = 'GT:GQ:AD:DP:PL'
-
-            PL = [ int(x) if not math.isinf(x) else x for x in genotype[2] ]
+
+            # Max GT likelihood of 1e6 (to avoid inf values)
+            PL = [ str(int(x)) if not math.isinf(x) else '1e6' for x in genotype[2] ]
             GT = '%s:%i:%s:%i:%s' % (genotype[0], genotype[1], AD, DP, ','.join(PL))
 
             outline = [CHROM, POS, '.', REF, ALT, '.', 'PASS', INFO, FORMAT, GT]
@@ -788,8 +787,7 @@ def __init__(self, parameters):
 
             else:
                 inputs.append((te, parameters))
-
-
+
         with mp.Pool(parameters.cpus) as pool:
             self.candidates = pool.map(self.process_tap, inputs)
 
@@ -1133,7 +1131,9 @@ def extract_deviant_reads(te,  region, parameters):
         Returns a dictionary d[read name] = [isize, AlignedSegment objects]
         
         IMPROVE: 
-            Only look up reads until start position of TE
+            - only look up reads until start position of TE
+            - alternative strategy for SE reads: count hard clipped reads as 
+            evidence for a TAP?
 
         """