diff --git a/tcrdb/resources/external/scRNAseq/Seurat3.rmd b/tcrdb/resources/external/scRNAseq/Seurat3.rmd index a3911a6c7..9fb58d5aa 100644 --- a/tcrdb/resources/external/scRNAseq/Seurat3.rmd +++ b/tcrdb/resources/external/scRNAseq/Seurat3.rmd @@ -1,13 +1,16 @@ --- title: 'Seurat scRNA-seq Analysis' +output: html_document + --- ```{r Setup} -knitr::opts_chunk$set(message=FALSE, warning=FALSE,echo=TRUE,error = FALSE) library(knitr) library(OOSAP) +knitr::opts_chunk$set(message=FALSE, warning=FALSE, echo=TRUE, error = TRUE) + cores <- Sys.getenv('SEQUENCEANALYSIS_MAX_THREADS') if (cores != ''){ print(paste0('Setting future::plan to ', cores, ' cores')) diff --git a/tcrdb/resources/web/tcrdb/panel/LibraryExportPanel.js b/tcrdb/resources/web/tcrdb/panel/LibraryExportPanel.js index 4988c82e5..e3d3b9934 100644 --- a/tcrdb/resources/web/tcrdb/panel/LibraryExportPanel.js +++ b/tcrdb/resources/web/tcrdb/panel/LibraryExportPanel.js @@ -88,7 +88,7 @@ Ext4.define('TCRdb.panel.LibraryExportPanel', { border: false }, items: [{ - html: 'Add an ordered list of plates, using tab-delimited columns. The first column(s) are plate ID and library type (GEX, VDJ, CITE, or HTO). These can either be one column (i.e. G234-1, C234-1, H234-1, or T234-1), or as two columns (234-1 GEX or 234-1 HTO). An optional next column is the lane assignment (i.e. Novaseq1, HiSeq1, HiSeq2). Finally, an optional final column can be used to provide the alias for this pool. This is mostly used for CITE-Seq/HTOs, where multiple libraries are pre-pooled. See these examples:
' + + html: 'Add an ordered list of plates, using tab-delimited columns. The first column(s) are plate ID and library type (GEX, VDJ, CITE, or HTO). These can either be one column (i.e. G234-1, C234-1, H234-1, or T234-1), or as two columns (234-1 GEX or 234-1 HTO). An optional next column is the lane assignment (i.e. Novaseq1, HiSeq1, HiSeq2). Finally, an optional final column can be used to provide the alias for this pool. This is mostly used for CITE-Seq/HTOs, where multiple libraries are pre-pooled. Note, a wildcard can be used to specify all plates beginning with that prefix. See these examples:
' + '
' +
                                                 '234-2\tGEX
' + '234-2\tVDJ
' + @@ -101,6 +101,7 @@ Ext4.define('TCRdb.panel.LibraryExportPanel', { '235-2\tHTO\tHiSeq2\tBNB-HTO-1
' + 'H235-2\tHiSeq1\tBNB-HTO-1
' + 'C235-2\tHiSeq1\tBNB-HTO-1' + + 'C235-*\tHiSeq2\tBNB-HTO-2' + '
', border: false },{ @@ -171,6 +172,7 @@ Ext4.define('TCRdb.panel.LibraryExportPanel', { }, this); var hadError = false; + var wildcards = {}; Ext4.Array.forEach(text, function(r){ if (r.length < 2){ hadError = true; @@ -186,6 +188,12 @@ Ext4.define('TCRdb.panel.LibraryExportPanel', { Ext4.Array.forEach(r, function(val, idx){ r[idx] = Ext4.String.trim(val); }, this); + + if (r[0].match('\\*$')) { + var m = r[0].match('\\*$'); + var val = r[0].substr(0, m.index); + wildcards[val] = r; + } }, this); if (hadError) { @@ -193,7 +201,58 @@ Ext4.define('TCRdb.panel.LibraryExportPanel', { return; } - this.onSubmit(btn, text); + if (!Ext4.Object.isEmpty(wildcards)) { + LABKEY.Query.selectRows({ + method: 'POST', + containerPath: Laboratory.Utils.getQueryContainerPath(), + schemaName: 'tcrdb', + queryName: 'cdnas', + columns: 'rowid,plateId', + filterArray: [LABKEY.Filter.create('plateId', Ext4.Object.getKeys(wildcards).join(';'), LABKEY.Filter.Types.CONTAINS_ONE_OF)], + scope: this, + failure: LDK.Utils.getErrorCallback(), + success: function (results) { + if (results.rows.length) { + var prefixToPlate = {}; + Ext4.Array.forEach(results.rows, function (row) { + Ext4.Array.forEach(Ext4.Object.getKeys(wildcards), function (prefix) { + if (row.plateId && row.plateId.includes(prefix)) { + prefix = prefix + '*'; + prefixToPlate[prefix] = prefixToPlate[prefix] || []; + prefixToPlate[prefix].push(row.plateId); + } + }, this); + }, this); + + Ext4.Array.forEach(Ext4.Object.getKeys(prefixToPlate), function (prefix) { + prefixToPlate[prefix] = Ext4.unique(prefixToPlate[prefix]); + }, this); + + var updatedText = []; + var prefixes = Ext4.Object.getKeys(prefixToPlate); + Ext4.Array.forEach(text, function (r, idx) { + var plateId = r[0]; + if (prefixes.indexOf(plateId) == -1) { + updatedText.push(r); + } + else { + Ext4.Array.forEach(prefixToPlate[plateId], function(newPlate){ + var r2 = [].concat(r); + r2[0] = newPlate; + updatedText.push(r2); + }, this); + } + }, this); + + text = updatedText; + } + + this.onSubmit(btn, text); + } + }); + } else { + this.onSubmit(btn, text); + } } }] }); @@ -265,11 +324,11 @@ Ext4.define('TCRdb.panel.LibraryExportPanel', { var instrument = btn.up('tcrdb-libraryexportpanel').down('#instrument').getValue(); var plateId = btn.up('tcrdb-libraryexportpanel').down('#sourcePlates').getValue(); var delim = 'TAB'; - var extention = 'txt'; + var extension = 'txt'; var split = '\t'; if (instrument !== 'NextSeq (MPSSR)') { delim = 'COMMA'; - extention = 'csv'; + extension = 'csv'; split = ','; } @@ -277,7 +336,7 @@ Ext4.define('TCRdb.panel.LibraryExportPanel', { var rows = LDK.Utils.CSVToArray(Ext4.String.trim(val), split); LABKEY.Utils.convertToTable({ - fileName: plateId + '.' + extention, + fileName: plateId + '.' + extension, rows: rows, delim: delim }); @@ -757,10 +816,10 @@ Ext4.define('TCRdb.panel.LibraryExportPanel', { var delim = instrument === 'Novogene' ? '\t' : ','; Ext4.Array.forEach(sortedRows, function (r) { - processType(readsetIds, rows, r, 'readsetId', 'GEX', 500, 1, 'G', null, false); - processType(readsetIds, rows, r, 'enrichedReadsetId', 'TCR', 700, 1, 'T', null, false); - processType(readsetIds, rows, r, 'hashingReadsetId', 'HTO', 182, 5, 'H', 'Cell hashing, 190bp amplicon. Please QC individually and pool in equal amounts per lane', true); - processType(readsetIds, rows, r, 'citeseqReadsetId', 'CITE', 182, 5, 'C', 'CITE-Seq, 190bp amplicon. Please QC individually and pool in equal amounts per lane', false); + processType(readsetIds, rows, r, 'readsetId', 'GEX', 500, 0.01, 'G', null, false); + processType(readsetIds, rows, r, 'enrichedReadsetId', 'TCR', 700, 0.01, 'T', null, false); + processType(readsetIds, rows, r, 'hashingReadsetId', 'HTO', 182, 0.05, 'H', 'Cell hashing, 190bp amplicon. Please QC individually and pool in equal amounts per lane', true); + processType(readsetIds, rows, r, 'citeseqReadsetId', 'CITE', 182, 0.05, 'C', 'CITE-Seq, 190bp amplicon. Please QC individually and pool in equal amounts per lane', false); }, this); //add missing barcodes: diff --git a/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerCellHashingHandler.java b/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerCellHashingHandler.java index 00be83f8c..3cd2ecf90 100644 --- a/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerCellHashingHandler.java +++ b/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerCellHashingHandler.java @@ -63,7 +63,7 @@ public static List getDefaultHashingParams(boolean incl ToolParameterDescriptor.create("scanEditDistances", "Scan Edit Distances", "If checked, CITE-seq-count will be run using edit distances from 0-3 and the iteration with the highest singlets will be used.", "checkbox", new JSONObject(){{ put("checked", false); }}, false), - ToolParameterDescriptor.create("editDistance", "Edit Distance", null, "ldk-integerfield", null, 3), + ToolParameterDescriptor.create("editDistance", "Edit Distance", null, "ldk-integerfield", null, 2), ToolParameterDescriptor.create("minCountPerCell", "Min Reads/Cell", null, "ldk-integerfield", null, 5), ToolParameterDescriptor.create("useSeurat", "Use Seurat Calling", "If checked, the seurat HTO calling algorithm will be used.", "checkbox", null, true), ToolParameterDescriptor.create("useMultiSeq", "Use MultiSeq Calling", "If checked, the MultiSeq HTO calling algorithm will be used.", "checkbox", null, true) @@ -248,7 +248,7 @@ public static File processBarcodeFile(SequenceOutputHandler.JobContext ctx, File //prepare whitelist of cell indexes File cellBarcodeWhitelist = utils.getValidCellIndexFile(); Set uniqueBarcodes = new HashSet<>(); - ctx.getLogger().debug("writing cell barcodes"); + ctx.getLogger().debug("writing cell barcodes, using file: " + perCellTsv.getPath()); try (CSVWriter writer = new CSVWriter(PrintWriters.getPrintWriter(cellBarcodeWhitelist), ',', CSVWriter.NO_QUOTE_CHARACTER);CSVReader reader = new CSVReader(IOUtil.openFileForBufferedUtf8Reading(perCellTsv), '\t')) { int rowIdx = 0; diff --git a/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerSeuratHandler.java b/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerSeuratHandler.java index 7d1fd5930..b927c37f3 100644 --- a/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerSeuratHandler.java +++ b/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerSeuratHandler.java @@ -60,6 +60,7 @@ public class CellRangerSeuratHandler extends AbstractParameterizedOutputHandler< { private FileType _fileType = new FileType("cloupe", false); public static final String SEURAT_MAX_THREADS = "seuratMaxThreads"; + private static final String GTF_FILE_ID = "gtfFileId"; public CellRangerSeuratHandler() { @@ -103,11 +104,11 @@ private static List getDefaultParams() put("storeValues", "simple;cca"); }}, "simple"), ToolParameterDescriptor.create(SEURAT_MAX_THREADS, "Seurat Max Threads", "Because seurat can behave badly with multiple threads, this allows a separate cap to be used from the main job. This will allow CITE-Seq-Count and other tools to run with more threads.", "ldk-integerfield", null, 1), - ToolParameterDescriptor.createExpDataParam("gtfFile", "Gene File", "This is the ID of a GTF file containing genes from this genome.", "sequenceanalysis-genomefileselectorfield", new JSONObject() + ToolParameterDescriptor.createExpDataParam(GTF_FILE_ID, "Gene File", "This is the ID of a GTF file containing genes from this genome.", "sequenceanalysis-genomefileselectorfield", new JSONObject() {{ put("extensions", Arrays.asList("gtf")); put("width", 400); - put("allowBlank", false); + put("allowBlank", true); }}, null) )); @@ -165,8 +166,6 @@ public boolean doSplitJobs() public class Processor implements SequenceOutputProcessor { - private static final String GTF_FILE_ID = "gtfFileIf"; - @Override public void init(PipelineJob job, SequenceAnalysisJobSupport support, List inputFiles, JSONObject params, File outputDir, List actions, List outputsToCreate) throws UnsupportedOperationException, PipelineJobException { @@ -184,19 +183,17 @@ public void init(PipelineJob job, SequenceAnalysisJobSupport support, List gtfIds = new HashSet<>(); - for (SequenceOutputFile so : inputFiles) + if (params.get(GTF_FILE_ID) == null) { - ExpData gtf = null; - ExpRun run = ExperimentService.get().getExpRun(so.getRunId()); - if (run != null) + job.getLogger().info("attempting to infer GTF:"); + + //TODO: collapse by filepath + Set gtfIds = new HashSet<>(); + for (SequenceOutputFile so : inputFiles) { - List gtfDatas = run.getInputDatas("GTF File", null); - if (!gtfDatas.isEmpty()) - { - gtf = gtfDatas.get(0); - } - else + ExpData gtf = null; + ExpRun run = ExperimentService.get().getExpRun(so.getRunId()); + if (run != null) { //Because existing runs didnt explicitly track GTF as an input, try to infer: PipelineStatusFile sf = PipelineService.get().getStatusFile(run.getJobId()); @@ -226,23 +223,23 @@ public void init(PipelineJob job, SequenceAnalysisJobSupport support, List inputFiles, JobContext c RecordedAction action = new RecordedAction(getName()); ctx.addActions(action); - int gtfId = ctx.getSequenceSupport().getCachedObject(GTF_FILE_ID, Integer.class); + int gtfId = ctx.getParams().optInt(GTF_FILE_ID, -1); + if (gtfId == -1) + { + ctx.getLogger().debug("GTF file was not specified, defaulting to inferred file"); + gtfId = ctx.getSequenceSupport().getCachedObject(GTF_FILE_ID, Integer.class); + } + File gtfFile = ctx.getSequenceSupport().getCachedData(gtfId); if (!gtfFile.exists()) { @@ -447,7 +450,7 @@ public void processFilesRemote(List inputFiles, JobContext c writer.println(); writer.println("setwd('/work')"); - writer.println("rmarkdown::render('" + rmdScript.getName() + "', clean=TRUE, output_file='" + outHtml.getName() + "')"); + writer.println("rmarkdown::render('" + rmdScript.getName() + "', clean=TRUE, output_format = 'html_document', output_file='" + outHtml.getName() + "')"); } catch (IOException e) { @@ -689,11 +692,15 @@ else if (rs.getReadsetId() == null) throw new PipelineJobException(e); } - if (htosForReadset > 0) + if (htosForReadset > 1) { ctx.getLogger().info("Total HTOs for readset: " + htosForReadset); finalCalls.put(barcodePrefix, CellRangerCellHashingHandler.processBarcodeFile(ctx, barcodes, rs, htoReadset, so.getLibrary_id(), action, getClientCommandArgs(ctx.getParams()), false, SeuratCellHashingHandler.CATEGORY, true, perReadsetHtos, true)); } + else if (htosForReadset == 1) + { + ctx.getLogger().info("Only single HTO used for lane, skipping cell hashing calling"); + } else { ctx.getLogger().info("No HTOs found for readset"); diff --git a/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerVDJUtils.java b/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerVDJUtils.java index a71cf5bac..4c376e533 100644 --- a/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerVDJUtils.java +++ b/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerVDJUtils.java @@ -376,7 +376,7 @@ public File runRemoteVdjCellHashingTasks(PipelineStepOutput output, String outpu File cellBarcodeWhitelist = getValidCellIndexFile(); Set uniqueBarcodes = new HashSet<>(); Set uniqueBarcodesIncludingNoCDR3 = new HashSet<>(); - _log.debug("writing cell barcodes"); + _log.debug("writing cell barcodes, using file: " + perCellTsv.getPath()); try (CSVWriter writer = new CSVWriter(PrintWriters.getPrintWriter(cellBarcodeWhitelist), ',', CSVWriter.NO_QUOTE_CHARACTER); CSVReader reader = new CSVReader(Readers.getReader(perCellTsv), ',')) { int rowIdx = 0; diff --git a/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerVDJWrapper.java b/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerVDJWrapper.java index b63bdedf5..19fe214ca 100644 --- a/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerVDJWrapper.java +++ b/tcrdb/src/org/labkey/tcrdb/pipeline/CellRangerVDJWrapper.java @@ -68,6 +68,7 @@ public CellRangerVDJWrapper(@Nullable Logger logger) public static final String TARGET_ASSAY = "targetAssay"; public static final String DELETE_EXISTING_ASSAY_DATA = "deleteExistingAssayData"; + public static final String INNER_ENRICHMENT_PRIMERS = "innerEnrichmentPrimers"; public static class VDJProvider extends AbstractAlignmentStepProvider { @@ -82,6 +83,13 @@ public VDJProvider() ToolParameterDescriptor.createCommandLineParam(CommandLineParam.create("--force-cells"), "force-cells", "Force Cells", "Force pipeline to use this number of cells, bypassing the cell detection algorithm. Use this if the number of cells estimated by Cell Ranger is not consistent with the barcode rank plot.", "ldk-integerfield", new JSONObject(){{ put("minValue", 0); }}, null), + ToolParameterDescriptor.createCommandLineParam(CommandLineParam.createSwitch("--disable-ui"), "disable-ui", "Disable UI", "If checked, this will run cellranger with the optional web-based UI disabled.", "checkbox", new JSONObject(){{ + put("checked", true); + }}, true), + ToolParameterDescriptor.create(INNER_ENRICHMENT_PRIMERS, "Inner Enrichment Primers", "An option comma-separated list of the inner primers used for TCR enrichment. These will be used for trimming.", "textarea", new JSONObject(){{ + put("height", 100); + put("width", 400); + }}, null), ToolParameterDescriptor.create(TARGET_ASSAY, "Target Assay", "Results will be loaded into this assay. If no assay is selected, a table will be created with nothing in the DB.", "tcr-assayselectorfield", new JSONObject(){{ put("autoSelectAssay", false); }}, null), @@ -296,6 +304,32 @@ public AlignmentStep.AlignmentOutput performAlignment(Readset rs, File inputFast File indexDir = AlignerIndexUtil.getIndexDir(referenceGenome, getIndexCachedDirName(getPipelineCtx().getJob())); args.add("--reference=" + indexDir.getPath()); + String primers = StringUtils.trimToNull(getProvider().getParameterByName(INNER_ENRICHMENT_PRIMERS).extractValue(getPipelineCtx().getJob(), getProvider(), getStepIdx(), String.class, null)); + if (primers != null) + { + primers = primers.replaceAll("\\s+", ","); + primers = primers.replaceAll(",+", ","); + + File primerFile = new File(outputDirectory, "primers.txt"); + try (PrintWriter writer = PrintWriters.getPrintWriter(primerFile)) + { + Arrays.stream(primers.split(",")).forEach(x -> { + x = StringUtils.trimToNull(x); + if (x != null) + { + writer.println(x); + } + }); + } + catch (IOException e) + { + throw new PipelineJobException(e); + } + + output.addIntermediateFile(primerFile); + args.add("--inner-enrichment-primers=" + primerFile.getPath()); + } + args.addAll(getClientCommandArgs("=")); Integer maxThreads = SequencePipelineService.get().getMaxThreads(getPipelineCtx().getLogger()); diff --git a/tcrdb/src/org/labkey/tcrdb/pipeline/SeuratCiteSeqHandler.java b/tcrdb/src/org/labkey/tcrdb/pipeline/SeuratCiteSeqHandler.java index 853a68a0b..582b5dd92 100644 --- a/tcrdb/src/org/labkey/tcrdb/pipeline/SeuratCiteSeqHandler.java +++ b/tcrdb/src/org/labkey/tcrdb/pipeline/SeuratCiteSeqHandler.java @@ -27,7 +27,7 @@ public class SeuratCiteSeqHandler extends AbstractParameterizedOutputHandler