From 181976dfea551220571cd66ebaa6a2ca47505865 Mon Sep 17 00:00:00 2001 From: Kimberly Meechan <kimberly.meechan@embl.de> Date: Tue, 28 Jan 2020 10:19:08 +0100 Subject: [PATCH] split main function into two --- .../extension/attributes/morphology_impl.py | 154 +++++++++++++++++- 1 file changed, 147 insertions(+), 7 deletions(-) diff --git a/scripts/extension/attributes/morphology_impl.py b/scripts/extension/attributes/morphology_impl.py index 0298a9d..d9e33fe 100644 --- a/scripts/extension/attributes/morphology_impl.py +++ b/scripts/extension/attributes/morphology_impl.py @@ -99,7 +99,7 @@ def load_data(ds, row, scale): # load the data from the bounding box return ds[bb] - +# @profile def morphology_row_features(mask, scale): # Calculate stats from skimage ski_morph = regionprops(mask.astype('uint8')) @@ -131,7 +131,7 @@ def morphology_row_features(mask, scale): return (volume_in_microns, extent, equiv_diameter, major_axis, minor_axis, surface_area, sphericity, max_radius) - +# @profile def intensity_row_features(raw, mask): intensity_vals_in_mask = raw[mask] # mean and stdev - use float64 to avoid silent overflow errors @@ -146,7 +146,7 @@ def intensity_row_features(raw, mask): return mean_intensity, st_dev, median_intensity, interquartile_range_intensity, total - +# @profile def radial_intensity_row_features(raw, mask, scale, stops=(0.0, 0.25, 0.5, 0.75, 1.0)): result = () @@ -163,7 +163,7 @@ def radial_intensity_row_features(raw, mask, scale, stops=(0.0, 0.25, 0.5, 0.75, return result - +# @profile def texture_row_features(raw, mask): # errors if there are small, isolated spots (because I'm using ignore zeros as true) # so here remove components that are < 10 pixels @@ -187,7 +187,7 @@ def texture_row_features(raw, mask): return tuple(hara) - +# @profile def radial_distribution(edt, mask, stops=(0.0, 0.25, 0.5, 0.75, 1.0)): result = () @@ -201,7 +201,7 @@ def radial_distribution(edt, mask, stops=(0.0, 0.25, 0.5, 0.75, 1.0)): return result - +# @profile def chromatin_row_features(chromatin, edt, raw, scale_chromatin): result = () @@ -224,7 +224,7 @@ def chromatin_row_features(chromatin, edt, raw, scale_chromatin): return result - +# @profile # compute morphology (and intensity features) for label range def morphology_features_for_label_range(table, ds, ds_raw, ds_chromatin, @@ -404,6 +404,145 @@ def compute_morphology_features(table, segmentation_path, raw_path, return stats + +def morphology_impl_nucleus (nucleus_segmentation_path, raw_path, chromatin_path, + table, + min_size, max_size, + max_bb, + nucleus_resolution, chromatin_resolution, + nucleus_seg_scale, raw_scale, chromatin_scale, + label_starts, label_stops): + + # keys for the different scales + nucleus_seg_key_full = 't00000/s00/0/cells' + nucleus_seg_key = 't00000/s00/%i/cells' % nucleus_seg_scale + raw_key_full = 't00000/s00/0/cells' + raw_key = 't00000/s00/%i/cells' % raw_scale + chromatin_key_full = 't00000/s00/0/cells' + chromatin_key = 't00000/s00/%i/cells' % chromatin_scale + + # get scale factors + scale_factor_nucleus_seg = get_scale_factor(nucleus_segmentation_path, nucleus_seg_key_full, nucleus_seg_key, + nucleus_resolution) + + if raw_path is not None: + log("Have raw path; compute intensity features") + # NOTE for now we can hard-code the resolution for the raw data here, + # but we might need to change this if we get additional dataset(s) + raw_resolution = [0.025, 0.01, 0.01] + scale_factor_raw = get_scale_factor(raw_path, raw_key_full, raw_key, raw_resolution) + else: + log("Don't have raw path; do not compute intensity features") + scale_factor_raw = None + + if chromatin_path is not None: + log("Have chromatin path; compute chromatin features") + scale_factor_chromatin = get_scale_factor(chromatin_path, chromatin_key_full, chromatin_key, + chromatin_resolution) + else: + log("Don't have chromatin path; do not compute chromatin features") + scale_factor_chromatin = None + + # remove zero label if it exists + table = table.loc[table['label_id'] != 0, :] + + # filter by size + if min_size is not None or max_size is not None: + table = filter_table(table, min_size, max_size) + log("Number of labels after size filter: %i" % table.shape[0]) + + # filter by bounding box size + if max_bb is not None: + table = filter_table_bb(table, max_bb) + log("Number of labels after bounding box size filter %i" % table.shape[0]) + + log("Computing morphology features") + stats = compute_morphology_features(table, segmentation_path, raw_path, + chromatin_path, exclude_nuc_path, + seg_key, raw_key, chromatin_key, + exclude_key, + scale_factor_seg, scale_factor_raw, + scale_factor_chromatin, + scale_factor_exclude, + label_starts, label_stops) + return stats + +def morphology_impl_cell (cell_segmentation_path, raw_path, + nucleus_segmentation_path, + table, mapping_path, + region_mapping_path, + min_size, max_size, + max_bb, + cell_resolution, nucleus_resolution, + cell_seg_scale, raw_scale, nucleus_seg_scale, + label_starts, label_stops): + + # keys for the different scales + cell_seg_key_full = 't00000/s00/0/cells' + cell_seg_key = 't00000/s00/%i/cells' % cell_seg_scale + raw_key_full = 't00000/s00/0/cells' + raw_key = 't00000/s00/%i/cells' % raw_scale + nucleus_seg_key_full = 't00000/s00/0/cells' + nucleus_seg_key = 't00000/s00/%i/cells' % nucleus_seg_scale + + # get scale factors + cell_seg_scale_factor = get_scale_factor(cell_segmentation_path, cell_seg_key_full, cell_seg_key, cell_resolution) + + if raw_path is not None: + log("Have raw path; compute intensity features") + # NOTE for now we can hard-code the resolution for the raw data here, + # but we might need to change this if we get additional dataset(s) + raw_resolution = [0.025, 0.01, 0.01] + scale_factor_raw = get_scale_factor(raw_path, raw_key_full, raw_key, raw_resolution) + else: + log("Don't have raw path; do not compute intensity features") + scale_factor_raw = None + + if nucleus_segmentation_path is not None: + log("Have nucleus path; exclude nucleus for intensity measures") + scale_factor_nucleus = get_scale_factor(nucleus_segmentation_path, nucleus_seg_key_full, nucleus_seg_key, + nucleus_resolution) + else: + log("Don't have exclude path; don't exclude nucleus area for intensity measures") + scale_factor_nucleus = None + + # remove zero label if it exists + table = table.loc[table['label_id'] != 0, :] + + # filter to only keep cells with assigned nuclei) + if mapping_path is not None: + log("Have mapping path %s" % mapping_path) + table = filter_table_from_mapping(table, mapping_path) + log("Number of labels after filter with mapping: %i" % table.shape[0]) + + # filter to exclude certain regions + if region_mapping_path is not None: + log("Have region mapping path %s" % region_mapping_path) + table = filter_table_region(table, region_mapping_path) + log("Number of labels after region filter: %i" % table.shape[0]) + + # filter by size of object (no. pixels) + if min_size is not None or max_size is not None: + table = filter_table(table, min_size, max_size) + log("Number of labels after size filter: %i" % table.shape[0]) + + # filter by bounding box size + if max_bb is not None: + table = filter_table_bb(table, max_bb) + log("Number of labels after bounding box size filter %i" % table.shape[0]) + + log("Computing morphology features") + stats = compute_morphology_features(table, segmentation_path, raw_path, + chromatin_path, exclude_nuc_path, + seg_key, raw_key, chromatin_key, + exclude_key, + scale_factor_seg, scale_factor_raw, + scale_factor_chromatin, + scale_factor_exclude, + label_starts, label_stops) + return stats + + def morphology_impl(segmentation_path, raw_path, chromatin_path, exclude_nuc_path, table, mapping_path, @@ -530,5 +669,6 @@ def morphology_impl(segmentation_path, raw_path, chromatin_path, return stats + if __name__ == '__main__': pass -- GitLab