Enable CITE-seq processing

Author: bbimber

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

@@ -415,11 +415,11 @@ protected Chunk createFinalChunk() throws PipelineJobException
 
     public static class Chunk
     {
-        String header;
-        String extraText;
-        String chunkName;
-        String chunkOpts;
-        List<String> bodyLines;
+        public String header;
+        public String extraText;
+        public String chunkName;
+        public String chunkOpts;
+        public List<String> bodyLines;
 
         public Chunk(String chunkName, @Nullable String header, @Nullable String extraText, List<String> bodyLines)
         {

singlecell/resources/chunks/AppendCiteSeq.R

@@ -5,10 +5,20 @@ for (datasetId in names(seuratObjects)) {
         stop(paste0('No CITE-seq information found for datasetId: ', datasetId))
     }
 
-    # callFile <- featureData[[datasetId]]
-    # if (!is.null(callFile)) {
-    #     seuratObj <- CellMembrane::AppendCiteSeq(seuratObj, barcodeCallFile = callFile, barcodePrefix = datasetId)
-    # }
+    adtWhitelist <- NULL
+    featureMetadata <- NULL
+    if (datasetId %in% names(featureMetadataFiles) && !is.null(featureMetadataFiles[[datasetId]])) {
+        featureMetadata <- read.table(featureMetadataFiles[[datasetId]], sep = '\t', header = T)
+        featureMetadata$rowname <- paste0(featureMetadata$tagname, '-', featureMetadata$sequence)
+        adtWhitelist <- featureMetadata$rowname
+    }
+
+    matrixDir <- featureData[[datasetId]]
+    if (!is.null(matrixDir)) {
+        seuratObj <- CellMembrane::AppendCiteSeq(seuratObj, unfilteredMatrixDir = matrixDir, normalizeMethod = normalizeMethod, datasetId = datasetId, featureMetadata = featureMetadata, adtWhitelist = adtWhitelist)
+    } else {
+        print('matrixDir was NULL, skipping CITE-seq')
+    }
 
     newSeuratObjects[[datasetId]] <- seuratObj
 

singlecell/src/org/labkey/singlecell/CellHashingServiceImpl.java

@@ -1360,6 +1360,11 @@ public Map<String, Object> executeCiteSeqCount(File outputDir, String basename,
                 baseArgs.add("-cells");
                 baseArgs.add(parameters.cells == null ? "0" : String.valueOf(parameters.cells));
             }
+            else if (parameters.cells != null)
+            {
+                baseArgs.add("-cells");
+                baseArgs.add(String.valueOf(parameters.cells));
+            }
 
             Integer cores = SequencePipelineService.get().getMaxThreads(log);
             if (cores != null)
@@ -1370,7 +1375,7 @@ public Map<String, Object> executeCiteSeqCount(File outputDir, String basename,
 
             for (ToolParameterDescriptor param : CellHashingHandler.getDefaultParams(parameters.type))
             {
-                if (parameters.cellBarcodeWhitelistFile != null && param.getName().equals("cells"))
+                if ((parameters.cellBarcodeWhitelistFile != null || parameters.cells != null) && param.getName().equals("cells"))
                 {
                     continue;
                 }

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

@@ -1,5 +1,6 @@
 package org.labkey.singlecell.pipeline.singlecell;
 
+import org.apache.commons.io.FileUtils;
 import org.json.JSONObject;
 import org.labkey.api.pipeline.PipelineJobException;
 import org.labkey.api.sequenceanalysis.SequenceOutputFile;
@@ -13,11 +14,13 @@
 import org.labkey.api.singlecell.pipeline.SingleCellOutput;
 import org.labkey.api.singlecell.pipeline.SingleCellStep;
 import org.labkey.api.util.PageFlowUtil;
+import org.labkey.singlecell.CellHashingServiceImpl;
 import org.labkey.singlecell.analysis.CellRangerSeuratHandler;
 import org.labkey.singlecell.analysis.SeuratCellHashingHandler;
 import org.labkey.singlecell.analysis.SeuratCiteSeqHandler;
 
 import java.io.File;
+import java.io.IOException;
 import java.util.ArrayList;
 import java.util.Collection;
 import java.util.Collections;
@@ -106,22 +109,75 @@ protected Map<Integer, File> prepareCountData(SingleCellOutput output, SequenceO
                 params.outputCategory = SeuratCiteSeqHandler.CATEGORY;
                 params.createOutputFiles = true;
                 params.genomeId = wrapper.getSequenceOutputFile().getLibrary_id();
-                params.cellBarcodeWhitelistFile = cellBarcodesParsed;
+                //params.cellBarcodeWhitelistFile = cellBarcodesParsed;
                 params.cells = 250000;
 
                 finalOutput = CellHashingService.get().processCellHashingOrCiteSeqForParent(parentReadset, output, ctx, params);
+
+                File validAdt = CellHashingServiceImpl.get().getValidCiteSeqBarcodeMetadataFile(ctx.getSourceDirectory(), parentReadset.getReadsetId());
+                if (!validAdt.exists())
+                {
+                    throw new PipelineJobException("Unable to find ADT metadata. expected: " + validAdt.getPath());
+                }
+
+                try
+                {
+                    FileUtils.copyFile(validAdt, getAdtMetadata(finalOutput));
+                }
+                catch (IOException e)
+                {
+                    throw new PipelineJobException(e);
+                }
             }
             else
             {
                 ctx.getLogger().info("CITE-seq not used, skipping: " + parentReadset.getName());
             }
 
-            dataIdToCalls.put(wrapper.getSequenceOutputFileId(), finalOutput);
+            dataIdToCalls.put(wrapper.getSequenceOutputFileId(), finalOutput.getParentFile());
         }
 
         return dataIdToCalls;
     }
 
+    public File getAdtMetadata(File countMatrix)
+    {
+        return new File(countMatrix.getParentFile(), "adtMetadata.txt");
+    }
+
+    @Override
+    protected Chunk createDataChunk(Map<Integer, File> hashingData)
+    {
+        Chunk ret = super.createDataChunk(hashingData);
+
+        List<String> lines = new ArrayList<>();
+
+        lines.add("featureMetadataFiles <- list(");
+        for (Integer key : hashingData.keySet())
+        {
+            if (hashingData.get(key) == null)
+            {
+                lines.add("\t'" + key + "' = NULL,");
+            }
+            else
+            {
+                File meta = getAdtMetadata(hashingData.get(key));
+                lines.add("\t'" + key + "' = '" + meta.getName() + "',");
+            }
+        }
+
+        // Remove trailing comma:
+        int lastIdx = lines.size() - 1;
+        lines.set(lastIdx, lines.get(lastIdx).replaceAll(",$", ""));
+
+        lines.add(")");
+        lines.add("");
+
+        ret.bodyLines.addAll(lines);
+
+        return ret;
+    }
+
     @Override
     public boolean isIncluded(SequenceOutputHandler.JobContext ctx, List<SequenceOutputFile> inputs) throws PipelineJobException
     {