Add support to calculate per-cell saturation

Author: bbimber

singlecell/api-src/org/labkey/api/singlecell/pipeline/AbstractSingleCellPipelineStep.java

@@ -148,7 +148,7 @@ protected File createRmd(SingleCellOutput output, SequenceOutputHandler.JobConte
             markdown.setup = new SetupChunk(getRLibraries());
             markdown.chunks = new ArrayList<>();
             markdown.chunks.add(createParamChunk(inputObjects, outputPrefix));
-            markdown.chunks.addAll(getChunks());
+            markdown.chunks.addAll(getChunks(ctx));
             markdown.chunks.add(createFinalChunk());
 
             markdown.print(out);
@@ -161,7 +161,7 @@ protected File createRmd(SingleCellOutput output, SequenceOutputHandler.JobConte
         return outfile;
     }
 
-    protected List<Chunk> getChunks() throws PipelineJobException
+    protected List<Chunk> getChunks(SequenceOutputHandler.JobContext ctx) throws PipelineJobException
     {
         List<Chunk> ret = new ArrayList<>();
         ret.add(new Chunk(getProvider().getName(), getProvider().getLabel(), null, loadChunkFromFile()));

singlecell/resources/chunks/AppendSaturation.R

@@ -0,0 +1,16 @@
+for (datasetId in names(seuratObjects)) {
+    seuratObj <- seuratObjects[[datasetId]]
+    seuratObjects[[datasetId]] <- NULL
+    molInfoFile <- molInfoFiles[[datasetId]]
+    if (is.null(molInfoFile)) {
+        stop(paste0('Unable to find molInfo file for: ', datasetId))
+    }
+
+    seuratObj <- CellMembrane::AppendPerCellSaturation(seuratObj, molInfoFile = molInfoFile)
+
+    newSeuratObjects[[datasetId]] <- seuratObj
+
+    # Cleanup
+    rm(seuratObj)
+    gc()
+}

singlecell/src/org/labkey/singlecell/SingleCellModule.java

@@ -160,6 +160,7 @@ public static void registerPipelineSteps()
         SequencePipelineService.get().registerPipelineStep(new CiteSeqPlots.Provider());
         SequencePipelineService.get().registerPipelineStep(new PhenotypePlots.Provider());
         SequencePipelineService.get().registerPipelineStep(new AppendMetadata.Provider());
+        SequencePipelineService.get().registerPipelineStep(new AppendSaturation.Provider());
     }
 
     @Override

singlecell/src/org/labkey/singlecell/analysis/ProcessSingleCellHandler.java

@@ -12,7 +12,7 @@
 
 public class ProcessSingleCellHandler extends AbstractSingleCellHandler
 {
-    private static FileType LOUPE_TYPE = new FileType("cloupe", false);
+    public static FileType LOUPE_TYPE = new FileType("cloupe", false);
 
     @Override
     public String getName()

singlecell/src/org/labkey/singlecell/pipeline/singlecell/AbstractCellHashingCiteseqStep.java

@@ -40,7 +40,7 @@ protected File createRmd(SingleCellOutput output, SequenceOutputHandler.JobConte
             markdown.chunks.add(createDataChunk(countData, ctx.getOutputDir()));
             markdown.chunks.addAll(addAdditionalChunks(ctx, inputObjects, countData));
 
-            markdown.chunks.addAll(getChunks());
+            markdown.chunks.addAll(getChunks(ctx));
             markdown.chunks.add(createFinalChunk());
 
             markdown.print(out);

singlecell/src/org/labkey/singlecell/pipeline/singlecell/AppendSaturation.java

@@ -0,0 +1,130 @@
+package org.labkey.singlecell.pipeline.singlecell;
+
+import au.com.bytecode.opencsv.CSVReader;
+import au.com.bytecode.opencsv.CSVWriter;
+import org.apache.commons.io.FileUtils;
+import org.labkey.api.pipeline.PipelineJobException;
+import org.labkey.api.reader.Readers;
+import org.labkey.api.sequenceanalysis.SequenceOutputFile;
+import org.labkey.api.sequenceanalysis.pipeline.AbstractPipelineStepProvider;
+import org.labkey.api.sequenceanalysis.pipeline.PipelineContext;
+import org.labkey.api.sequenceanalysis.pipeline.SequenceOutputHandler;
+import org.labkey.api.singlecell.pipeline.SingleCellStep;
+import org.labkey.api.writer.PrintWriters;
+
+import java.io.File;
+import java.io.IOException;
+import java.util.ArrayList;
+import java.util.Arrays;
+import java.util.List;
+
+import static org.labkey.singlecell.analysis.ProcessSingleCellHandler.LOUPE_TYPE;
+
+public class AppendSaturation extends AbstractCellMembraneStep
+{
+    public AppendSaturation(PipelineContext ctx, AppendSaturation.Provider provider)
+    {
+        super(provider, ctx);
+    }
+
+    public static class Provider extends AbstractPipelineStepProvider<SingleCellStep>
+    {
+        public Provider()
+        {
+            super("AppendSaturation", "Append Saturation", "CellMembrane", "This will calculate sequencing saturation per-cell and append as metadata.", Arrays.asList(
+
+            ), null, null);
+        }
+
+        @Override
+        public AppendSaturation create(PipelineContext ctx)
+        {
+            return new AppendSaturation(ctx, this);
+        }
+    }
+
+    @Override
+    public void init(SequenceOutputHandler.JobContext ctx, List<SequenceOutputFile> inputFiles) throws PipelineJobException
+    {
+        try (CSVWriter writer = new CSVWriter(PrintWriters.getPrintWriter(getMolInfoTable(ctx))))
+        {
+            for (SequenceOutputFile so : inputFiles)
+            {
+                if (!LOUPE_TYPE.isType(so.getFile()))
+                {
+                    throw new PipelineJobException("All input files must be loupe files to use sequence saturation");
+                }
+
+                File molInfo = new File(so.getFile().getParentFile(), "molecule_info.h5");
+                if (!molInfo.exists())
+                {
+                    throw new PipelineJobException("Cannot find file: " + molInfo.getPath());
+                }
+
+                writer.writeNext(new String[]{String.valueOf(so.getRowid()), molInfo.getPath()});
+            }
+        }
+        catch (IOException e)
+        {
+            throw new PipelineJobException(e);
+        }
+    }
+
+    private File getMolInfoTable(SequenceOutputHandler.JobContext ctx)
+    {
+        return new File(ctx.getSourceDirectory(), "molInfo.txt");
+    }
+
+    @Override
+    protected List<Chunk> getChunks(SequenceOutputHandler.JobContext ctx) throws PipelineJobException
+    {
+        List<String> lines = new ArrayList<>();
+
+        lines.add("molInfoFiles <- list(");
+        File table = getMolInfoTable(ctx);
+        if (!table.exists())
+        {
+            throw new PipelineJobException("Unable to find molInfo table: " + table.getPath());
+        }
+
+        try (CSVReader reader = new CSVReader(Readers.getReader(table)))
+        {
+            String[] line;
+            while ((line = reader.readNext()) != null)
+            {
+                File source = new File(line[1]);
+                File dest = new File(ctx.getWorkingDirectory(), line[0] + ".molInfo.h5");
+                if (dest.exists())
+                {
+                    dest.delete();
+                }
+                FileUtils.copyFile(source, dest);
+
+                lines.add("\t'" + line[0] + "' = '" + dest.getName() + "',");
+                ctx.getFileManager().addIntermediateFile(dest);
+            }
+        }
+        catch (IOException e)
+        {
+            throw new PipelineJobException(e);
+        }
+
+        // Remove trailing comma:
+        int lastIdx = lines.size() - 1;
+        lines.set(lastIdx, lines.get(lastIdx).replaceAll(",$", ""));
+
+        lines.add(")");
+        lines.add("");
+
+        List<Chunk> ret = super.getChunks(ctx);
+        ret.add(new Chunk("molInfoFiles", null, null, lines, null));
+
+        return ret;
+    }
+
+    @Override
+    public String getFileSuffix()
+    {
+        return "saturation";
+    }
+}