diff --git a/README.md b/README.md index 040dd22..2065e2a 100644 --- a/README.md +++ b/README.md @@ -92,4 +92,3 @@ pytest ## General outline of SPACEc analysis ![SPACEc](https://github.com/yuqiyuqitan/SPACEc/tree/master/docs/overview.png?raw=true "") - diff --git a/notebooks/0_tissue_extractor.ipynb b/notebooks/0_tissue_extractor.ipynb index ae356a0..c118cfe 100644 --- a/notebooks/0_tissue_extractor.ipynb +++ b/notebooks/0_tissue_extractor.ipynb @@ -74,7 +74,7 @@ "metadata": {}, "outputs": [], "source": [ - "root_path = \"/home/tim/Dokumente/GitHub/SPACEc/\" # replace with your path\n", + "root_path = \"/home/user/path/SPACEc/\" # replace with your path\n", "data_path = root_path + 'example_data/raw/' # where the data is stored\n", "\n", "# where you want to store the output\n", diff --git a/notebooks/1_cell_segmentation.ipynb b/notebooks/1_cell_segmentation.ipynb index 6e62227..2fe6ec3 100644 --- a/notebooks/1_cell_segmentation.ipynb +++ b/notebooks/1_cell_segmentation.ipynb @@ -77,11 +77,11 @@ "metadata": {}, "outputs": [], "source": [ - "root_path = \"/home/tim/Dokumente/GitHub/SPACEc/\"\n", + "root_path = \"/home/user/path/SPACEc/\" # inset your own path\n", "data_path = root_path + 'example_data/raw/' # where the data is stored\n", "\n", "# where you want to store the output\n", - "output_dir = root_path + 'example_data/output/'\n", + "output_dir = root_path + 'example_data/output/' # inset your own path\n", "os.makedirs(output_dir, exist_ok=True)" ] }, diff --git a/notebooks/2_preprocessing.ipynb b/notebooks/2_preprocessing.ipynb index d7e0fba..75fe2e2 100755 --- a/notebooks/2_preprocessing.ipynb +++ b/notebooks/2_preprocessing.ipynb @@ -57,7 +57,7 @@ "metadata": {}, "outputs": [], "source": [ - "root_path = \"/home/tim/Dokumente/GitHub/SPACEc/\" # replace with your path\n", + "root_path = \"/home/user/path/SPACEc/\" # inset your own path\n", "data_path = root_path + 'example_data/raw/' # where the data is stored\n", "\n", "# where you want to store the output\n", @@ -399,7 +399,7 @@ "# Identify the lowest 1% for cell size and nuclear marker intensity to get a better idea of potential segmentation artifacts.\n", "df_filt = sp.pp.filter_data(\n", " df_seg, \n", - " nuc_thres=1, # remove cells with DAPI intensity below threshold\n", + " nuc_thres=one_percent_nuc, # remove cells with DAPI intensity below threshold\n", " size_thres=one_percent_area, # remove cells with area below threshold\n", " nuc_marker=\"DAPI\", # name of nuclear marker\n", " cell_size = \"area\", # name of cell size column\n", diff --git a/notebooks/3_cell_annotation_STELLAR.ipynb b/notebooks/3_cell_annotation_STELLAR.ipynb index d059070..b3d45e0 100644 --- a/notebooks/3_cell_annotation_STELLAR.ipynb +++ b/notebooks/3_cell_annotation_STELLAR.ipynb @@ -67,7 +67,7 @@ "metadata": {}, "outputs": [], "source": [ - "root_path = \"/home/tim/Dokumente/GitHub/SPACEc/\" # replace with your path\n", + "root_path = \"/home/user/path/SPACEc/\" # inset your own path\n", "data_path = root_path + 'example_data/raw/' # where the data is stored\n", "\n", "# where you want to store the output\n", diff --git a/notebooks/3_cell_annotation_ml.ipynb b/notebooks/3_cell_annotation_ml.ipynb index 9fb1703..ab2d1e8 100644 --- a/notebooks/3_cell_annotation_ml.ipynb +++ b/notebooks/3_cell_annotation_ml.ipynb @@ -50,7 +50,7 @@ "metadata": {}, "outputs": [], "source": [ - "root_path = \"/home/tim/Dokumente/GitHub/SPACEc/\" # replace with your path\n", + "root_path = \"/home/user/path/SPACEc/\" # inset your own path\n", "data_path = root_path + 'example_data/raw/' # where the data is stored\n", "\n", "# where you want to store the output\n", diff --git a/notebooks/3_clustering.ipynb b/notebooks/3_clustering.ipynb index 9d56667..4f7f117 100644 --- a/notebooks/3_clustering.ipynb +++ b/notebooks/3_clustering.ipynb @@ -65,7 +65,7 @@ "metadata": {}, "outputs": [], "source": [ - "root_path = \"/home/tim/Dokumente/GitHub/SPACEc/\" # replace with your path\n", + "root_path = \"/home/user/path/SPACEc/\" # inset your own path\n", "data_path = root_path + 'example_data/raw/' # where the data is stored\n", "\n", "# where you want to store the output\n", diff --git a/notebooks/4_cell_neighborhood_analysis.ipynb b/notebooks/4_cell_neighborhood_analysis.ipynb index a000969..3b3a61c 100755 --- a/notebooks/4_cell_neighborhood_analysis.ipynb +++ b/notebooks/4_cell_neighborhood_analysis.ipynb @@ -55,7 +55,7 @@ "metadata": {}, "outputs": [], "source": [ - "root_path = \"/home/tim/Dokumente/SPACEc_Apr_2024/\"\n", + "root_path = \"/home/user/path/SPACEc/\" # inset your own path\n", "\n", "data_path = root_path + 'data/' # where the data is stored\n", "\n", diff --git a/notebooks/5_distance_permutation_analysis.ipynb b/notebooks/5_distance_permutation_analysis.ipynb index c7d9e89..49134d3 100644 --- a/notebooks/5_distance_permutation_analysis.ipynb +++ b/notebooks/5_distance_permutation_analysis.ipynb @@ -62,7 +62,7 @@ "metadata": {}, "outputs": [], "source": [ - "root_path = \"/home/tim/Dokumente/SPACEc_Apr_2024/\"\n", + "root_path = \"/home/user/path/SPACEc/\" # inset your own path\n", "\n", "data_path = root_path + 'data/' # where the data is stored\n", "\n", diff --git a/notebooks/6_patch_proximity_analysis.ipynb b/notebooks/6_patch_proximity_analysis.ipynb index 4f145c6..45ccb6d 100644 --- a/notebooks/6_patch_proximity_analysis.ipynb +++ b/notebooks/6_patch_proximity_analysis.ipynb @@ -59,7 +59,7 @@ "metadata": {}, "outputs": [], "source": [ - "root_path = \"/home/tim/Dokumente/SPACEc_Apr_2024/\" \n", + "root_path = \"/home/user/path/SPACEc/\" # inset your own path\n", "\n", "data_path = root_path + 'data/' # where the data is stored\n", "\n", @@ -99,7 +99,7 @@ } ], "source": [ - "adata = sc.read('/home/tim/Dokumente/GitHub/SPACEc/tests/data/processed/tonsil/1/adata_nn_demo_annotated_long.h5ad')\n", + "adata = sc.read(output_dir + 'adata_nn_demo_annotated.h5ad')\n", "adata" ] }, diff --git a/notebooks/7_TissUUmaps.ipynb b/notebooks/7_TissUUmaps.ipynb index 8e91de2..f5358a8 100644 --- a/notebooks/7_TissUUmaps.ipynb +++ b/notebooks/7_TissUUmaps.ipynb @@ -56,7 +56,7 @@ "metadata": {}, "outputs": [], "source": [ - "root_path = \"/home/tim/Dokumente/SPACEc_Apr_2024/\"\n", + "root_path = \"/home/user/path/SPACEc/\" # inset your own path\n", "\n", "data_path = root_path + 'data/' # where the data is stored\n", "\n", @@ -353,7 +353,7 @@ "image_list, csv_paths = sp.tl.tm_viewer(\n", " adata,\n", " images_pickle_path= output_dir + 'seg_output_tonsil2.pickle',\n", - " directory = \"/home/tim/Dokumente/SPACEc_Apr_2024/cache\",\n", + " directory = \"/home/user/path/SPACEc/cache\" # inset your own path where you want to cache your images for TM visualization (you can delete this once you are done with TM)\n", " region_column = \"unique_region\",\n", " region = \"reg002\",\n", " xSelector = \"y\",\n", diff --git a/src/spacec/_shared/segmentation.py b/src/spacec/_shared/segmentation.py index 9d28ed8..30d2678 100644 --- a/src/spacec/_shared/segmentation.py +++ b/src/spacec/_shared/segmentation.py @@ -5,6 +5,23 @@ def create_multichannel_tiff(input_dir, output_dir, output_filename): + """ + Create a multi-channel TIFF image from individual TIFF files. + + Parameters + ---------- + input_dir : str + Directory containing the input TIFF files. + output_dir : str + Directory to save the output TIFF file. + output_filename : str + Name of the output TIFF file. + + Returns + ------- + list of str + List of channel names. + """ # Get a list of all TIFF files in the input directory tiff_files = [f for f in os.listdir(input_dir) if f.endswith((".tiff", ".tif"))] @@ -35,18 +52,36 @@ def create_multichannel_tiff(input_dir, output_dir, output_filename): # combine multiple channels in one image and add as new image to image_dict with the name segmentation_channel def combine_channels(image_dict, channel_list, new_channel_name): + """ + Combine multiple channels into a single channel. + + Parameters + ---------- + image_dict : dict + Dictionary with channel names as keys and images as values. + channel_list : list of str + List of channel names to combine. + new_channel_name : str + Name of the new channel. + + Returns + ------- + dict + Updated dictionary with the new channel added. + """ + # Determine bit depth of input images + bit_depth = image_dict[channel_list[0]].dtype + # Create empty image new_image = np.zeros( - (image_dict[channel_list[0]].shape[0], image_dict[channel_list[0]].shape[1]) + (image_dict[channel_list[0]].shape[0], image_dict[channel_list[0]].shape[1]), + dtype=bit_depth, ) # Add channels to image as maximum projection for channel in channel_list: new_image = np.maximum(new_image, image_dict[channel]) - # generate greyscale image - new_image = np.uint8(new_image) - # Add image to image_dict image_dict[new_channel_name] = new_image @@ -62,6 +97,31 @@ def format_CODEX( stack=True, input_format="Multichannel", ): + """ + Format images based on the input format. + + Parameters + ---------- + image : ndarray or str + Input image or directory containing images. + channel_names : list of str, optional + List of channel names. + number_cycles : int, optional + Number of cycles in the CODEX format. + images_per_cycle : int, optional + Number of images per cycle in the CODEX format. + stack : bool, default=True + If True, stack the images in the list. + input_format : str, default="Multichannel" + Format of the input images. Options are "CODEX", "Multichannel", and "Channels". + + Returns + ------- + dict + Dictionary with channel names as keys and images as values. If `stack` is True and `input_format` is "CODEX", + also returns a stacked image as a numpy array. + """ + if input_format == "CODEX": total_images = number_cycles * images_per_cycle image_list = [None] * total_images # pre-allocated list diff --git a/src/spacec/helperfunctions/_general.py b/src/spacec/helperfunctions/_general.py index 3c0321f..de7d9a4 100644 --- a/src/spacec/helperfunctions/_general.py +++ b/src/spacec/helperfunctions/_general.py @@ -1053,6 +1053,18 @@ def is_dark(color): def check_for_gpu(): + """ + Check if a GPU is available for use by TensorFlow and PyTorch. + + This function checks if a GPU is available for use by TensorFlow and PyTorch. + It prints a message indicating whether a GPU is available for each library, + and returns a boolean indicating whether a GPU is available for PyTorch. + + Returns + ------- + bool + True if a GPU is available for PyTorch, False otherwise. + """ if tf.config.list_physical_devices("GPU"): print("GPU is available to Tensorflow") else: @@ -1061,3 +1073,4 @@ def check_for_gpu(): use_GPU = use_gpu() yn = ["GPU is not available to Pytorch", "GPU is available to Pytorch"] print(f"{yn[use_GPU]}") + return use_GPU diff --git a/src/spacec/plotting/_general.py b/src/spacec/plotting/_general.py index a1c1132..816b979 100644 --- a/src/spacec/plotting/_general.py +++ b/src/spacec/plotting/_general.py @@ -3816,6 +3816,34 @@ def cn_map( output_dir="./", rand_seed=1, ): + """ + Generates a CNMap plot using the provided data and parameters. + + Parameters + ---------- + adata : anndata.AnnData + Annotated data matrix. + cnmap_dict : dict + Dictionary containing graph, tops, e0, e1, and simp_freqs. + cn_col : str + Column name in adata to be used for color coding. + palette : dict, optional + Color palette to use for the plot. If None, a random color palette is generated. + figsize : tuple, optional + Size of the figure. Defaults to (40, 20). + savefig : bool, optional + Whether to save the figure or not. Defaults to False. + output_fname : str, optional + The filename for the saved figure. Required if savefig is True. Defaults to "". + output_dir : str, optional + The directory where the figure will be saved. Defaults to "./". + rand_seed : int, optional + Seed for random number generator. Defaults to 1. + + Returns + ------- + None + """ graph = cnmap_dict["g"] tops = cnmap_dict["tops"] e0 = cnmap_dict["e0"] @@ -3854,6 +3882,46 @@ def cn_map( c=col, zorder=-1, ) + # Dummy scatter plots for legend + freqs = simp_freqs * 10000 + max_size = max(freqs) + sizes = [ + round(max_size) / 4, + round(max_size) / 2, + round(max_size), + ] # Replace with the sizes you want in the legend + labels = [ + str(round(max_size / 100) / 4) + "%", + str(round(max_size / 100) / 2) + "%", + str(round(max_size / 100)) + "%", + ] # Replace with the labels you want in the legend + + # Add legend + legend_elements = [ + plt.Line2D( + [0], + [0], + marker="o", + color="w", + label=label, + markerfacecolor="black", + markersize=size**0.5, + ) + for size, label in zip(sizes, labels) + ] + + # Add first legend + legend1 = plt.legend( + handles=legend_elements, + loc="lower right", + title="Total frequency", + title_fontsize=30, + fontsize=30, + handlelength=6, + handletextpad=1, + bbox_to_anchor=(0.0, -0.15, 1.0, 0.102), + ) + if n in tops: plt.text( pos[n][0], @@ -3903,10 +3971,11 @@ def cn_map( ] # Add legend to bottom of plot + plt.gca().add_artist(legend1) plt.legend( handles=legend_patches, bbox_to_anchor=(0.0, -0.15, 1.0, 0.102), - loc="lower center", + loc="lower left", ncol=3, borderaxespad=0.0, fontsize=35, diff --git a/src/spacec/plotting/_segmentation.py b/src/spacec/plotting/_segmentation.py index 5c96a39..f3eed62 100644 --- a/src/spacec/plotting/_segmentation.py +++ b/src/spacec/plotting/_segmentation.py @@ -13,7 +13,9 @@ def segmentation_ch( file_name, # image for segmentation channel_file, # all channels used for staining - output_dir, # + output_dir, + savefig=False, # new + output_fname="", # new extra_seg_ch_list=None, # channels used for membrane segmentation nuclei_channel="DAPI", input_format="Multichannel", # CODEX or Phenocycler --> This depends on the machine you are using and the resulting file format (see documentation above) @@ -68,7 +70,18 @@ def segmentation_ch( ax[1].imshow(image_dict["segmentation_channel"]) ax[0].set_title("nuclei") ax[1].set_title("membrane") - plt.show() + + # save or plot figure + if savefig: + plt.savefig( + output_dir + output_fname + ".pdf", + format="pdf", + dpi=300, + transparent=True, + bbox_inches="tight", + ) + else: + plt.show() def show_masks( diff --git a/src/spacec/tools/_segmentation.py b/src/spacec/tools/_segmentation.py index c1b8cf3..0fa3ea9 100644 --- a/src/spacec/tools/_segmentation.py +++ b/src/spacec/tools/_segmentation.py @@ -9,6 +9,7 @@ import pandas as pd import requests import skimage.io +import tensorflow as tf from cellpose import io, models from deepcell.applications import Mesmer from deepcell.utils.plot_utils import create_rgb_image, make_outline_overlay @@ -43,6 +44,7 @@ def cell_segmentation( model_path="./models", resize_factor=1, custom_model=False, + differentiate_nucleus_cytoplasm=False, # experimental ): """ Perform cell segmentation on an image. @@ -82,11 +84,21 @@ def cell_segmentation( Whether to save the segmentation mask as a PNG file. Default is False. model_path : str, optional The path to the model. Default is './models'. + differentiate_nucleus_cytoplasm : bool, optional + Whether to differentiate between nucleus and cytoplasm. Default is False. Returns ------- dict A dictionary containing the original image ('img'), the segmentation masks ('masks'), and the image dictionary ('image_dict'). """ + if use_gpu == True: + gpus = tf.config.list_physical_devices("GPU") + + if gpus: + print("GPU(s) available") + for gpu in gpus: + tf.config.experimental.set_memory_growth(gpu, True) + print("Create image channels!") # check input format @@ -139,72 +151,290 @@ def cell_segmentation( ) # Replace the original image with the resized image in the dictionary segmentation_image_dict[channel] = resized_img - if seg_method == "mesmer": - print("Segmenting with Mesmer!") - if membrane_channel_list is None: + + if differentiate_nucleus_cytoplasm == True: + if membrane_channel_list == None: print( - "Mesmer expects two-channel images as input, where the first channel must be a nuclear channel (e.g. DAPI) and the second channel must be a membrane or cytoplasmic channel (e.g. E-Cadherin)." + "Provide membrane channel for differentiation between nucleus and cytoplasm" ) - sys.exit("Please provide any membrane or cytoplasm channel!") + return else: - masks = mesmer_segmentation( - nuclei_image=segmentation_image_dict[nuclei_channel], - membrane_image=segmentation_image_dict["segmentation_channel"], - plot_predictions=plot_predictions, # plot segmentation results - compartment=compartment, - model_path=model_path, - ) # segment whole cells or nuclei only - else: - print("Segmenting with Cellpose!") - if membrane_channel_list is None: - masks, flows, styles, input_image, rgb_channels = cellpose_segmentation( - image_dict=segmentation_image_dict, - output_dir=output_dir, - membrane_channel=None, - cytoplasm_channel=cytoplasm_channel_list, - nucleus_channel=nuclei_channel, - use_gpu=use_gpu, - model=model, - custom_model=custom_model, - diameter=diameter, - save_mask_as_png=save_mask_as_png, + if seg_method == "mesmer": + print("Segmenting with Mesmer!") + + masks_nuclei = mesmer_segmentation( + nuclei_image=segmentation_image_dict[nuclei_channel], + membrane_image=None, + plot_predictions=plot_predictions, # plot segmentation results + compartment="nuclear", + model_path=model_path, + ) # segment whole cells or nuclei only + + masks_whole_cell = mesmer_segmentation( + nuclei_image=segmentation_image_dict[nuclei_channel], + membrane_image=segmentation_image_dict["segmentation_channel"], + plot_predictions=plot_predictions, # plot segmentation results + compartment=compartment, + model_path=model_path, + ) # segment whole cells or nuclei only + else: + print("Segmenting with Cellpose!") + + ( + masks_nuclei, + flows, + styles, + input_image, + rgb_channels, + ) = cellpose_segmentation( + image_dict=segmentation_image_dict, + output_dir=output_dir, + membrane_channel=None, + cytoplasm_channel=cytoplasm_channel_list, + nucleus_channel=nuclei_channel, + use_gpu=use_gpu, + model=model, + custom_model=custom_model, + diameter=diameter, + save_mask_as_png=save_mask_as_png, + ) + + ( + masks_whole_cell, + flows, + styles, + input_image, + rgb_channels, + ) = cellpose_segmentation( + image_dict=segmentation_image_dict, + output_dir=output_dir, + membrane_channel="segmentation_channel", + cytoplasm_channel=cytoplasm_channel_list, + nucleus_channel=nuclei_channel, + use_gpu=use_gpu, + model=model, + custom_model=custom_model, + diameter=diameter, + save_mask_as_png=save_mask_as_png, + ) + + # Remove single-dimensional entries from the shape of segmentation_masks + masks_whole_cell = masks_whole_cell.squeeze() + # Get the original dimensions of any one of the images + original_height, original_width = image_dict[ + nuclei_channel + ].shape # or any other channel + # Resize the masks back to the original size + masks_whole_cell = cv2.resize( + masks_whole_cell, + (original_width, original_height), + interpolation=cv2.INTER_NEAREST, ) - else: - masks, flows, styles, input_image, rgb_channels = cellpose_segmentation( - image_dict=segmentation_image_dict, - output_dir=output_dir, - membrane_channel="segmentation_channel", - cytoplasm_channel=cytoplasm_channel_list, - nucleus_channel=nuclei_channel, - use_gpu=use_gpu, - model=model, - custom_model=custom_model, - diameter=diameter, - save_mask_as_png=save_mask_as_png, + + # Remove single-dimensional entries from the shape of segmentation_masks + masks_nuclei = masks_nuclei.squeeze() + # Get the original dimensions of any one of the images + original_height, original_width = image_dict[ + nuclei_channel + ].shape # or any other channel + # Resize the masks back to the original size + masks_nuclei = cv2.resize( + masks_nuclei, + (original_width, original_height), + interpolation=cv2.INTER_NEAREST, + ) + + # Create binary masks + binary_masks_nuclei = masks_nuclei > 0 + binary_masks_whole_cell = masks_whole_cell > 0 + + # Subtract the binary nuclei mask from the binary whole cell mask + binary_masks_cytoplasm = binary_masks_whole_cell & ~binary_masks_nuclei + + # Now, if you want to get a labeled mask for the cytoplasm, you can use a function like `label` from `scipy.ndimage` + from scipy.ndimage import label + + masks_cytoplasm, num_labels = label(binary_masks_cytoplasm) + + print("Quantifying features after segmentation!") + print("Quantifying features nuclei") + nuc = extract_features( + image_dict=image_dict, # image dictionary + segmentation_masks=masks_nuclei, # segmentation masks generated by cellpose + channels_to_quantify=channel_names, # list of channels to quantify (here: all channels) + output_file=pathlib.Path(output_dir) + / ( + output_fname + "_" + seg_method + "_nuclei_result.csv" + ), # output path to store results as csv + size_cutoff=size_cutoff, + ) # size cutoff for segmentation masks (default = 0) + print("Quantifying features cytoplasm") + cyto = extract_features( + image_dict=image_dict, # image dictionary + segmentation_masks=masks_cytoplasm, # segmentation masks generated by cellpose + channels_to_quantify=channel_names, # list of channels to quantify (here: all channels) + output_file=pathlib.Path(output_dir) + / ( + output_fname + "_" + seg_method + "_cytoplasm_result.csv" + ), # output path to store results as csv + size_cutoff=size_cutoff, + ) # size cutoff for segmentation masks (default = 0) + + print("Quantifying features whole cell") + whole = extract_features( + image_dict=image_dict, # image dictionary + segmentation_masks=masks_whole_cell, # segmentation masks generated by cellpose + channels_to_quantify=channel_names, # list of channels to quantify (here: all channels) + output_file=pathlib.Path(output_dir) + / ( + output_fname + "_" + seg_method + "_whole_cell_result.csv" + ), # output path to store results as csv + size_cutoff=size_cutoff, + ) # size cutoff for segmentation masks (default = 0) + print("Done!") + + # remove + out = [ + "x", + "y", + "eccentricity", + "perimeter", + "convex_area", + "area", + "axis_major_length", + "axis_minor_length", + "label", + ] + + # keep metadata + whole_meta = whole[out] + + # remove from nuc + nuc = nuc.drop(out, axis=1) + # add whole metadata to cyto + nuc_save = pd.concat([nuc, whole_meta], axis=1) + nuc_save.to_csv( + output_dir + + output_fname + + "_" + + seg_method + + "_nuclei_intensities_result.csv" + ) + # remove from cyto + cyto = cyto.drop(out, axis=1) + # add whole metadata to cyto + cyto_save = pd.concat([cyto, whole_meta], axis=1) + cyto_save.to_csv( + output_dir + + output_fname + + "_" + + seg_method + + "_cytoplasm_intensities_result.csv" ) - # Remove single-dimensional entries from the shape of segmentation_masks - masks = masks.squeeze() - # Get the original dimensions of any one of the images - original_height, original_width = image_dict[ - nuclei_channel - ].shape # or any other channel - # Resize the masks back to the original size - masks = cv2.resize( - masks, (original_width, original_height), interpolation=cv2.INTER_NEAREST - ) - print("Quantifying features after segmentation!") - extract_features( - image_dict=image_dict, # image dictionary - segmentation_masks=masks, # segmentation masks generated by cellpose - channels_to_quantify=channel_names, # list of channels to quantify (here: all channels) - output_file=pathlib.Path(output_dir) - / ( - output_fname + "_" + seg_method + "_result.csv" - ), # output path to store results as csv - size_cutoff=size_cutoff, - ) # size cutoff for segmentation masks (default = 0) - print("Done!") - return {"img": img, "masks": masks, "image_dict": image_dict} + + whole.to_csv( + output_dir + + output_fname + + "_" + + seg_method + + "_whole_cell_intensities_result.csv" + ) + whole = whole.drop(out, axis=1) + + # add identifier to each column name + nuc.columns = [str(col) + "_nuc" for col in nuc.columns] + cyto.columns = [str(col) + "_cyto" for col in cyto.columns] + whole.columns = [str(col) + "_whole" for col in whole.columns] + + # combine the dataframes and save as csv + result = pd.concat([nuc, cyto, whole, whole_meta], axis=1) + result = result.loc[:, ~result.columns.str.contains("Unnamed: 0")] + result.to_csv( + output_dir + + output_fname + + "_" + + seg_method + + "_segmentation_results_combined.csv" + ) + + return { + "img": img, + "masks": masks_whole_cell, + "image_dict": image_dict, + "masks_cytoplasm": masks_cytoplasm, + "masks_nuclei": masks_nuclei, + } + + else: + if seg_method == "mesmer": + print("Segmenting with Mesmer!") + if membrane_channel_list is None: + masks = mesmer_segmentation( + nuclei_image=segmentation_image_dict[nuclei_channel], + membrane_image=None, + plot_predictions=plot_predictions, # plot segmentation results + compartment="nuclear", + model_path=model_path, + ) # segment whole cells or nuclei only + else: + masks = mesmer_segmentation( + nuclei_image=segmentation_image_dict[nuclei_channel], + membrane_image=segmentation_image_dict["segmentation_channel"], + plot_predictions=plot_predictions, # plot segmentation results + compartment=compartment, + model_path=model_path, + ) # segment whole cells or nuclei only + else: + print("Segmenting with Cellpose!") + if membrane_channel_list is None: + masks, flows, styles, input_image, rgb_channels = cellpose_segmentation( + image_dict=segmentation_image_dict, + output_dir=output_dir, + membrane_channel=None, + cytoplasm_channel=cytoplasm_channel_list, + nucleus_channel=nuclei_channel, + use_gpu=use_gpu, + model=model, + custom_model=custom_model, + diameter=diameter, + save_mask_as_png=save_mask_as_png, + ) + else: + masks, flows, styles, input_image, rgb_channels = cellpose_segmentation( + image_dict=segmentation_image_dict, + output_dir=output_dir, + membrane_channel="segmentation_channel", + cytoplasm_channel=cytoplasm_channel_list, + nucleus_channel=nuclei_channel, + use_gpu=use_gpu, + model=model, + custom_model=custom_model, + diameter=diameter, + save_mask_as_png=save_mask_as_png, + ) + # Remove single-dimensional entries from the shape of segmentation_masks + masks = masks.squeeze() + # Get the original dimensions of any one of the images + original_height, original_width = image_dict[ + nuclei_channel + ].shape # or any other channel + # Resize the masks back to the original size + masks = cv2.resize( + masks, (original_width, original_height), interpolation=cv2.INTER_NEAREST + ) + print("Quantifying features after segmentation!") + extract_features( + image_dict=image_dict, # image dictionary + segmentation_masks=masks, # segmentation masks generated by cellpose + channels_to_quantify=channel_names, # list of channels to quantify (here: all channels) + output_file=pathlib.Path(output_dir) + / ( + output_fname + "_" + seg_method + "_result.csv" + ), # output path to store results as csv + size_cutoff=size_cutoff, + ) # size cutoff for segmentation masks (default = 0) + print("Done!") + return {"img": img, "masks": masks, "image_dict": image_dict} def extract_features( @@ -292,6 +522,8 @@ def extract_features( # Export to CSV markers.to_csv(output_file) + return markers + def cellpose_segmentation( image_dict, @@ -622,7 +854,13 @@ def mesmer_segmentation( # Create a combined image stack # Assumes nuclei_image and membrane_image are numpy arrays of the same shape - combined_image = np.stack([nuclei_image, membrane_image], axis=-1) + if membrane_image is None: + # generate empty membrane image + print("No membrane image provided. Nuclear segmentation only.") + membrane_image = np.zeros_like(nuclei_image) + combined_image = np.stack([nuclei_image, membrane_image], axis=-1) + else: + combined_image = np.stack([nuclei_image, membrane_image], axis=-1) # Add an extra dimension to make it compatible with Mesmer's input requirements # Changes shape from (height, width, channels) to (1, height, width, channels)