Xml 20260826

.gitignore

@@ -22,4 +22,5 @@ dockers
 smrtlink_dockers
 slim_smrtlink_dockers
 bin/fix_rg_ids.py
-bin/fix_rg_ids.sh
+bin/fix_rg_ids.sh
+docker

README.md

@@ -13,20 +13,18 @@ The pipeline starts with HiFi BAM files and has the following steps:
 2. The `LIST_HYBRIDS` and `REMOVE_HYBRIDS` processes identify and remove any reads with mismatched i7 and i5 seqWell barcode sequences in the remaining non-demultiplexed reads.
 3. The second Lima process, `LIMA_EITHER_END`, demultiplexes reads with only an i7 or i5 seqWell barcode sequence.
 4. The BAM files for each sample within each pool are merged in the `MERGE_READS` process and merged FASTQ files and bam files are created.
-5. The `CREATE_XML` process combines the indexed, merged BAMs for each pool into a single `ConsensusReadSet` XML file using PacBio SMRT Tools' `dataset create` CLI. This step only runs when the pipeline is executed with the default `main.nf` (docker); it is not available with `main_conda.nf`, since `dataset create` is only provided in the `seqwell/smrttools` docker image and is not installable via conda. See [Pipeline Entry Points](#pipeline-entry-points-mainnf-vs-main_condanf) below.
+5. The `CREATE_XML` process combines the indexed, merged BAMs for each pool into a single `ConsensusReadSet` XML file using the `pbcore` Python API. This step runs under both the `docker` and `conda` profiles.
 6. The `DEMUX_STATS` process generates a summary of the demultiplexing steps.
 7. If a `rename_map` is provided, the `RENAME_DEMUX_STATS` process renames the sample identifiers in the demultiplexing summary to match the user-defined sample names.
 8. `NANOSTAT` and `MULTIQC` are used to generate summary metrics for the reads assigned to each sample in the pool.
 9. `NANOSTAT_UNBARCODED` generates sequencing metrics for the unbarcoded reads remaining after both lima steps.
    Because the unbarcoded BAM is unaligned, reads are first converted to FASTQ via pysam before being passed to NanoStat.
 10. `DEMUX_QC` combines lima barcode statistics, per-sample NanoStat results, and unbarcoded NanoStat results into two final output tables per pool: a per-well stats table and a per-pool summary table.
 
-The final output from this pipeline includes Lima output files, demultiplexed BAM and FASTQ files, a demultiplexing summary, a MultiQC report collating NanoStat results, comprehensive per-pool and per-well demux statistics, and — when run with `main.nf` (docker) — a combined `ConsensusReadSet` XML file per pool.
+The final output from this pipeline includes Lima output files, demultiplexed BAM and FASTQ files, a demultiplexing summary, a MultiQC report collating NanoStat results, comprehensive per-pool and per-well demux statistics, and a combined `ConsensusReadSet` XML file per pool.
 
 ![Fig1. LongPlex Workflow](./docs/LongPlex_Workflow.png)
 
-*Dashed node: `CREATE_XML` only runs in the docker pipeline (`main.nf`); it is skipped when running `main_conda.nf`.*
-
 ## Dependencies
 
 This pipeline requires installation of [Nextflow](https://www.nextflow.io/docs/latest/install.html) version 22.04.5 through 25.10.4. This pipeline is not compatible with Nextflow 26.04.0 and later.
@@ -45,20 +43,21 @@ All docker containers used in this pipeline are publicly available.
 - *multiqc*: quay.io/biocontainers/multiqc:1.21--pyhdfd78af_0
 - *python*: python:3.12-bookworm
 - *pandas*: quay.io/biocontainers/pandas:1.5.2
-- *smrttools*: seqwell/smrttools:26.1.0-slim (docker-only; used by `CREATE_XML` to run PacBio SMRT Tools' `dataset create` CLI)
+- *pbcore*: seqwell/pbcore-2.6.0-2026:v1.0 (provides pbindex and pbcore for the CREATE_XML step)
 
 ## Conda Environment
 
-The conda environment is defined in `environment-pipeline.yml` and will be built automatically if the pipeline is run with `-profile conda`. Note that this profile is only supported on Linux systems, as **lima (v2.13.0)** is only available for Linux.
-
-## Pipeline Entry Points: `main.nf` vs `main_conda.nf`
+The conda environment is defined in `environment-pipeline.yml` and will be built automatically if the pipeline is run with `-profile conda`. This environment includes `pbtk` and `pbcore`, so the `CREATE_XML` step runs under conda as well as docker. Note that this profile is only supported on Linux systems, as **lima (v2.13.0)** is only available for Linux.
 
-This pipeline ships with two entry points:
+## Profiles:
 
-- **`main.nf`** (default): supports the `docker`, `apptainer`, and `singularity` profiles. This is the only entry point that includes the `CREATE_XML` step, since that step requires the `seqwell/smrttools` docker image.
-- **`main_conda.nf`**: used with `-profile conda`. This entry point is identical to `main.nf` except that it omits the `CREATE_XML` step, as SMRT Tools' `dataset create` CLI is not available via conda.
+Several profiles are available and can be selected with the `-profile` option at the command line.
 
-Be sure to point `nextflow run` at the entry point that matches your chosen profile, e.g. `-profile docker main.nf` or `-profile conda main_conda.nf`.
+- `apptainer`
+- `aws`
+- `conda`
+- `docker`
+- `singularity`
 
 # How to run the pipeline:
 
@@ -88,7 +87,6 @@ There are four required columns:
 The output directory path can be a local absolute path or an AWS S3 URI.
 If it is an AWS S3 URI, please make sure to [set your security credentials appropriately](https://www.nextflow.io/docs/latest/amazons3.html#security-credentials).
 
-
 ## Optional Parameters
 
 ### `rename_map`
@@ -121,16 +119,6 @@ When `rename_map` is provided:
 - When multiple pools are present in the `pool_sheet`, the `rename_map` may contain entries for all pools.
   Each pool's `DEMUX_QC` run will automatically filter the map to only its own entries using the `pool_ID` prefix, ensuring no cross-pool mixing.
 
-## Profiles:
-
-Several profiles are available and can be selected with the `-profile` option at the command line.
-
-- `apptainer`
-- `aws`
-- `conda` (use with `main_conda.nf`; skips the `CREATE_XML` step)
-- `docker` (use with `main.nf`)
-- `singularity`
-
 ## Example Command
 
 A minimal execution might look like:
@@ -181,7 +169,7 @@ nextflow run \
 ```bash
 nextflow run \
     -profile conda \
-    main_conda.nf \
+    main.nf \
     -c nextflow.config \
     --pool_sheet "${PWD}/tests/pool_sheet.csv" \
     --output "${PWD}/test_output" \
@@ -227,7 +215,7 @@ test_output/
 │   ├── multiqc/
 │   │   └── bc1015_multiqc_report.html                       # MultiQC report including NanoStat results
 │   └── xml/
-│       └── bc1015.combined.consensusreadset.xml             # Combined ConsensusReadSet XML for pool bc1015 (docker/main.nf only)
+│       └── bc1015.combined.consensusreadset.xml             # Combined ConsensusReadSet XML for pool bc1015
 └── logs/
     ├── execution_report_[DATE-TIME-STAMP].html              # Nextflow execution report
     ├── execution_timeline_[DATE-TIME-STAMP].html            # Nextflow execution timeline
@@ -268,3 +256,4 @@ One summary table per pool covering the full run.
 | Barcoded HiFi Read Length (mean, kb) | Weighted mean read length across all barcoded wells |
 | Unbarcoded HiFi Read Length (mean, kb) | Mean read length from `NANOSTAT_UNBARCODED` |
 
+

bin/create_combined_xml.py

@@ -0,0 +1,44 @@
+#!/usr/bin/env python3
+"""Build a combined ConsensusReadSet XML from a pool's demultiplexed BAMs.
+
+Equivalent to:
+    dataset create --generateIndices --force --type ConsensusReadSet \
+        <pool_id>.combined.consensusreadset.xml <bams...>
+"""
+import argparse
+import os
+import subprocess
+
+from pbcore.io import ConsensusReadSet
+
+
+def build_combined_dataset(inputs, pool_id, threads=2):
+    # inputs may be a mixed list of .bam and .bam.pbi paths (MERGE_READS
+    # emits both together as bam_w_index) -- separate them out.
+    bams = [f for f in inputs if f.endswith(".bam")]
+
+    for bam in bams:
+        if not os.path.exists(bam + ".pbi"):
+            subprocess.run(["pbindex", "-j", str(threads), bam], check=True)
+
+    ds = ConsensusReadSet(*bams, strict=True)
+    ds.name = f"{pool_id} (combined)"
+
+    xml_path = f"{pool_id}.combined.consensusreadset.xml"
+    ds.write(xml_path)
+    return xml_path
+
+
+def main():
+    parser = argparse.ArgumentParser(description=__doc__)
+    parser.add_argument("pool_id")
+    parser.add_argument("inputs", nargs="+", help="bam and/or bam.pbi files")
+    parser.add_argument("-j", "--threads", type=int, default=2)
+    args = parser.parse_args()
+
+    xml_path = build_combined_dataset(args.inputs, args.pool_id, args.threads)
+    print(f"Combined dataset: {xml_path}")
+
+
+if __name__ == "__main__":
+    main()

docs/flow.txt

@@ -32,13 +32,13 @@ digraph workflow {
     MERGE_READS -> NANOSTAT;
     NANOSTAT -> MULTIQC;
 
-    // Docker-only path (main.nf); skipped by main_conda.nf
-    CREATE_XML [style="rounded,filled,dashed", fillcolor="#cfe2f3", color="#6c8ebf"];
-    MERGE_READS -> CREATE_XML [style=dashed];
+    // Runs under both docker and conda profiles
+    CREATE_XML [style="rounded,filled", fillcolor="#cfe2f3", color="#6c8ebf"];
+    MERGE_READS -> CREATE_XML;
 
     // Compact column grouping
     { rank=same; LIST_HYBRIDS; REMOVE_HYBRIDS }
     { rank=same; DEMUX_STATS; LIMA_EITHER_END }
     { rank=same; NANOSTAT; MERGE_READS }
     { rank=same; MULTIQC; CREATE_XML }
-}
+}

environment-pipeline.yml

@@ -8,6 +8,7 @@ dependencies:
   - bioconda::multiqc=1.21
   - bioconda::picard=3.2.0
   - bioconda::samtools=1.19.2
+  - bioconda::pbtk=3.5.0
   - conda-forge::python=3.11.0
   - conda-forge::pip=26.1.2
   - conda-forge::r-base=4.3.3
@@ -17,6 +18,10 @@ dependencies:
   - conda-forge::r-stringr=1.5.1
   - conda-forge::r-tidyr=1.3.1
   - conda-forge::defopt=6.4.0
+  - conda-forge::setuptools<82
+  - conda-forge::git
+
 
   - pip:
     - git+https://github.com/seqwell/longplexpy@0.2.2
+    - git+https://github.com/PacificBiosciences/pbcore

modules/create_xml.nf

@@ -1,21 +1,15 @@
 process CREATE_XML {
     tag "${meta.pool_ID}"
-
+    
     input:
     tuple val(meta), path(bams)
-
     output:
     tuple val(meta), path("*.combined.consensusreadset.xml"), emit: xml
-
     script:
     """
-    dataset create \\
-        --generateIndices --force --type ConsensusReadSet \\
-        ${meta.pool_ID}.combined.consensusreadset.xml \\
-        ${bams}
-
-    sed -i -E 's#(ResourceId=")[^"]*/([^/"]+)"#\\1${params.output}/${meta.pool_ID}/\\2"#g' \
-    ${meta.pool_ID}.combined.consensusreadset.xml
+    create_combined_xml.py ${meta.pool_ID} ${bams}
 
+    sed -i -E 's#(ResourceId=")[^"]*/([^/"]+)"#\\1${params.output}/${meta.pool_ID}/\\2"#g' \\
+    ${meta.pool_ID}.combined.consensusreadset.xml
     """
 }

nextflow.config

@@ -155,7 +155,7 @@ process {
     withName:  CREATE_XML { 
         memory = '8.0G'
         cpus = 2
-        container = 'seqwell/smrttools:26.1.0-slim' 
+        container = 'seqwell/pbcore-2.6.0-2026:v1.0' 
         }
       publishDir = [
             path: { "${params.output}/${meta.pool_ID}/xml" },

tests/run_longplex_integration.sh

@@ -4,6 +4,6 @@ nextflow run \
     -profile docker \
     -c "${PWD}/nextflow.config" \
     -c "${PWD}/tests/limit-resources.config" \
-    main_conda.nf \
+    main.nf \
     --pool_sheet "${PWD}/tests/pool_sheet.csv" \
     --output  "${PWD}/output/"