Density fig updates

nuc_morph_analysis/analyses/dataset_images_for_figures/figure_helper.py

@@ -98,7 +98,7 @@ def get_save_dir_and_fig_panel_str(figure, panel):
     return savedir, fig_panel_str
 
 
-def return_glasbey_on_dark(N=255, cell_id=None, cell_color=None):
+def return_glasbey_on_dark(N=255, cell_id=None, cell_color=None, from_list=False):
     """
     Publication
     The Glasbey LUT is based on the publication:
@@ -109,8 +109,9 @@ def return_glasbey_on_dark(N=255, cell_id=None, cell_color=None):
     N: int, the number of colors to return
     cell_id: int, the cell id to change the color of
     cell_color: array, the RGB color to change the cell to
+    from_list: bool, whether colormap should be returned from ListedColormap or not
     """
-    from matplotlib.colors import LinearSegmentedColormap
+    from matplotlib.colors import LinearSegmentedColormap, ListedColormap
     from nuc_morph_analysis.analyses.dataset_images_for_figures.glasbey_on_dark import (
         glasbey_on_dark_func,
     )
@@ -136,7 +137,8 @@ def return_glasbey_on_dark(N=255, cell_id=None, cell_color=None):
     cmap_name = "glasbey_on_dark"
     rgb_array0_1 = [tuple(np.asarray(x) / 255) for x in glasbey_on_dark]
     cmap = LinearSegmentedColormap.from_list(cmap_name, rgb_array0_1, N=N + 1)
-
+    if from_list:
+        cmap = ListedColormap(rgb_array0_1, cmap_name, N=N)
     return np.asarray(rgb_array0_255), cmap, rgb_array0_255
 
 

nuc_morph_analysis/analyses/density/pseudocell_based_density_analysis.py → nuc_morph_analysis/analyses/density/extra_checks/check_for_size_dependence.py

@@ -1,4 +1,5 @@
 #%%
+# this code is not used in the final analysis, but was used to confirm that there is no correlation between density and nuclear area
 from nuc_morph_analysis.lib.visualization.reference_points import COLONY_COLORS, COLONY_LABELS
 import matplotlib.pyplot as plt
 import numpy as np
@@ -8,6 +9,9 @@
 from nuc_morph_analysis.lib.preprocessing import filter_data
 from pathlib import Path
 from sklearn.linear_model import LinearRegression
+from matplotlib.ticker import MaxNLocator
+import matplotlib
+
 
 #%%
 # set figure directory
@@ -21,57 +25,23 @@
 dfm = filter_data.all_timepoints_minimal_filtering(df)
 
 #%% # plot density over time for each colony  along colony time 
-x_col = "colony_time"
-column_val = 'label_img'
-feature_list = ['2d_area_nuc_cell_ratio','density','2d_area_nucleus','2d_area_pseudo_cell',
-                '2d_area_cyto','inv_cyto_density','2d_intensity_min_edge','2d_intensity_max_edge','2d_intensity_mean_edge'
-                ]
-for y_col in feature_list:
-    fig,ax = plt.subplots(figsize=(4,3))
-
-    for colony in ['small','medium','large']:
-
-        dfsub = dfm[dfm['colony']==colony].copy()
-        dfsub.dropna(subset=[y_col],inplace=True)
-
-        # create a pivot of the dataframe to get a 2d array of track_id x timepoint with each value being the density
-        pivot = dfsub.pivot(index=x_col, columns=column_val, values=y_col)
-        pivot.head()
+matplotlib.rcParams.update({'font.size': 8})
+matplotlib.rcParams.update({'axes.titlesize': 8})
+matplotlib.rcParams.update({'axes.labelsize': 8})
+matplotlib.rcParams.update({'xtick.labelsize': 8})
+matplotlib.rcParams.update({'ytick.labelsize': 8})
+matplotlib.rcParams.update({'legend.fontsize': 7})
+matplotlib.rcParams.update({'legend.title_fontsize': 8})
+matplotlib.rcParams.update({'figure.titlesize': 8})
+# make axis linewidth thinner
+matplotlib.rcParams.update({'axes.linewidth': 0.5})
 
-        mean = pivot.median(axis=1)
-        lower = pivot.quantile(0.05,axis=1)
-        upper = pivot.quantile(0.95,axis=1)
-
-        xscale_factor, xlabel, xunit, xlimit = get_plot_labels_for_metric(x_col)
-        x = mean.index * xscale_factor
-        yscale_factor, ylabel, yunit, ylimit = get_plot_labels_for_metric(y_col)
-        y = mean.values * yscale_factor
-        yl = lower.values * yscale_factor
-        yu = upper.values * yscale_factor
-
-        ax.plot(x, y, label=COLONY_LABELS[colony], color=COLONY_COLORS[colony])
-        ax.fill_between(x, yl, yu, alpha=0.2, color=COLONY_COLORS[colony],
-                        edgecolor='none')
-        ax.set_xlabel(f"{xlabel} {xunit}")
-        ax.set_ylabel(f"{ylabel} {yunit}\n(90% interpercentile range)")
-
-
-    ax.legend(loc="upper right", handletextpad=0.7, frameon=False)
-    plt.tight_layout()
-    for ext in ['.png','.pdf']:
-        save_and_show_plot(
-            f"{figdir}/{y_col}_vs_{x_col}_by_colony",
-            file_extension=ext,
-            dpi=300,
-            transparent=False,
-        )
-    plt.show()
 
 #%%
-# plot density as a function of nucleus size (and compare to old density metric)
+# plot density as a function of nucleus size 
 colony='medium'
 x_col = '2d_area_nucleus'
-for yi,y_col in enumerate(['2d_area_nuc_cell_ratio','density','inv_cyto_density']):
+for yi,y_col in enumerate(['2d_area_nuc_cell_ratio']):
 
     dfsub = dfm[dfm['colony']==colony].copy()
     dfsub.dropna(subset=[y_col],inplace=True)
@@ -88,6 +58,7 @@
     ax.set_xlabel(f"{xlabel} {xunit}")
     ax.set_ylabel(f"{ylabel} {yunit}")
 
+
     # add best fit line
     reg = LinearRegression().fit(x.reshape(-1,1), y)
     y_pred = reg.predict(x.reshape(-1,1))

nuc_morph_analysis/analyses/density/validate_new_filters.py → nuc_morph_analysis/analyses/density/extra_checks/visualize_filter_for_uncaught_pseudo_cell_artifacts.py

@@ -1,4 +1,7 @@
 #%%
+# this code visualizes cells identified in this workflow `remove_uncaught_pseudo_cell_artifacts(df, apply_to_nucleus_too, verbose)`
+# and highlights the cells that were removed
+
 from nuc_morph_analysis.lib.preprocessing.twoD_zMIP_area import watershed_workflow
 from pathlib import Path
 from nuc_morph_analysis.lib.preprocessing import global_dataset_filtering

nuc_morph_analysis/analyses/density/watershed_validate.py → nuc_morph_analysis/analyses/density/extra_checks/visually_validate_watershed_psuedo_cell_seg_workflow.py

@@ -14,13 +14,13 @@
 
 def get_contours_from_pair_of_2d_seg_image(nuc_mip,cell_mip,dft=None):
         contour_list = [] # (label, nucleus_contour, cell_contour, color)
-        rgb_array0_255, _, _ = return_glasbey_on_dark()
-        for label_img in range(0,nuc_mip.max()+1):
+        rgb_array0_255, cmapper, _ = return_glasbey_on_dark()
+        for label_img in range(1,nuc_mip.max()+1):
             if dft is not None:
                 #ask if the label_img is in the dataframe
                 if label_img not in dft['label_img'].values:
                     continue
-            color = np.float64(rgb_array0_255[label_img % len(rgb_array0_255)])
+            color = np.float64(cmapper(label_img))
 
             # get the nucleus boundary
             nucleus_boundary = nuc_mip == label_img
@@ -36,15 +36,39 @@ def get_contours_from_pair_of_2d_seg_image(nuc_mip,cell_mip,dft=None):
 
         return contour_list
 
-def draw_contours_on_image(axlist,contour_list,new_color=None):
+def draw_contours_on_image(axlist,contour_list,new_color=None,filled=False,colorize=False,dft=None,linewidth=1,colorfeat='2d_area_nuc_cell_ratio',cmapstr='PiYG',usedfcolors=False,cmap=None):
     # draw contours
+    if dft is not None:
+        minval = dft[colorfeat].min()
+        maxval = dft[colorfeat].max()
     for label, nuc_contours, cell_contours, color in contour_list:
         if new_color is not None:
             color = new_color
+
+        if usedfcolors:
+            if label not in dft['label_img'].values:
+                color = np.asarray((0.4,0.4,0.4,1))
+            else:
+                color = np.float64(cmap(dft.loc[dft['label_img']==label,colorfeat].values[0]))
+
+        if colorize:
+            if label not in dft['label_img'].values:
+                continue
+            value = dft[dft['label_img']==label][colorfeat].values[0]
+            #rescale between 0 and 255
+            new_value = (value - minval) / (maxval - minval)
+            cmap = cm.get_cmap(cmapstr)
+            color = np.float64(cmap(new_value))
+
         for contour in nuc_contours:
-            axlist.plot(contour[:, 1], contour[:, 0], linewidth=1, color=color/255)
+            axlist.plot(contour[:, 1], contour[:, 0], linewidth=linewidth, color=color)
+            if filled: # now draw as filled
+                axlist.fill(contour[:, 1], contour[:, 0], color=color, alpha=0.5)
         for contour in cell_contours:
-            axlist.plot(contour[:, 1], contour[:, 0], linewidth=1, color=color/255)
+            axlist.plot(contour[:, 1], contour[:, 0], linewidth=linewidth, color=color)
+            if filled: # now draw as filled
+                axlist.fill(contour[:, 1], contour[:, 0], color=color, alpha=0.5)
+
     return axlist
 
 def plot_colorized_image_with_contours(img_dict,dft,feature,cmapstr,colony='test',TIMEPOINT=None,RESOLUTION_LEVEL=None,categorical=False,draw_contours=True):
@@ -71,8 +95,6 @@ def plot_colorized_image_with_contours(img_dict,dft,feature,cmapstr,colony='test
         cmap = cm.get_cmap(cmapstr)
         if categorical:
             cimg = np.round(cimg).astype('uint16')
-
-        # rgb = np.take(np.uint16(cmaparr*255),cimg.astype('uint16'),axis=0)
 
         # create the figure
         fig, axlist = plt.subplots(1, 1, figsize=(6, 4))
@@ -88,7 +110,7 @@ def plot_colorized_image_with_contours(img_dict,dft,feature,cmapstr,colony='test
         if draw_contours:
             # create the contours
             contour_list = get_contours_from_pair_of_2d_seg_image(nuc_mip,cell_mip,dft)
-            draw_contours_on_image(axlist,contour_list,new_color=np.asarray((200,0,200)))
+            draw_contours_on_image(axlist,contour_list,new_color=np.asarray((200,0,200))/255)
         # remove the axis
         plt.axis('off')
         plt.xlim([x1,x1+w])
@@ -219,6 +241,5 @@ def run_validation_and_plot(TIMEPOINT=48,colony='medium',RESOLUTION_LEVEL=1,plot
 
 if __name__ == '__main__':
     # set the details
-    # dft_test = run_validation_and_plot(testing=True)
-    # dft0 = run_validation_and_plot(plot_everything=True)
-    dft0 = run_validation_and_plot(247,colony='small',RESOLUTION_LEVEL=1,plot_everything=True)
+    dft_test = run_validation_and_plot(testing=True,plot_everything=False)
+    dft0 = run_validation_and_plot(testing=False,plot_everything=False)