Pre-sort BAM by tag for velocyto

Author: bbimber

SequenceAnalysis/api-src/org/labkey/api/sequenceanalysis/pipeline/SamSorter.java

@@ -4,7 +4,6 @@
 import org.apache.commons.io.FileUtils;
 import org.apache.commons.lang3.StringUtils;
 import org.apache.logging.log4j.Logger;
-import org.apache.logging.log4j.LogManager;
 import org.jetbrains.annotations.Nullable;
 import org.labkey.api.pipeline.PipelineJobException;
 import org.labkey.api.pipeline.PipelineJobService;
@@ -30,6 +29,11 @@ public SamSorter(@Nullable Logger logger)
     }
 
     public File execute(File input, @Nullable File output, SAMFileHeader.SortOrder sortOrder) throws PipelineJobException
+    {
+        return execute(input, output, sortOrder, null);
+    }
+
+    public File execute(File input, @Nullable File output, SAMFileHeader.SortOrder sortOrder, @Nullable List<String> extraArgs) throws PipelineJobException
     {
         getLogger().info("Sorting SAM/BAM: " + input.getPath());
 
@@ -50,7 +54,7 @@ public File execute(File input, @Nullable File output, SAMFileHeader.SortOrder s
         }
         else if (sortOrder != SAMFileHeader.SortOrder.coordinate)
         {
-            throw new PipelineJobException("Improper sort order: " + sortOrder.name());
+            getLogger().debug("Using sort order: " + sortOrder.name());
         }
 
         Integer threads = SequencePipelineService.get().getMaxThreads(getLogger());
@@ -71,6 +75,11 @@ else if (sortOrder != SAMFileHeader.SortOrder.coordinate)
         params.add("-T");
         params.add(tmpDir);
 
+        if (extraArgs != null)
+        {
+            params.addAll(extraArgs);
+        }
+
         params.add("-o");
         params.add(output.getPath());
 

singlecell/src/org/labkey/singlecell/run/VelocytoAlignmentStep.java

@@ -2,6 +2,7 @@
 
 import au.com.bytecode.opencsv.CSVReader;
 import au.com.bytecode.opencsv.CSVWriter;
+import htsjdk.samtools.SAMFileHeader;
 import org.apache.logging.log4j.Logger;
 import org.jetbrains.annotations.Nullable;
 import org.json.JSONObject;
@@ -14,6 +15,8 @@
 import org.labkey.api.sequenceanalysis.pipeline.AlignmentStepProvider;
 import org.labkey.api.sequenceanalysis.pipeline.PipelineContext;
 import org.labkey.api.sequenceanalysis.pipeline.ReferenceGenome;
+import org.labkey.api.sequenceanalysis.pipeline.SamSorter;
+import org.labkey.api.sequenceanalysis.pipeline.SamtoolsRunner;
 import org.labkey.api.sequenceanalysis.pipeline.SequencePipelineService;
 import org.labkey.api.sequenceanalysis.pipeline.ToolParameterDescriptor;
 import org.labkey.api.sequenceanalysis.run.AbstractCommandWrapper;
@@ -149,6 +152,13 @@ public File runVelocytoFor10x(File localBam, File gtf, File outputFolder, @Nulla
 
             getLogger().debug("Using BAM: " + localBam.getPath());
 
+            //Pre-sort since velocyto dies when doing this itself:
+            //samtools ``sort -l [compression] -m [mb_to_use]M -t [tagname] -O BAM -@ [threads_to_use] -o cellsorted_[bamfile] [bamfile]
+            List<String> sortArgs = Arrays.asList("-t", "CB");
+            File sortedBam = new File(localBam.getParentFile(), "cellsorted_" + localBam.getName());
+            new SamSorter(getLogger()).execute(localBam, sortedBam, SAMFileHeader.SortOrder.unknown, sortArgs);
+
+
             SimpleScriptWrapper wrapper = new SimpleScriptWrapper(getLogger());
             List<String> args = new ArrayList<>();
             args.add(SequencePipelineService.get().getExeForPackage("VELOCYTOPATH", "velocyto").getPath());
@@ -205,6 +215,8 @@ public File runVelocytoFor10x(File localBam, File gtf, File outputFolder, @Nulla
                 throw new PipelineJobException("Missing expected file: " + loom.getPath());
             }
 
+            sortedBam.delete();
+
             return loom;
         }
     }