Merge pull request #41 from bcgsc/tigmint-long-bugs

Fix barcode naming bug

Author: janetxinli

README.md

@@ -121,7 +121,7 @@ tigmint-make metrics draft=myassembly reads=myreads ref=GRCh38 G=3088269832
 To run Tigmint with long reads in fasta or fastq format (`reads.fa.gz` or `reads.fq.gz`) on the draft assembly `draft.fa`:
 
 ```sh
-tigmint-make tigmint-long draft=myassembly reads=myreads span=auto gsize=genomesize
+tigmint-make tigmint-long draft=myassembly reads=myreads span=auto G=genomesize
 ```
 
 - `minimap2 map-ont` is used to align long reads from the Oxford Nanopore Technologies (ONT) platform, which is the default input for Tigmint. To use PacBio long reads specify the parameter `longmap=pb`
@@ -131,7 +131,7 @@ tigmint-make tigmint-long draft=myassembly reads=myreads span=auto gsize=genomes
 
 + `tigmint-make` is a Makefile script, and so any `make` options may also be used with `tigmint-make`, such as `-n` (`--dry-run`).
 + The file extension of the assembly must be `.fa` and the reads `.fq.gz` (or `.fa.gz` for long reads), and the extension is not included in the parameters `draft` and `reads`. These specific file name requirements result from implementing the pipeline in GNU Make.
-+ The minimum spanning molecules parameter (`span`) for `tigmint-cut` is heavily dependent on the sequence coverage of the linked or long reads provided. When running Tigmint with long reads, use `span=auto` and set `gsize` to your assembly organism's haploid genome size for this parameter to be calculated automatically, or explicitly set `span` to a specific number if you are interested in adjusting it. See [Tips](https://github.com/bcgsc/tigmint/tree/long-to-linked/#Note) for more details. 
++ The minimum spanning molecules parameter (`span`) for `tigmint-cut` is heavily dependent on the sequence coverage of the linked or long reads provided. When running Tigmint with long reads, use `span=auto` and set `G` to your assembly organism's haploid genome size for this parameter to be calculated automatically, or explicitly set `span` to a specific number if you are interested in adjusting it. See [Tips](https://github.com/bcgsc/tigmint/tree/long-to-linked/#Note) for more details. 
 
 # tigmint-make commands
 
@@ -178,7 +178,7 @@ tigmint-make tigmint-long draft=myassembly reads=myreads span=auto gsize=genomes
 - Sort by BX tag using `samtools sort -tBX`.
 - Merge multiple BAM files using `samtools merge -tBX`.
 - When aligning long reads with Minimap2, use the `-y` option to include the barcode in the BX tag of the alignments.
-- When using long reads, the minimum spanning molecule thresholds (`span`) should be no greater than 1/4 of the sequence coverage. Setting the parameter `span=auto` allows the appropriate parameter value to be selected automatically (this setting requires the parameter `gsize` as well).
+- When using long reads, the minimum spanning molecule thresholds (`span`) should be no greater than 1/4 of the sequence coverage. Setting the parameter `span=auto` allows the appropriate parameter value to be selected automatically (this setting requires the parameter `G` as well).
 - When using long reads, the edit distance threshold (`nm`) is automatically set to the cut length (`cut`) to compensate for the higher error rate and length. This parameter should be kept relatively high to include as many alignments as possible.
 
 # Using stLFR linked reads

azure-pipelines.yml

@@ -27,7 +27,8 @@ jobs:
       echo path:
       echo $PATH
       sudo /home/linuxbrew/.linuxbrew/bin/pip3 install pytest
-      sudo HOMEBREW_NO_AUTOUPDATE=1 /home/linuxbrew/.linuxbrew/bin/brew install bedtools bwa samtools gnu-time
+      sudo HOMEBREW_NO_AUTOUPDATE=1 /home/linuxbrew/.linuxbrew/bin/brew tap brewsci/bio 
+      sudo HOMEBREW_NO_AUTOUPDATE=1 /home/linuxbrew/.linuxbrew/bin/brew install bedtools bwa samtools gnu-time minimap2
       cd tests/
       pytest -vs tigmint_test.py
     env:

bin/long-to-linked

@@ -23,15 +23,17 @@ def cut_reads(length, min_size, auto_span, genome_size, span_out):
     total_bases = 0
 
     with open("/dev/stdin", 'rt') as reads:
+        barcode = 0
         if not reads.isatty():  # Check that reads are being piped in from stdin
             for header, seq, _, _ in read_fasta(reads):
                 read_length = len(seq)
                 total_bases += read_length
                 if read_length >= min_size:
                     split_reads = split_long_read(seq, length)
                     for i, read in enumerate(split_reads):
-                        print(">" + header + "_" + str(i) + " BX:Z:" + str(i), file=sys.stdout)
+                        print(">" + header + "_" + str(i) + " BX:Z:" + str(barcode), file=sys.stdout)
                         print(read, file=sys.stdout)
+                    barcode += 1
             if auto_span:
                 with open(span_out, "w") as span_file:
                     span = int(total_bases / genome_size * cov_to_span)
@@ -46,7 +48,7 @@ def main():
     parser = argparse.ArgumentParser(description="Segment long reads into short reads and calculate span.")
     parser.add_argument("-v", "--version",
                         action="version",
-                        version="long-to-linked 1.2.0")
+                        version="long-to-linked 1.2.1")
     parser.add_argument("-l", "--length",
                         type=int,
                         default=500,

bin/tigmint-cut

@@ -238,7 +238,7 @@ def get_span(filename):
 
 def main():
     parser = argparse.ArgumentParser(description="Find misassembled regions in assembly using Chromium molecule extents")
-    parser.add_argument("--version", action="version", version="tigmint-cut 1.2.0")
+    parser.add_argument("--version", action="version", version="tigmint-cut 1.2.1")
     parser.add_argument("fasta", type=str, help="Reference genome fasta file (must have FAI index generated)")
     parser.add_argument("bed", type=str, help="Sorted bed file of molecule extents")
     parser.add_argument("-o", "--fastaout", type=str, help="The output FASTA file.", required=True)

bin/tigmint-make

@@ -123,7 +123,7 @@ help:
 	@echo 'For more information see https://bcgsc.github.io/tigmint/'
 
 version:
-	@echo "Tigmint 1.2.0"
+	@echo "Tigmint 1.2.1"
 	@echo "Written by Shaun Jackman @sjackman."
 
 all: tigmint arcs

bin/tigmint-molecule

@@ -230,7 +230,7 @@ class MolecIdentifier:
             "Read a SAM/BAM file and output a TSV file. "
             "The SAM/BAM file must be sorted by BX tag and then by position.")
         parser.add_argument(
-            '--version', action='version', version='tigmint-molecule 1.2.0')
+            '--version', action='version', version='tigmint-molecule 1.2.1')
         parser.add_argument(
             metavar="BAM", dest="in_bam_filename",
             help="Input BAM file sorted by BX tag then position, - for stdin")

setup.py

@@ -5,7 +5,7 @@
 
 setuptools.setup(
     name="tigmint",
-    version="1.2.0",
+    version="1.2.1",
     author="Shaun Jackman",
     author_email="[email protected]",
     description="Correct misassemblies using linked or long reads",

tests/span_gsize_100000 tests/span_G_100000tests/span_gsize_100000 renamed to tests/span_G_100000

@@ -1,26 +1,24 @@
-test	0	4998	14	4
-test	3	4998	2	7
-test	4	4003	16	6
-test	58	3827	26	5
-test	80	4998	7	5
-test	169	4737	8	5
-test	306	3950	31	5
-test	348	4993	13	5
-test	393	4993	1	9
-test	398	4067	11	6
-test	434	4998	0	8
-test	478	4998	18	6
-test	556	4995	21	6
-test	590	4874	6	4
-test	856	4994	12	7
-test	859	4998	29	6
-test	987	3538	10	4
-test	1004	3993	19	5
-test	1057	4993	22	5
-test	1126	4614	5	5
-test	1142	4776	15	4
-test	1583	4510	20	5
-test	1634	4997	4	6
-test	1998	4632	3	6
-test	2685	4998	23	5
-test	2687	4788	24	4
+test	0	2393	13	5
+test	4	2453	31	5
+test	4	2607	294	5
+test	6	2186	21	4
+test	80	2149	48	4
+test	169	2624	7	4
+test	479	2585	35	4
+test	2685	4993	11	5
+test	2685	4994	36	6
+test	2685	4998	140	6
+test	2685	4998	16	6
+test	2685	4998	18	5
+test	2685	4998	46	6
+test	2686	4995	256	5
+test	2686	4997	26	5
+test	2686	4998	10	5
+test	2687	4993	6	5
+test	2689	4998	30	5
+test	2700	4776	38	4
+test	2717	4737	14	4
+test	2735	4874	45	5
+test	2786	4788	55	4
+test	2800	4829	51	4
+test	2912	4998	27	4