<jupyter_start><jupyter_text>Nuclear Morphology and Chromatin Organization FeaturesHere we aim to compute a library of features that exhaustively describe the nuclear morphology and chromatin organization for each segmented nucleus in a given image.<jupyter_code># import libraries
%load_ext autoreload
import sys
sys.path.append("../")

from tifffile import imread
import pandas as pd 
from skimage import measure
import numpy as np
import matplotlib.pyplot as plt
import cv2 as cv
import src.nuclear_features.Boundary_local_curvature as BLC
import src.nuclear_features.Boundary_global as BG
import src.nuclear_features.Int_dist_features as IDF
import src.nuclear_features.Img_texture as IT
import os

# initialising paths
labelled_image_path = os.path.join(os.path.dirname(os.getcwd()),'example_data/images/TMA_nuc_labels.tif')
raw_image_path = os.path.join(os.path.dirname(os.getcwd()),'example_data/images/TMA_DAPI.tif')

feature_path = os.path.join(os.path.dirname(os.getcwd()),'example_data/')<jupyter_output><empty_output><jupyter_text>Below is an example of the data that can be used.<jupyter_code>#Read in Images
labelled_image = imread(labelled_image_path)
raw_image = imread(raw_image_path)

#Subset the image
labelled_image=labelled_image[4000:6000,5000:7000]
raw_image=raw_image[4000:6000,5000:7000]

# normalize images
raw_image = ((raw_image-np.min(raw_image))/(np.max(raw_image)-np.min(raw_image)))*255
raw_image = raw_image.astype(int)
#Visulaise the data

#save plots to show clusters
fig = plt.figure(figsize=(8, 4))
ax0 = fig.add_subplot(121)
ax1 = fig.add_subplot(122)
#show raw image 
ax0.imshow(raw_image,aspect='auto',cmap='inferno') 
ax0.axis('off')
ax0.title.set_text('Image')
#show segmented image
ax1.imshow(labelled_image,aspect='auto',cmap='viridis') 
ax1.axis('off')
ax1.title.set_text('Nuclear Labels')<jupyter_output><empty_output><jupyter_text>One can now access each nucleus in the labelled image as well as the raw image.<jupyter_code>#Get indexing for the individual nuclei in the image
props = measure.regionprops(labelled_image,raw_image)

fig = plt.figure(figsize=(12, 6))
ax0 = fig.add_subplot(241)
ax1 = fig.add_subplot(242)
ax2 = fig.add_subplot(243)
ax3 = fig.add_subplot(244)
ax4 = fig.add_subplot(245)
ax5 = fig.add_subplot(246)
ax6 = fig.add_subplot(247)
ax7 = fig.add_subplot(248)

#show raw image 
ax0.imshow(props[77].intensity_image,aspect='auto',cmap='inferno') 
ax0.title.set_text('Nucleus 1')
ax0.axis('off')
ax1.imshow(props[452].intensity_image,aspect='auto',cmap='inferno') 
ax1.title.set_text('Nucleus 2')
ax1.axis('off')
ax2.imshow(props[567].intensity_image,aspect='auto',cmap='inferno') 
ax2.title.set_text('Nucleus 3')
ax2.axis('off')
ax3.imshow(props[114].intensity_image,aspect='auto',cmap='inferno') 
ax3.title.set_text('Nucleus 4')
ax3.axis('off')

#show segmented image
ax4.imshow(props[77].image,aspect='auto',cmap='viridis') 
ax4.title.set_text('Label 1')
ax4.axis('off')
ax5.imshow(props[452].image,aspect='auto',cmap='viridis') 
ax5.title.set_text('Label 2')
ax5.axis('off')
ax6.imshow(props[567].image,aspect='auto',cmap='viridis') 
ax6.title.set_text('Label 3')
ax6.axis('off')
ax7.imshow(props[114].image,aspect='auto',cmap='viridis') 
ax7.title.set_text('Label 4')
ax7.axis('off')<jupyter_output><empty_output><jupyter_text>Basic FeaturesScikit provides several informative features that describe "region properties". First we extract such built in features. For more information on how thhe features were computed check out the documentation(https://scikit-image.org/docs/dev/api/skimage.measure.htmlskimage.measure.regionprops_table).<jupyter_code>#Measure scikit's built in features
propstable = pd.DataFrame(measure.regionprops_table(labelled_image,raw_image,cache=True,
                               properties=['label', 'area','perimeter','bbox_area','convex_area',
                                           'equivalent_diameter','major_axis_length','minor_axis_length',
                                           'eccentricity','orientation',
                                            'centroid','weighted_centroid',
                                           'weighted_moments','weighted_moments_normalized',
                                           'weighted_moments_central','weighted_moments_hu',
                                            'moments','moments_normalized','moments_central','moments_hu']))
propstable.iloc[[77,452,567,114]]<jupyter_output><empty_output><jupyter_text>Global Boundary featuresHere we compute the features that describe morphology of a given object. These include 1. Calliper distances2. Distribution features of radii(centroid to boundary) distancesBelow are features computed for 4 nuclei<jupyter_code>BG_feat = pd.concat([BG.boundary_features(props[77].image,centroids=props[77].local_centroid),
           BG.boundary_features(props[452].image,centroids=props[452].local_centroid),
           BG.boundary_features(props[567].image,centroids=props[567].local_centroid),
           BG.boundary_features(props[114].image,centroids=props[114].local_centroid)])
BG_feat<jupyter_output><empty_output><jupyter_text>Local Boundary FeaturesHere we compute the features that describe local curvature of a given object.Approach:For a given object we obtain the edge pixels and compute the local curvature of each point on the curve +/- a given stepsize. Larger steps give a smoother curvature. We define the local curvature between 3 points as the inverse of the radius of their circumcircle and if the circumcenter is inside the object then the sign of curvature is positive. Below is the radius of curvature for Nucleus 1<jupyter_code>r_c= BLC.local_radius_curvature(props[77].image,step=5,show_boundary=True)
#calculate local curvature features
local_curvature=[np.divide(1,r_c[x]) if r_c[x]!=0 else 0 for x in range(len(r_c))]<jupyter_output><empty_output><jupyter_text>Now that we have the local curvature for all points on the boundary, we compute features that describe it such as Average and Standard Deviation of curature (positive and negative), number of time the polarity changes etc. Feature names are self-descriptive.<jupyter_code>#compute local and global features
global_features = [BLC.global_curvature_features(np.array(local_curvature))]
global_features = pd.DataFrame([o.__dict__ for o in global_features])
global_features<jupyter_output><empty_output><jupyter_text>We also check to see if there are any prominant jumps in curvature the image.<jupyter_code>prominant_features = [BLC.prominant_curvature_features(local_curvature,show_plot=True)]
prominant_features = pd.DataFrame([o.__dict__ for o in prominant_features])
prominant_features<jupyter_output><empty_output><jupyter_text>Below are the features computed for 4 nuclei.<jupyter_code>BLC_feat= pd.concat([BLC.curvature_features(props[77].image,step=5),
           BLC.curvature_features(props[452].image,step=5),
           BLC.curvature_features(props[567].image,step=5),
           BLC.curvature_features(props[114].image,step=5)])
BLC_feat<jupyter_output><empty_output><jupyter_text>Intensity FeaturesHere we compute features that describe the intensity distribution.These include features that describe the intensity distribution, entropy and heterocromatin ratios Below are the features computed for 4 nuclei.<jupyter_code>Int_feat= pd.concat([IDF.intensity_features(props[77].image,props[77].intensity_image),
                     IDF.intensity_features(props[452].image,props[452].intensity_image),
                     IDF.intensity_features(props[567].image,props[567].intensity_image),
                     IDF.intensity_features(props[114].image,props[114].intensity_image)])
Int_feat<jupyter_output><empty_output><jupyter_text>Image TexturesHere we compute features that describe the texture of the image.These include the GCLM features. Below are the features computed for 4 nuclei.<jupyter_code>Int_Text= pd.concat([IT.texture_features(props[77].image,props[77].intensity_image,props[77].local_centroid),
                     IT.texture_features(props[452].image,props[452].intensity_image,props[452].local_centroid),
                     IT.texture_features(props[567].image,props[567].intensity_image,props[567].local_centroid),
                     IT.texture_features(props[114].image,props[114].intensity_image,props[114].local_centroid)])
Int_Text<jupyter_output><empty_output><jupyter_text>Misc. featuresWe merge all features and compute some related features.<jupyter_code>features = pd.concat([propstable.iloc[[77,452,567,114]].reset_index(drop=True), 
           BG_feat.reset_index(drop=True), 
           BLC_feat.reset_index(drop=True), 
           Int_feat.reset_index(drop=True),
           Int_Text.reset_index(drop=True)], axis=1)

features['Concavity']=(features['convex_area']-features['area'])/features['convex_area']
features['Solidity']=features['area']/features['convex_area']
features['A_R']=features['minor_axis_length']/features['major_axis_length']
features['Shape_Factor']=(features['perimeter']**2)/(4*np.pi*features['area'])
features['Area_bbArea']=features['area']/features['bbox_area']
features['Center_Mismatch']=np.sqrt((features['weighted_centroid-0']-features['centroid-0'])**2+
                             (features['weighted_centroid-1']-features['centroid-1'])**2)
features['Smallest_largest_Calliper']=features['Min_Calliper']/features['Max_Calliper']
features['Frac_Peri_w_posi_curvature']=features['Len_posi_Curvature']/features['perimeter']
features['Frac_Peri_w_neg_curvature']=features['Len_neg_Curvature'].replace(to_replace ="NA",value =0)/features['perimeter']
features['Frac_Peri_w_polarity_changes']=features['nPolarity_changes']/features['perimeter']

features<jupyter_output><empty_output><jupyter_text>For a quick extraction of all features given a segmented image use the following code:<jupyter_code>from src.utlis.Run_nuclear_feature_extraction import run_nuclear_chromatin_feat_ext

features = run_nuclear_chromatin_feat_ext(raw_image_path,labelled_image_path,feature_path)
features_1 = features.replace('NA',0, regex=True)
features_1 = features_1.replace('NaN',0, regex=True)
features_1<jupyter_output><empty_output><jupyter_text>Tissue level summary:In order to characterise the nuclear density/crowding in a given tissue, we compute the distribution characteristics of each of the above features.The measures available are: Median, Min, Max, Standard Deviation (SD) Coefficient of Variation (CV) and Coefficient of Dispersion (CD), Inter_Quartile_Range(IQR) and Quartile Coeffient of Dispersrion (QCD).<jupyter_code>from src.utlis.summarising_features import summarise_feature_table
features_1 = features.replace('NA',0, regex=True)
features_1 = features_1.replace('NaN',0, regex=True)
summarise_feature_table(features_1.drop(['Image'],axis=1))<jupyter_output>/home/pathy_s/anaconda3/lib/python3.8/site-packages/numpy/lib/nanfunctions.py:995: RuntimeWarning: All-NaN slice encountered
  result = np.apply_along_axis(_nanmedian1d, axis, a, overwrite_input)
/home/pathy_s/anaconda3/lib/python3.8/site-packages/numpy/lib/nanfunctions.py:1391: RuntimeWarning: All-NaN slice encountered
  result = np.apply_along_axis(_nanquantile_1d, axis, a, q,
{'fraction_non_alphanumeric': 0.11587867975022302, 'fraction_numerical': 0.027386262265834075, 'mean_word_length': 6.0910815939278935, 'pattern_counts': {'":': 0, '<': 56, '<?xml version=': 0, '>': 56, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 12, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '14', 'dir': 'jupyter-structured-clean-dedup', 'id': '495469', 'license': 'MIT', 'n_tokens_mistral': 3865, 'n_tokens_neox': 3539, 'n_words': 682, 'path': 'notes_on_feature_extraction/Nuclear_Features.ipynb', 'repo_name': 'GVS-Lab/genomic-scoring-breast-cancer-progression', 'size': '11210'}