-
Kai Schleicher authoredKai Schleicher authored
2d_spots_in_fibers.py 15.87 KiB
#@ OpService ops
#@ CommandService command
#@ RoiManager rm
#@ File (label="select image", description="select your input image", style=file) path_to_image
#@ Integer (label="select image file series", description="leave 1 if not needed", value=1, min=1, max=20, stepSize=1, persist=false, style=slider) series_number
#@ File (label="Fiber Segmentation", description="select the Fiber Segmentation RoiSet .zip file", style=file) fiber_segmentation_roiset
#@ Integer (label="DAPI channel", description="select the DAPI channel", min=1, max=20, stepSize=1, persist=true, style=slider) dapi_channel_number
#@ Integer (label="DAPI threshold", description="0 = Auto", min=0) threshold
#@ String (label="processing channels", description="comma separated list of channels, e.g. 2,3,6,7", value="1,2,3,4") processing_channels_string
#@ String (label="Spot detection quality threshold", description="comma separated list of values, one per channel", value="80,80,80,80" ) quality_thresholds_string
#@ Double (label="Minimum spot diameter", description="smaller spots will be discarded", value=0.3, min=1, max=10, stepSize=0.1) spot_diameter
# trackmate imports
from fiji.plugin.trackmate import Settings
from fiji.plugin.trackmate import TrackMate
from fiji.plugin.trackmate.detection import DogDetectorFactory
from fiji.plugin.trackmate.action import LabelImgExporter
# ij imports
from ij import IJ
from ij.plugin import ImageCalculator
from ij.measure import ResultsTable
from ij.measure import Measurements as M
from ij.plugin.filter import Analyzer
# Bio-Formats imports
from loci.plugins import BF
from loci.common import Region
from loci.plugins.in import ImporterOptions
from loci.formats import ImageReader, MetadataTools, Memoizer
# MorpholibJ imports
from inra.ijpb.binary import BinaryImages
# BIOP imports
from ch.epfl.biop.ij2command import Labels2Rois
# python imports
import os
def get_parentdir_filename_ext_from_path(path):
""" returns parent directory, filename, and file extension from a given path
Parameters
----------
path : String
the full path to a file
Returns
-------
String
the parent directory, the filename without extension and the extension split after the dot
"""
parent_dir = os.path.dirname(path)
base = os.path.basename(path)
(filename, ext) = os.path.splitext(base)
return (parent_dir, filename, ext)
def write_bf_memoryfile(path_to_file):
"""Write a BF memo-file so subsequent access to the same file is faster.
The Bio-Formats memo-file is written next to the image file (i.e. in the
same folder as the given file).
Parameters
----------
string
The full path to the image file.
"""
reader = Memoizer(ImageReader())
reader.setId(path_to_file)
reader.close()
def get_ome_metadata(path_to_image, series_number):
"""Get some metadata fields using BF from an image on disk from a given series.
Gets number of channels, the image width, image height and bit depth.
Parameters
----------
path_to_image : str
full path to the image
series_number : int
the Bio-Formats series number
Returns
-------
int
n_channels : the number of channels
int
image_width : the image width in px
int
image_height : the image height in px
int
bit_depth : the image bit depth
"""
series_number -= 1 # python starts counting from 0
reader = ImageReader()
ome_meta = MetadataTools.createOMEXMLMetadata()
reader.setMetadataStore(ome_meta)
reader.setId(path_to_image)
reader.setSeries(series_number)
n_channels = reader.getSizeC()
image_width = reader.getSizeX()
image_height = reader.getSizeY()
bit_depth = reader.getBitsPerPixel()
reader.close()
return n_channels, image_width, image_height, bit_depth
def BFopen_image(path_to_image, channel_number, series_number, region=None, z_slice_number=None):
"""Use Bio-Formats to open a given image or subset.
Parameters
----------
path_to_image : str
full path to the image file on disk
channel_number : int
the number of the channel to open
series_number : int
the Bio-Formats series number to open
region : list, optional
Bio-Formats crop region, by default None. format = [start_x, start_y, width in px, height in px]
z_slice_number : int, optional
the number of the z-slice to open, by default None
Returns
-------
ImagePlus
the requested image
"""
series_number -= 1 # python starts counting from 0
channel_number -= 1 # python starts counting from 0
options = ImporterOptions()
options.setId(path_to_image)
options.setColorMode(ImporterOptions.COLOR_MODE_GRAYSCALE)
options.setSplitChannels(False)
options.setSplitFocalPlanes(False)
options.setSplitTimepoints(False)
options.setSeriesOn(series_number, True)
if region is not None:
options.setCrop(True)
options.setCropRegion(series_number, Region(region[0], region[1], region[2], region[3])) # x,y,w,h
options.setSpecifyRanges(True)
options.setCBegin(series_number, channel_number)
options.setCEnd(series_number, channel_number)
options.setCStep(series_number, 1)
if z_slice_number is not None:
z_slice_number -= 1 # python starts counting from 0
options.setZBegin(series_number, z_slice_number)
options.setZEnd(series_number, z_slice_number)
options.setZStep(series_number, 1)
imps = BF.openImagePlus(options)
return imps[0]
def run_trackmate_dog_spot_detector(imp, spot_diameter, quality_threshold):
"""Run TrackMates DoG detector with a given quality threshold
on a target image
Parameters
----------
imp : ImagePlus
the input image
quality_threshold : float or double
the spot detection quality threshold
Returns
-------
ImagePlus
a label image of the detected spots.
"""
cal = imp.getCalibration()
settings = Settings(imp)
settings.detectorFactory = DogDetectorFactory()
settings.detectorSettings['DO_SUBPIXEL_LOCALIZATION'] = True
settings.detectorSettings['RADIUS'] = spot_diameter / 2 # type = double
settings.detectorSettings['TARGET_CHANNEL'] = 0
settings.detectorSettings['THRESHOLD'] = quality_threshold # type = double
settings.detectorSettings['DO_MEDIAN_FILTERING'] = False
trackmate = TrackMate(settings)
trackmate.computeSpotFeatures(False)
trackmate.computeTrackFeatures(False)
if not trackmate.checkInput():
print( str( trackmate.getErrorMessage() ) )
if not trackmate.process():
print( str( trackmate.getErrorMessage() ) )
exportSpotsAsDots = True
exportTracksOnly = False
labelIdPainting = LabelImgExporter.LabelIdPainting.LABEL_IS_INDEX_MOVIE_UNIQUE
label_imp = LabelImgExporter.createLabelImagePlus(trackmate, exportSpotsAsDots, exportTracksOnly, labelIdPainting)
label_imp.setCalibration(cal)
imp.close()
return label_imp
def get_threshold_from_method(imp, method):
"""returns the threshold value of chosen IJ AutoThreshold method from a stack (hyperstack?)
Parameters
----------
imp : ImagePlus (maybe also accepts ij2 dataset, img?)
the image from which to get the threshold value
method : string
the AutoThreshold method to use, e.g. "otsu", "yen", ...
Returns
-------
int
"""
histogram = ops.run("image.histogram", imp)
threshold_value = ops.run("threshold.%s" % method, histogram)
threshold_value = int(round(threshold_value.get()))
return threshold_value
def convert_to_binary(imp, threshold=1, scale_binary=True):
"""Convert grayscale image to a mask / binary image
Parameters
----------
imp : ImagePlus
the input image
threshold : int, optional
the threshold above which pixels are considered foreground, by default 1
scale_binary : bool, optional
scale the mask image from 0/255 to 0/1, by default True
Returns
-------
ImagePlus
the mask / binary image
"""
mask_imp = imp.duplicate()
bitdepth = mask_imp.getBitDepth()
if bitdepth > 16.0:
IJ.run("Conversions...", " ")
IJ.run(mask_imp, "16-bit", "")
IJ.run("Conversions...", "scale")
IJ.setRawThreshold(mask_imp, threshold, (2**bitdepth)-1, None)
IJ.run(mask_imp, "Convert to Mask", "")
if scale_binary == True:
IJ.run(mask_imp, "Divide...", "value=255")
IJ.run(mask_imp, "Enhance Contrast", "saturated=0.35")
mask_imp.changes = False
return mask_imp
def save_image_as_IJtif(imp, filename, suffix, target):
"""Save the image as ImageJ-Tiff
Parameters
----------
imp : ImagePlus
the input Image to be saved
filename : str
the image filename without extension
suffix : int, float or str
a suffix that will be added to the filename as a string
target : str
the target directory in with to save the image
"""
savename = filename + "_" + str(suffix) + ".tif"
savepath = os.path.join(target, savename)
IJ.log("now saving: " + str(savepath))
IJ.saveAs(imp, "Tiff", savepath)
def convert_labelimage_to_binary(label_imp, scale_binary=True):
"""Convert a label image to a mask / binary image
Parameters
----------
label_imp : ImagePlus
the input label image
scale_binary : bool, optional
scale the mask image from 0/255 to 0/1, by default True
Returns
-------
ImagePlus
the mask / binary image
"""
binary_imp = BinaryImages.binarize(label_imp)
if scale_binary == True:
IJ.run(binary_imp, "Divide...", "value=255")
return binary_imp
def measure_intensity_sum(imp, rm):
"""
Measure theraw integrated intensity for all rois in target imp.
Parameters
----------
imp : ImagePlus
The image from which the intensity will be measured.
rm : RoiManager
The ROI Manager containing the regions of interest to be analyzed.
Returns
-------
label_id : list of str
A list of labels corresponding to each ROI.
sum_intensity : list of int
A list of summed integrated intensities for each ROI.
Notes
-----
The results are stored in a `ij.ResultsTable`,
from which the raw integrated density values are extracted.
Example
-------
>>> labels, intensities = measure_intensity_sum(imp, roi_manager)
"""
rt_ = ResultsTable()
options = M.INTEGRATED_DENSITY
an = Analyzer (imp, options, rt_)
an.setPrecision(0)
label_id = []
sum_intensity = []
for index, roi in enumerate(rm.getRoisAsArray()):
imp.setRoi(roi)
an.measure()
label_id.append(roi.getName())
sum_intensity.append(int(rt_.getColumn("RawIntDen")[index]))
return label_id, sum_intensity
def add_results_to_resultstable(results_table, column, values):
"""Add values to the ResultsTable starting from row 0 of a given column.
Parameters
----------
results_table : ij.measure.ResultsTable
a reference of the IJ-ResultsTable
column : string
the column in which to add the values
values : list(int, double or float)
array with values to be added
"""
for index, value in enumerate(values):
results_table.setValue(column, index, value)
def save_results_table(results_table, filename, suffix, target):
"""Save an IJ-ResultsTable as txt to a given location
Parameters
----------
results_table : IJ-ResultsTable
an IJ-ResultsTable
filename : str
filename to include in the tables name
suffix : str, int or float
a suffix to be included in the tables name
target : str
the target directory in which to save the table
"""
savename = filename + "_" + str(suffix) + ".txt"
savepath = os.path.join(target, savename)
results_table.save(savepath)
def save_labelimage_as_ijroiset(label_imp, rm, filename, suffix, target):
"""Save all labels i a label Images as IJ-ROIset to target location.
Only 2D, requires PTBIOP update site.
Parameters
----------
label_imp : ImagePlus
the input label image
filename : str
the filename to include in the ROIset name
suffix : str, int or float
a suffix to be included in the ROIset name
target : str
the target directory in which to save the table
"""
command.run( Labels2Rois , False , 'imp' , label_imp , 'rm', rm).get()
savename = filename + "_" + str(suffix) + ".zip"
savepath = os.path.join(target, savename)
rm.runCommand("Save", savepath)
rm.reset()
def load_rois_from_zip(path, rm):
"""Load ROIs from the given zip file and add them to the RoiManager.
Parameters
----------
path : string
Path to the ROI zip file.
rm : RoiManager
The ROI Manager.
"""
rm.runCommand("Open", path)
def close_images(list_of_imps):
"""Close given ImagePlus images
Parameters
----------
list_of_imps : list of ImagePlus
list of ImagePlus images to be closed
"""
for imp in list_of_imps:
imp.changes = False
imp.close()
path_to_image = str(path_to_image).replace("\\","/")
fiber_segmentation_roiset = str(fiber_segmentation_roiset).replace("\\","/")
parent_dir, filename, ext = get_parentdir_filename_ext_from_path(path_to_image)
write_bf_memoryfile(path_to_image)
n_channels, image_width, image_height, bit_depth = get_ome_metadata(path_to_image, series_number)
processing_channels_string = processing_channels_string.replace(" ", "")
processing_channels = processing_channels_string.split(",")
quality_thresholds_string = quality_thresholds_string.replace(" ", "")
quality_thresholds = quality_thresholds_string.split(",")
# threshold DAPI channel and convert to binary
dapi_channel = BFopen_image(path_to_image, dapi_channel_number, series_number)
if threshold <= 0:
threshold = get_threshold_from_method(dapi_channel, "otsu")
dapi_binary = convert_to_binary(dapi_channel, threshold)
# detect spots and count them per fiber
results_table = ResultsTable()
for index, channel in enumerate(processing_channels):
channel = int(channel)
quality_threshold = float(quality_thresholds[index])
spots_channel = BFopen_image(path_to_image, channel, series_number)
spots_label_imp = run_trackmate_dog_spot_detector(spots_channel, spot_diameter, quality_threshold)
save_labelimage_as_ijroiset(spots_label_imp, filename, "spots_ch" + str(channel), parent_dir)
spots_binary_imp = convert_labelimage_to_binary(spots_label_imp)
dapi_positive_spots_binary = ImageCalculator.run(spots_binary_imp, dapi_binary, "Multiply create")
dapi_negative_spots_binary = ImageCalculator.run(spots_binary_imp, dapi_positive_spots_binary, "Subtract create")
load_rois_from_zip(fiber_segmentation_roiset)
fiber_label_IDs, n_spots_per_fiber = measure_intensity_sum(spots_binary_imp, rm)
_, n_dapi_positive_spots_per_fiber = measure_intensity_sum(dapi_positive_spots_binary, rm)
_, n_dapi_negative_spots_per_fiber = measure_intensity_sum(dapi_negative_spots_binary, rm)
# TODO: the column fiber label ID needs to be added only one time. Maybe check before if column exists.
add_results_to_resultstable(results_table, "fiber label ID", fiber_label_IDs)
add_results_to_resultstable(results_table, "n spots channel " + str(channel), n_spots_per_fiber)
add_results_to_resultstable(results_table, "n spots dapi positive channel " + str(channel), n_dapi_positive_spots_per_fiber)
add_results_to_resultstable(results_table, "n spots dapi negative channel " + str(channel), n_dapi_negative_spots_per_fiber)
close_images([spots_channel, spots_label_imp, spots_binary_imp, spots_binary_imp, dapi_binary, dapi_positive_spots_binary, dapi_negative_spots_binary])
results_table.show("Spots")
save_results_table(results_table, filename, "Spots", parent_dir)
IJ.run("Collect Garbage", "")
IJ.log("DONE")