diff --git a/Dockerfile b/Dockerfile
index 91574fc..db580d2 100644
--- a/Dockerfile
+++ b/Dockerfile
@@ -1,13 +1,23 @@
-FROM continuumio/miniconda3
-
-COPY environment.yml .
-RUN apt-get update -qq && apt-get install -y \
- build-essential \
- ffmpeg \
- libsm6 \
- libxext6
-
-RUN conda env create -f environment.yml
-ENV PATH="/opt/conda/envs/backsub/bin:$PATH"
-WORKDIR /background_subtraction
-COPY . .
\ No newline at end of file
+FROM mambaorg/micromamba:1.5.10-noble
+
+# Copy conda environment file
+COPY --chown=$MAMBA_USER:$MAMBA_USER ./environment.yml /tmp/conda.yml
+
+# Install environment
+RUN micromamba install -y -n base -f /tmp/conda.yml \
+ && micromamba install -y -n base conda-forge::procps-ng \
+ && micromamba env export --name base --explicit > environment.lock \
+ && echo ">> CONDA_LOCK_START" \
+ && cat environment.lock \
+ && echo "<< CONDA_LOCK_END" \
+ && micromamba clean -a -y
+
+# Switch to root to copy everything
+USER root
+
+# Ensure micromamba binaries are in PATH
+ENV PATH="$MAMBA_ROOT_PREFIX/bin:$PATH"
+
+# Copy the rest of the current directory into /app inside the container
+WORKDIR /app
+COPY ./backsub .
\ No newline at end of file
diff --git a/README.md b/README.md
index 03b7f18..f970494 100644
--- a/README.md
+++ b/README.md
@@ -13,7 +13,17 @@ Marker*corrected* = Marker*raw* - Background / Exposure1:
+ print(f"""Warning: Background channel with name {markers.background[channel]}
+ appears several times in the column "marker_name".
+ Only the first occurrence will be used for the subtraction.""" )
+ bg_idx=bg_idx[0]
+ scaling_factor[channel] = markers.exposure[channel] / markers.exposure[bg_idx]
+ background_idx[channel] = bg_idx
+
+ markers.insert(markers.shape[1], "factor", scaling_factor)
+ markers.insert(markers.shape[1], "bg_idx", background_idx)
+
+ return markers
+
+
+def extract_img_props(img_path,pixel_size=None):
+
+ #Extract data_type, pyramidal specs, height, width
+ with tifff.TiffFile(img_path) as tif:
+ pyr_levels=len(tif.series[0].levels)
+ is_pyramid=pyr_levels > 1
+ data_type=tif.series[0].dtype.name
+ height,width=tif.series[0].shape[-2::]
+ dask_chunksize=da.from_array(tif.pages[0].asarray(), chunks='auto').chunksize
+ if is_pyramid:
+ #dimensions of the reduced resolution layers(subresolution_dimensions)
+ subres_dims=[tif.series[0].levels[lvl].pages[0].shape
+ for lvl in range(1,pyr_levels)
+ ]
+ else:
+ subres_dims=None
+
+ #Try to extract pixel size from ome-xml
+ if pixel_size is None:
+ print('Pixel size not specified in the arguments (-mpp)')
+ try:
+ metadata = ome_types.from_tiff(img_path)
+ pixel_size = metadata.images[0].pixels.physical_size_x
+ pixel_size_unit = metadata.images[0].pixels.physical_size_x_unit
+ except Exception as err:
+ print(err)
+ print('Pixel size or pixel size unit detection using ome-types failed')
+ pixel_size = 1
+ pixel_size_unit="pixel"
+ else:
+ pixel_size=pixel_size
+ pixel_size_unit="µm"
+
+
+ img_props={"pixel_size":pixel_size,
+ "pixel_size_unit":pixel_size_unit,
+ "data_type":data_type,
+ "pyramid":is_pyramid,
+ "levels":pyr_levels,
+ "sub_levels_dims":subres_dims,
+ "size_x":width,
+ "size_y":height ,
+ "chunksize":dask_chunksize
+ }
+
+ return img_props
+
+
+def subtract_channels(src_img_path,
+ signal_index,
+ background_index,
+ factor,
+ ref_chunksize,
+ ref_dtype,
+ task_no
+ ):
+ """
+ This function executes the background substraction using generators, each element of the generator
+ is a tuple with 3 values, such that:
+
+ tuple=(img_with_backsub[array],calculate or extract pyramid [str],pyramid_from_index[int])
+
+ The second entry of the tuple indicates in the writing process if the pyramid should be calculated using
+ pyramid_gaussian from scikit image. If extract, the index given in the third entry will fetch all the pyramid
+ levels from the original image stack(src_img_path).
+ """
+ factor=np.float32(factor)#limiting precision to float32 saves memory
+ signal_as_zarr = zarr.open( tifff.imread( src_img_path, aszarr=True, series=0, level=0,key=int(signal_index) ) )
+ background_as_zarr =zarr.open( tifff.imread( src_img_path, aszarr=True, series=0, level=0,key=int(background_index) ) )
+ signal=da.from_zarr(signal_as_zarr, chunks=ref_chunksize )
+ background=da.from_zarr(background_as_zarr, chunks=ref_chunksize)
+ subtraction=da.clip(signal-factor*background,0,65535).astype(ref_dtype)
+ with ResourceProfiler(dt=0.25) as resources:
+ with ProgressBar():
+ result=subtraction.compute()
+ print(f"Resources used by dask during subtraction {task_no}:")
+ print(resources.results[0],"([sec],[MB],[% CPU usage])")
+ return result
+
+def extract_sublevels_from_tiff(path,ch,levs):
+ with tifff.TiffFile(path) as tif:
+ for l in range(1,levs):
+ yield tif.series[0].levels[l].pages[ch].asarray()
+
+def write_pyramid(src_img_path,
+ tasks_table,
+ outdir,
+ levels,
+ sub_lvls_dims,
+ file_name,
+ src_data_type,
+ is_src_pyramid=False,
+ save_ram=False
+ ):
+
+ outdir.mkdir(parents=True, exist_ok=True)
+ out_file_path=outdir / file_name
+ sub_levels=levels-1
+
+ total_operations=tasks_table['processed'].values.sum()#Count True values
+ count=1
+
+ with tifff.TiffWriter(out_file_path,bigtiff=True) as tif:
+ #write first the original resolution image,i.e. first layer
+ for _,channel in tasks_table.iterrows():
+ if channel.processed:
+ operation_count=f"({count}/{total_operations})"
+ print(f"\n {operation_count} Calculating subtraction of background {channel.background} from {channel.marker_name} signal:")
+ first_layer=subtract_channels(src_img_path, channel.ind, channel.bg_idx, channel.factor, (4096,4096), src_data_type,operation_count)
+ pyramid_action="calculate"
+ count+=1
+ else:
+ first_layer=tifff.imread(src_img_path,series=0,level=0,key=int(channel.ind))
+
+ if (save_ram or not is_src_pyramid):
+ pyramid_action="calculate"
+ else:
+ pyramid_action="extract"
+
+ tif.write(
+ first_layer,
+ subifds=sub_levels,
+ tile=(256, 256),
+ photometric='minisblack',
+ compression="lzw"
+ )
+
+ if pyramid_action=="calculate":
+ if save_ram:
+ pyramid=pyramid_save_ram(first_layer,sub_levels)
+ else:
+ pyramid=pyramidal_levels(first_layer,sub_levels,sub_lvls_dims)
+
+ elif pyramid_action=="extract":
+ pyramid=extract_sublevels_from_tiff(src_img_path,int(channel.ind),levels)
+
+
+ for sub_layer in pyramid:
+ tif.write(
+ sub_layer,
+ subfiletype=1,
+ tile=(256, 256),
+ photometric='minisblack',
+ compression="lzw"#lzw works better when saving channel-by-channel and jpeg 2000 when saving the whole stack at once
+ )
+
+ return out_file_path
+
+@memocron
+def main(version):
+ args=CLI.get_args()
+ in_path = args.root
+ out_path = args.output
+
+ # 0) Validate input_path is not the same as output_path,pixel data is read into RAM lazily, cannot overwrite input file
+ assert out_path != in_path
+
+ # 1) Extract image properties
+ src_props = extract_img_props(in_path, args.pixel_size,)
+ # 2) Modify pyramid_levels if required
+ if src_props["pyramid"]:
+ levels=src_props["levels"]
+ else:
+ levels=args.pyramid_levels
+
+ # 3) Read/Create markers table and update it to include the information of the processing tasks
+ if args.comet_metadata:
+ registration_marker="DAPI"
+ meta_table=meta_from_file(in_path,registration_marker)
+ markers = process_markers( assign_background(meta_table,rmv_ref=True,ref_marker=registration_marker) )
+
+ elif args.markers:
+ markers = process_markers(pd.read_csv(args.markers))
+
+ markers_updated=markers.loc[ markers.keep]
+ #4) Write updated markers.csv without appended columns. This file contains the markers information of the final image stack
+ markers_preview = markers_updated.drop(columns=['keep','ind','processed','factor','bg_idx'])
+ markers_preview["channel_number"] = np.arange(1, len(markers_preview)+1)
+ markers_preview.to_csv(args.markerout, index=False)
+
+ logger.info("\nTASKS PREVIEW:\n{}",markers_updated)
+ tasks=1
+ for _,channel in markers_updated.iterrows():
+ if channel.processed:
+ print(f"\n(Task_{tasks}): background subtraction, Channel {channel.marker_name} (Background {channel.background})")
+ tasks+=1
+
+ #5) Calculate subtractions and write output file
+ out_file_name=f'{ (in_path.stem).split(".ome")[0] }_backsub.ome.tif'
+ logger.info(f"\nTASKS PROGRESS" )
+
+ print(f"\nCommencing writing of pyramidal ome.tif file into {out_path / out_file_name}")
+ print(f"\nCommencing subtraction tasks\n")
+
+ pyramid_abs_path=write_pyramid(
+ in_path,
+ markers_updated,
+ out_path,
+ levels,
+ src_props["sub_levels_dims"],
+ out_file_name,
+ src_props["data_type"],
+ is_src_pyramid=src_props["pyramid"],
+ save_ram=args.save_ram
+ )
+
+ #6) Write metadata in OME format into the pyramidal file
+ channel_names=markers_updated["marker_name"].tolist()
+ ome_xml=ome_writer.create_ome(channel_names,src_props,version)
+ tifff.tiffcomment(pyramid_abs_path, ome_xml.encode("utf-8"))
+
+
+
+ logger.info(f'\nSCRIPT FINISHED PROCESSING TASKS ')
+ print(f'\nPyramidal image with {levels} levels was successfully written ')
+
+
+
+
+if __name__ == '__main__':
+ _version = 'v0.5.0'
+ main(_version)
+
+
+
+
+
+
diff --git a/backsub/metadata2markers_table.py b/backsub/metadata2markers_table.py
new file mode 100644
index 0000000..2451cb8
--- /dev/null
+++ b/backsub/metadata2markers_table.py
@@ -0,0 +1,129 @@
+import pathlib
+from ome_types import from_tiff
+import argparse
+import pandas as pd
+import numpy as np
+
+#CLI
+def get_args():
+ parser=argparse.ArgumentParser()
+ parser.add_argument('-i',
+ '--input_img',
+ required=True,
+ type=pathlib.Path,
+ help='absolute path of the input image stack (.tif)'
+ )
+
+ parser.add_argument('-o',
+ '--output_dir',
+ required=True,
+ type=pathlib.Path,
+ help='absolute path of the directory where the output .csv file will be written'
+ )
+
+ parser.add_argument('-fn',
+ '--output_file_name',
+ required=False,
+ type=str,
+ default="markers.csv",
+ help='name of the csv file'
+ )
+
+ parser.add_argument('-rr',
+ '--remove_background_references',
+ required=False,
+ action='store_true',
+ help='setup the removal of all reference background channels and all DAPI except the first occurrence.'
+ )
+
+ parser.add_argument('-rm',
+ '--registration_marker',
+ required=False,
+ type=str,
+ default="DAPI",
+ help='name of the csv file'
+ )
+
+
+ args=parser.parse_args()
+ return args
+
+def meta_from_file(src_img_path,ref_marker_name):
+ #Fetch metadata object
+ ome=from_tiff(src_img_path)
+ #Fetch image attributes from ome
+ ch_names = [ element.name for element in ome.images[0].pixels.channels ]
+ exp_times = [ element.exposure_time for element in ome.images[0].pixels.planes ]
+ #cycles=[int(element.attributes["CycleID"])+1 for element in ome.structured_annotations[0].value.any_elements[0].children]
+ filters= [ element.attributes["FluorescenceChannel"] for element in ome.structured_annotations[0].value.any_elements[0].children ]
+ background=[None if ref_marker_name in element
+ else element for element in filters]
+
+ aux_dict={"channel_number":list(range(1,len(ch_names)+1)),
+ #"cycle_number":cycles,
+ "marker_name":ch_names,
+ "Filter":filters,
+ "background":background,
+ "exposure":exp_times
+ }
+
+ df=pd.DataFrame(aux_dict)
+ return df
+
+
+def assign_background(df,rmv_ref=False,ref_marker="DAPI"):
+ #Create column ["backsub_process"] indicating which rows will be processed with backsub
+ filters_=df.Filter.unique().tolist()
+ #Strings corresponding to filters/background names are set to False, since the are not processed
+ backsub_process=df["marker_name"].replace(filters_,value=False,regex=True)
+ #Marker_name corresponding to signal will be set to True for processing
+ backsub_process=np.where(backsub_process==False,False,True)
+ df.insert(df.shape[1],"backsub_process",backsub_process)
+
+ #Assign the latest mention of the autofluorescence channel to the correspondent row in the background column
+ rename_background=[]
+ #List with row indices of background channels to be removed
+
+ for idx,row in df.iterrows():
+
+ if row.backsub_process:
+ previous_channels=reversed(df.iloc[:idx].marker_name.to_list())
+ for element in previous_channels:
+ #Supposes background name is a subset of the channel/marker name
+ if row.background in element:
+ rename_background.append(element)
+ break
+ else:
+ rename_background.append(None)
+
+ df.drop(columns=["backsub_process"],inplace=True)
+ df.background=rename_background
+ if rmv_ref:
+ remove_val=len(df)*[""]
+ df.insert(df.shape[1],"remove",remove_val)
+ df.loc[ df["background"].isnull() , ["remove"] ]="TRUE"
+ first_ref_marker=df[df.marker_name== ref_marker].index[0]
+ df.loc[first_ref_marker,"remove"]=""
+ return df
+
+
+
+def main():
+ args=get_args()
+ img_path=args.input_img
+ out_dir=args.output_dir
+ file_name=args.output_file_name
+ global_ref_marker=args.registration_marker
+
+ df=meta_from_file(img_path,global_ref_marker)
+ df_updated=assign_background(df,args.remove_background_references,global_ref_marker)
+ df_updated.to_csv( out_dir/file_name ,index=False)
+
+
+if __name__ == '__main__':
+ main()
+
+
+
+
+
diff --git a/backsub/ome_schema.py b/backsub/ome_schema.py
new file mode 100644
index 0000000..4c19bb9
--- /dev/null
+++ b/backsub/ome_schema.py
@@ -0,0 +1,166 @@
+#!/usr/bin/python
+import ome_types
+from ome_types.model import OME,Image,Pixels,TiffData,Channel,Plane
+import platform
+
+
+def INPUTS(frame):
+ """
+ This function creates a dictionary with the metadata of the tiles.
+ Args:
+ frame (pd.DataFrame): dataframe containing the metadata of the tiles.
+ conformed_markers (list): list of tuples with the name of the markers and their corresponding fluorophore.
+ Returns:
+ dict: dictionary with the metadata of the tiles.
+ """
+ inputs=frame.to_dict('list')
+
+ return inputs
+
+
+def TIFF_array(no_of_channels, inputs={'offset':0}):
+ """
+ This function creates a list of TIFFData objects.
+ Args:
+ no_of_channels (int): number of channels.
+ inputs (dict): dictionary with the metadata of the tiles.
+ Returns:
+ list: list of TIFFData objects.
+ """
+ TIFF = [
+ TiffData(
+ first_c=ch,
+ ifd=n,
+ plane_count=1
+ )
+ for n,ch in enumerate(range(0,no_of_channels), start=inputs['offset'])
+ ]
+
+ return TIFF
+
+
+def PLANE_array(no_of_channels, inputs):
+ """
+ This function creates a list of Plane objects.
+ Args:
+ no_of_channels (int): number of channels.
+ inputs (dict): dictionary with the metadata of the tiles.
+ Returns:
+ list: list of Plane objects.
+ """
+
+ PLANE = [
+ Plane(
+ the_c=ch,
+ the_t=0,
+ the_z=0,
+ position_x= inputs['position_x'][ch] if 'position_x' in inputs.keys() else 0,
+ position_y= inputs['position_y'][ch] if 'position_y' in inputs.keys() else 0,
+ position_z=0,
+ position_x_unit= inputs['position_x_unit'][ch] if 'position_x_unit'in inputs.keys() else "pixel" ,
+ position_y_unit= inputs['position_y_unit'][ch] if 'position_y_unit'in inputs.keys() else "pixel"
+ )
+ for ch in range(0,no_of_channels)
+ ]
+
+ return PLANE
+
+
+def CHANN_array(no_of_channels, inputs):
+ """
+ This function creates a list of Channel objects.
+ Args:
+ no_of_channels (int): number of channels.
+ inputs (dict): dictionary with the metadata of the tiles.
+ Returns:
+ list: list of Channel objects.
+ """
+
+ CHANN = [
+ Channel(
+ id=f"Channel:{str(ch)}", # 'Channel:{y}:{x}:{marker_name}'.format(x=ch,y=100+int( inputs['tile'][ch] ) ,marker_name=inputs['marker'][ch] )
+ name=inputs["name"][ch],
+ color=(255,255,255)
+ )
+ for ch in range(0,no_of_channels)
+ ]
+
+ return CHANN
+
+
+def PIXELS_array(chann_block, plane_block, tiff_block, inputs):
+ """
+ This function creates a Pixels object.
+ Args:
+ chann_block (list): list of Channel objects.
+ plane_block (list): list of Plane objects.
+ tiff_block (list): list of TIFFData objects.
+ inputs (dict): dictionary with the metadata of the tiles.
+ Returns:
+ Pixels: Pixels object.
+ """
+
+ PIXELS = Pixels(
+ id=f"Pixels:{inputs['tile'][0]}",
+ dimension_order='XYCZT',
+ size_c=len(chann_block),
+ size_t=1,
+ size_x=inputs['size_x'][0],
+ size_y=inputs['size_y'][0],
+ size_z=1,
+ type=inputs['type'][0],#bit_depth
+ big_endian=False,
+ channels=chann_block,
+ interleaved=False,
+ physical_size_x=inputs['physical_size_x'][0],
+ physical_size_x_unit=inputs['physical_size_x_unit'][0],
+ physical_size_y=inputs['physical_size_y'][0],
+ physical_size_y_unit=inputs['physical_size_y_unit'][0],
+ physical_size_z=1.0,
+ planes=plane_block,
+ significant_bits=inputs['significant_bits'][0],
+ tiff_data_blocks=tiff_block
+ )
+
+ return PIXELS
+
+
+def IMAGE_array(pixels_block, imageID):
+ """
+ This function creates an Image object.
+ Args:
+ pixels_block (Pixels): Pixels object.
+ imageID (int): identifier of the image.
+ Returns:
+ Image: Image object.
+ """
+
+ IMAGE = Image(
+ id =f'Image:{imageID}',
+ pixels=pixels_block
+ )
+
+ return IMAGE
+
+
+def OME_metadata(image_block,software):
+ """
+ This function creates an OME object.
+ Args:
+ image_block (list): list of Image objects.
+ Returns:
+ OME: OME object.
+ """
+ ome = OME()
+ ome.creator = " ".join([software,
+ ome_types.__name__,
+ ome_types.__version__,
+ '/ python version-',
+ platform.python_version()
+ ]
+ )
+
+ ome.images = image_block
+ ome_xml = ome_types.to_xml(ome)
+
+ return ome, ome_xml
diff --git a/backsub/ome_writer.py b/backsub/ome_writer.py
new file mode 100644
index 0000000..ba698ff
--- /dev/null
+++ b/backsub/ome_writer.py
@@ -0,0 +1,44 @@
+import ome_schema as schema
+import pandas as pd
+
+
+def create_ome(conformed_markers,info,software_version):
+ """
+ This function creates an OME-XML file from a pandas dataframe containing the metadata of the tiles.
+ Args:
+ tile_info (pd.DataFrame): dataframe containing the metadata of the tiles.
+ conformed_markers (list): list with the name of the markers in the corresponding order of their appearance in the ome.tif file .
+ Returns:
+ str: OME-XML file.
+ """
+ software=f'backsub {software_version}'
+ no_of_channels = len(conformed_markers)
+ tile_info_dict={
+ "tile": no_of_channels *[1],
+ "name":conformed_markers ,
+ "type": no_of_channels *[ info["data_type"] ],
+ "size_x":no_of_channels *[ info["size_x"] ] ,
+ "size_y":no_of_channels *[info["size_y"]],
+ "physical_size_x": no_of_channels *[info["pixel_size"]],
+ "physical_size_x_unit": no_of_channels *[info["pixel_size_unit"]],
+ "physical_size_y": no_of_channels *[info["pixel_size"]],
+ "physical_size_y_unit": no_of_channels *[info["pixel_size_unit"]],
+ "significant_bits": no_of_channels *["16"]
+ }
+ tile_info=pd.DataFrame(tile_info_dict)
+
+ grouped_tiles = tile_info.groupby(['tile'])
+
+ tiles_counter = 0
+ image = []
+ for tileID, frame in grouped_tiles:
+ metadata = schema.INPUTS(frame)
+ tiff = schema.TIFF_array(no_of_channels, inputs={'offset': no_of_channels * tiles_counter})
+ plane = schema.PLANE_array(no_of_channels, metadata)
+ channel = schema.CHANN_array(no_of_channels, metadata)
+ pixels = schema.PIXELS_array(channel, plane, tiff, metadata)
+ image.append(schema.IMAGE_array (pixels, tiles_counter))
+ tiles_counter += 1
+ ome, ome_xml = schema.OME_metadata(image,software)
+
+ return ome_xml
diff --git a/environment.yml b/environment.yml
index 957a760..db53f01 100644
--- a/environment.yml
+++ b/environment.yml
@@ -2,17 +2,15 @@ name: backsub
channels:
- conda-forge
- defaults
- - anaconda
dependencies:
- - "python=3.9"
- - "openslide=3.4.1"
- - "scikit-image=0.19.2"
- - "numexpr=2.8.3"
- - "tifffile=2022.8.12"
- - "scipy=1.9.3"
- - "pandas=2.1.1"
- - "zarr=2.3.2"
- - procps-ng
- - pip
+ - python=3.12.11
+ - pandas=2.3.1
+ - tifffile=2025.6.11
+ - scikit-image=0.25.2
+ - dask=2025.7.0
+ - ome-types=0.6.0
+ - numpy=2.3.2
+ - loguru=0.7.3
+ - dask-image=2024.5.3
- pip:
- - palom
\ No newline at end of file
+ - "zarr==3.1.1"
\ No newline at end of file