Crock

CHANGELOG.md

@@ -3,14 +3,21 @@
 The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/)
 and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
-## v1.0.0dev - [date]
+## v1.1.0dev - [date]
 
-Initial release of nf-core/microscopy, created with the [nf-core](https://nf-co.re/) template.
+First update to nfcore/microscopy.
 
 ### `Added`
 
+* Selection of segmentation method as a parameter.
+* Altered workflow to handle pre-stitched or legacy HALO tiff inputs.
+* Added a downscaling step which reduces resolution to 1px/um for faster segmentation.
+* Expanded README and added a metro diagram.
+
 ### `Fixed`
 
+* Fixed a bug where warnings were exported into tiff metadata xml.
+
 ### `Dependencies`
 
 ### `Deprecated`

README.md

@@ -7,47 +7,52 @@
 
 ## Introduction
 
-**nf-core/microscopy** is a bioinformatics pipeline that ...
+**nf-core/microscopy** is a bioinformatics pipeline designed to streamline the analysis of multiplex immunohistochemsitry (mIHC) samples.
 
-<!-- TODO nf-core:
-   Complete this sentence with a 2-3 sentence summary of what types of data the pipeline ingests, a brief overview of the
-   major pipeline sections and the types of output it produces. You're giving an overview to someone new
-   to nf-core here, in 15-20 seconds. For an example, see https://github.com/nf-core/rnaseq/blob/master/README.md#introduction
--->
+The purpose of the pipeline is to convert tiled or stitched multi-channel microscopy images into clean single-cell data. It takes as input a samplesheet, and an optional list of markers present on the panel. It stitches images together, performs cell segmentation on the DAPI channel, quantifies outputs into single-cell format and then runs a few basic clustering analyses.
 
-<!-- TODO nf-core: Include a figure that guides the user through the major workflow steps. Many nf-core
-     workflows use the "tube map" design for that. See https://nf-co.re/docs/guidelines/graphic_design/workflow_diagrams#examples for examples.   -->
-<!-- TODO nf-core: Fill in short bullet-pointed list of the default steps in the pipeline -->
+The segmentation model to be used can be selected (currently mesmer is default, and cellpose is also available). The pipeline returns a stitched image file, segmented image mask, cell x feature spreadsheet, as well as UMAP representations of the data clustered using various methods.
+
+![nf-core/microscopy metro diagram](assets/microscopy_metro.png)
 
 ## Usage
 
 > [!NOTE]
 > If you are new to Nextflow and nf-core, please refer to [this page](https://nf-co.re/docs/usage/installation) on how to set-up Nextflow. Make sure to [test your setup](https://nf-co.re/docs/usage/introduction#how-to-run-a-pipeline) with `-profile test` before running the workflow on actual data.
 
-<!-- TODO nf-core: Describe the minimum required steps to execute the pipeline, e.g. how to prepare samplesheets.
-     Explain what rows and columns represent. For instance (please edit as appropriate):
-
-First, prepare a samplesheet with your input data that looks as follows:
+First, prepare a samplesheet that looks as follows:
 
 `samplesheet.csv`:
 
 ```csv
-sample,fastq_1,fastq_2
-CONTROL_REP1,AEG588A1_S1_L002_R1_001.fastq.gz,AEG588A1_S1_L002_R2_001.fastq.gz
+sample,tiffs
+SAMPLE_NAME,/path/to/tiff/directory
 ```
 
-Each row represents a fastq file (single-end) or a pair of fastq files (paired end).
+Each row represents a directory which contains several .tiff tiles.
 
--->
+Next, you'll need to prepare a list of your markers:
 
-Now, you can run the pipeline using:
+`markers.csv`:
+
+```csv
+marker_name
+DAPI
+PTPRC
+HMWCK
+CD3
+CD11c
+```
 
-<!-- TODO nf-core: update the following command to include all required parameters for a minimal example -->
+Each row represents a different channel in your images. Note that marker names cannot contain spaces!
+
+Now, you can run the pipeline using:
 
 ```bash
 nextflow run nf-core/microscopy \
    -profile <docker/singularity/.../institute> \
    --input samplesheet.csv \
+   --markers markers.csv \
    --outdir <OUTDIR>
 ```
 
@@ -67,8 +72,8 @@ For more details about the output files and reports, please refer to the
 nf-core/microscopy was originally written by Song Li.
 
 We thank the following people for their extensive assistance in the development of this pipeline:
+  * Patrick Crock
 
-<!-- TODO nf-core: If applicable, make list of people who have also contributed -->
 
 ## Contributions and Support
 

assets/schema_input.json

@@ -17,8 +17,16 @@
                 "type": "string",
                 "format": "directory-path",
                 "exists": true,
-                "pattern": "^\\S+$",
-                "errorMessage": "TIFF folder path must be provided and cannot contain spaces"
+                "pattern": "^.+$",
+                "errorMessage": "TIFF folder path must be provided, and if it contains spaces, it should be quoted like so: 'example/path to/tiffs'"
+            },
+            "format": {
+                "type": "string",
+                "pattern": "^(tiles|stitched|fused)$",
+                "default": "tiles",
+                "errorMessage": "Format must be one of: tiles, stitched, fused",
+                "meta": ["format"],
+                "help_text": "Input format: 'tiles' for tiled TIFFs (default), 'stitched' for .ome.tiff files, 'fused' for HALO-fused .tiff files"
             }
         },
         "required": ["sample", "tiffs"]

bin/convert_ome_tiff.py

@@ -1,7 +1,7 @@
 #!/usr/bin/env python
 
-# Written by Song Li
-# Version: 0.0.1
+# Written by Song Li & Patrick Crock
+# Version: 0.0.2 Changed default order to match upstream outputs.
 
 import argparse
 import tifffile
@@ -11,8 +11,8 @@ def main():
     parser = argparse.ArgumentParser(description="Extract channel from image")
     parser.add_argument("-o", "--output", type=str, required=True, help="Output .tif file for extracted channel")
     parser.add_argument("-i", "--image", type=str, required=True, help="Input .tif image")
-    parser.add_argument("--order", type=str, default='2,0,1', help="Transpose dimensions, assuming Y,X,C input dimension order")
-    
+    parser.add_argument("--order", type=str, default='0,1,2', help="Transpose dimensions, assuming X,Y,C input dimension order")
+
     args = parser.parse_args()
 
     img = tifffile.imread(args.image)
@@ -29,4 +29,4 @@ def main():
 
 
 if __name__ == "__main__":
-    main()
+    main()

bin/extract_image_channel.py

@@ -1,32 +1,52 @@
 #!/usr/bin/env python
 
-# Written by Song Li
-# Version 0.0.1
-
+# Written by Song Li & Patrick Crock
+# Version 0.0.3 - Made extraction robust to dimension order
 import sys, os
 import argparse
-import xmltodict
-import skimage
+import xml.etree.ElementTree as ET
+import skimage.io
+import tifffile
+import numpy as np
 
 def extract_channel(xml, channel_name):
-    with open(xml, 'r') as file:
-        xml_str = file.read()
+    tree = ET.parse(xml)
+    root = tree.getroot()
 
-    xml_dict = xmltodict.parse(xml_str, attr_prefix='', cdata_key='')
-    channels = xml_dict['OME']['Image']['Pixels']['Channel'] 
+    # Handle namespace
+    ns = {'ome': 'http://www.openmicroscopy.org/Schemas/OME/2016-06'}
+    channels = root.findall('.//ome:Channel', ns)
 
-    for channel in channels:
-        if channel_name in channel['Name'].upper():
-            ch_extracted = channel['ID']
-            ch_extracted = ch_extracted.split(':')[-1]
+    channel_name_upper = channel_name.upper()
 
-            return int(ch_extracted) 
+    for i, channel in enumerate(channels):
+        name = channel.get('Name', '').upper()
+        if channel_name_upper in name:
+            return i  # Return the channel index directly
 
-    return 
+    return None
 
 def save_channel_image(ch, img_path, out_path):
-    img = skimage.io.imread(img_path)
-    channel = img[:,:,3]
+    with tifffile.TiffFile(img_path) as tif:
+        arr = tif.asarray()
+        axes = tif.series[0].axes
+        print("Array shape:", arr.shape, "Axes:", axes)
+
+    if "C" in axes:
+        c_index = axes.index("C")
+    elif "S" in axes:
+        c_index = axes.index("S")
+    else:
+        raise RuntimeError(f"No channel axis in axes string {axes}")
+
+    # move channel axis to last for consistency
+    arr = np.moveaxis(arr, c_index, -1)
+
+    if ch >= arr.shape[-1]:
+        print(f"Error: Channel {ch} not found, image has {arr.shape[-1]} channels")
+        sys.exit(os.EX_SOFTWARE)
+
+    channel = arr[..., ch]
     skimage.io.imsave(out_path, channel)
 
 def main():
@@ -35,19 +55,18 @@ def main():
     parser.add_argument("-i", "--image", type=str, required=True, help="Input .tif image")
     parser.add_argument("-x", "--xml", type=str, required=True, help="Metadata .xml for .tif image")
     parser.add_argument("-c", "--channel", type=str, default='DAPI', help="Channel name to extract")
-
-    args = parser.parse_args()
 
+    args = parser.parse_args()
     channel_extracted = extract_channel(args.xml, args.channel)
 
-    if channel_extracted:
+    if channel_extracted is not None:
         save_channel_image(channel_extracted, args.image, args.output)
+        print(f"Successfully extracted channel {channel_extracted} ({args.channel}) to {args.output}")
     else:
         print(f"{args.channel} channel could not be found")
         sys.exit(os.EX_SOFTWARE)
 
     return
 
-
 if __name__ == "__main__":
-    main()
+    main()

bin/indicaTIFF_to_ome.py

@@ -0,0 +1,128 @@
+#!/usr/bin/env python
+
+## from: https://forum.image.sc/t/trouble-generating-ome-tiffs-in-the-right-shape/89191/7
+
+import argparse
+from xml.etree import ElementTree
+from tifffile import TiffFile, TiffWriter
+from xml.dom import minidom
+
+def tiles(series):
+    # yield raw tiles from all pages in TIFF series
+    fh = series.parent.filehandle
+    for page in series:
+        for offset, bytecount in zip(page.dataoffsets, page.databytecounts):
+            fh.seek(offset)
+            yield fh.read(bytecount)
+
+
+def main():
+    parser = argparse.ArgumentParser(description="Convert HALO/Indica tif to OME tif")
+    parser.add_argument("-o", "--output", type=str, required=True, help="Output .ome.tif file")
+    parser.add_argument("-i", "--image", type=str, required=True, help="Input .tif image")
+
+    args = parser.parse_args()
+
+    with TiffFile(args.image) as tif:
+        print(tif)
+
+        tree = ElementTree.fromstring(tif.pages.first.description)
+        channel_names = [
+            channel.attrib['name'] for channel in tree.iter('channel')
+        ]
+
+        objective_value = float(tree.find('.//objective').attrib['value'])
+
+        with TiffWriter(
+            args.output, bigtiff=True, ome=True, byteorder=tif.byteorder
+        ) as ome:
+            for series in tif.series:
+                print("SERIES:", series)
+                print("  series.axes:", series.axes, "series.shape:", series.shape)
+                assert series.axes == 'IYX'
+                assert series.shape[0] == len(channel_names)
+
+                for i, level in enumerate(series.levels):
+                    page = level.keyframe
+                    print("  LEVEL", i, "level.shape:", level.shape, "page.keyframe.shape:", getattr(page, 'shape', None))
+                    print("    page.imagewidth:", getattr(page, 'imagewidth', None),
+                          "page.imagelength:", getattr(page, 'imagelength', None))
+                    print("    page.tile:", getattr(page, 'tile', None))
+                    print("    TileWidth/TileLength tags:",
+                          page.tags.get('TileWidth'), page.tags.get('TileLength'))
+
+                    assert page.is_tiled
+
+                    # --- compute rows/cols and samples confidently from page tags ---
+                    rows = getattr(page, 'imagelength', None)
+                    cols = getattr(page, 'imagewidth', None)
+                    # number of channels/samples (fall back to level.shape if unsure)
+                    samples = series.shape[0] if series.shape else (level.shape[0] if len(level.shape) == 3 else 1)
+
+                    # If level.shape looks like (rows, cols, samples) convert to (samples, rows, cols)
+                    if len(level.shape) == 3:
+                        a, b, c = level.shape
+                        if (a == rows and b == cols and c == samples):
+                            inferred_shape = (samples, rows, cols)
+                        elif (a == samples and b == rows and c == cols):
+                            inferred_shape = (a, b, c)  # already (samples, rows, cols)
+                        else:
+                            # fallback to explicit (samples, rows, cols)
+                            inferred_shape = (samples, rows, cols)
+                    else:
+                        inferred_shape = (rows, cols)
+
+                    # get tile dims from tags if present (TileWidth, TileLength)
+                    try:
+                        tile_w = page.tags['TileWidth'].value
+                        tile_l = page.tags['TileLength'].value
+                        tile = (tile_w, tile_l)
+                    except Exception:
+                        tile = page.tile or None
+
+                    if i == 0:
+                        # base-level metadata construction (unchanged)
+                        if len(series.levels) > 1:
+                            subifds = len(series.levels) - 1
+                        else:
+                            subifds = None
+                        resx, resy = page.get_resolution('micrometer')
+
+                        # Correct resolution by objective magnification
+                        resx = resx / objective_value
+                        resy = resy / objective_value
+
+                        metadata = {
+                            'axes': 'CYX',
+                            'PhysicalSizeX': resx,
+                            'PhysicalSizeXUnit': 'µm',
+                            'PhysicalSizeY': resy,
+                            'PhysicalSizeYUnit': 'µm',
+                            'Channel': {'Name': channel_names},
+                        }
+                    else:
+                        subifds = None
+                        metadata = None
+
+                    dtype = 'float32' if level.dtype == 'uint32' else level.dtype
+
+                    print("    -> writing with shape:", inferred_shape, "tile:", tile)
+
+                    ome.write(
+                        tiles(level),
+                        shape=inferred_shape,
+                        dtype=dtype,
+                        photometric=page.photometric,
+                        compression=page.compression,
+                        resolution=page.resolution,
+                        resolutionunit=page.resolutionunit,
+                        tile=tile,
+                        subifds=subifds,
+                        metadata=metadata,
+                    )
+
+    return
+
+
+if __name__ == "__main__":
+    main()

main.nf

@@ -64,7 +64,7 @@ workflow {
         params.input
     )
 
-    
+
     // WORKFLOW: Run main workflow
     //
     NFCORE_MICROSCOPY (