<jupyter_start><jupyter_text>6. Functions In the previous chapter we developped a small procedure to segment our image of nuclei. If you develop such a routine you are going to re-use it multiple times, so it makes sense to package it into a re-usable unit.We will summarize here how to achieve that in this brief chapter.<jupyter_code>#importing packages
import numpy as np
import matplotlib.pyplot as plt
from skimage.external.tifffile import TiffFile

import course_functions
datapath = course_functions.define_data_path()

import skimage.morphology as skm
import skimage.filters as skf
#load the image to process
data = TiffFile(datapath+'Data/30567/30567.tif')
image = data.pages[3].asarray()
plt.imshow(image,cmap = 'gray')
plt.show()<jupyter_output><empty_output><jupyter_text>Let us summarize all the necessary steps within one code block<jupyter_code>#median filter
image_med = skf.rank.median(image,selem=np.ones((2,2)))
#otsu thresholding
image_local_threshold = skf.threshold_local(image_med,block_size=51)
image_local = image > image_local_threshold
#remove tiny features
image_local_eroded = skm.binary_erosion(image_local, selem= skm.disk(1))
#label image
image_labeled = label(image_local_eroded)
#analyze regions
our_regions = regionprops(image_labeled)
#create a new mask with constraints on the regions to keep
newimage = np.zeros(image.shape)
#fill in using region coordinates
for x in our_regions:
    if (x.area>200) and (x.eccentricity<0.8):
        newimage[x.coords[:,0],x.coords[:,1]] = 1<jupyter_output>/Users/gw18g940/miniconda3/envs/teaching_test/lib/python3.6/site-packages/skimage/filters/rank/generic.py:102: UserWarning: Bitdepth of 14 may result in bad rank filter performance due to large number of bins.
  "performance due to large number of bins." % bitdepth)<jupyter_text>We see that we have an error message. We didn't import the label() function from skimage. Let's add it to the code, as well as the regionprops() function that is missing as well.<jupyter_code>from skimage.measure import label, regionprops

#median filter
image_med = skf.rank.median(image,selem=np.ones((2,2)))
#otsu thresholding
image_local_threshold = skf.threshold_local(image_med,block_size=51)
image_local = image > image_local_threshold
#remove tiny features
image_local_eroded = skm.binary_erosion(image_local, selem= skm.disk(1))
#label image
image_labeled = label(image_local_eroded)
#analyze regions
our_regions = regionprops(image_labeled)
#create a new mask with constraints on the regions to keep
newimage = np.zeros(image.shape)
#fill in using region coordinates
for x in our_regions:
    if (x.area>200):# and (x.eccentricity<0.8):
        newimage[x.coords[:,0],x.coords[:,1]] = 1
plt.figure(figsize=(10,10))
plt.imshow(newimage, cmap = 'gray')
plt.show()<jupyter_output><empty_output><jupyter_text>We can now make a function out of it. You can choose the "level" of your function depending on your needs. For example you could pass a filename and a plane index to the function and make it import your data, or you can pass directly an image.In addition to the image, you coud pass other arguments if you want to make your function more general. For example, you might not always want to filter objects of the same size or shape, and so you can set those as parameters:<jupyter_code>from skimage.measure import label, regionprops

def detect_nuclei(image, size = 200, shape = 0.8):
    #median filter
    image_med = skf.rank.median(image,selem=np.ones((2,2)))
    #otsu thresholding
    image_local_threshold = skf.threshold_local(image_med,block_size=51)
    image_local = image > image_local_threshold
    #remove tiny features
    image_local_eroded = skm.binary_erosion(image_local, selem= skm.disk(1))
    #label image
    image_labeled = label(image_local_eroded)
    #analyze regions
    our_regions = regionprops(image_labeled)
    #create a new mask with constraints on the regions to keep
    newimage = np.zeros(image.shape)
    #fill in using region coordinates
    for x in our_regions:
        if (x.area>size) and (x.eccentricity<shape):
            newimage[x.coords[:,0],x.coords[:,1]] = 1
            
    return newimage<jupyter_output><empty_output><jupyter_text>And now we can test the function (which appears also now in autocompletion):<jupyter_code>nuclei = detect_nuclei(image, size = 400)
plt.imshow(nuclei, cmap = 'gray')
plt.show()<jupyter_output><empty_output><jupyter_text>In order to avoid cluttering your notebooks with function definitions and to be able to reuse your functions across multiple notebooks, I also strongly advice you to create your own module files. Those are .py files that group multipe functions and that can be called from any notebook.Let's create one, call it my_module.py and copy our function in it. Now we can use the function like this:<jupyter_code>import my_module
#or alternatively: from my_module import detect_nuclei
my_module.detect_nuclei()
nuclei2 = my_module.detect_nuclei(image)<jupyter_output>/Users/gw18g940/miniconda3/envs/teaching_test/lib/python3.6/site-packages/skimage/filters/rank/generic.py:102: UserWarning: Bitdepth of 14 may result in bad rank filter performance due to large number of bins.
  "performance due to large number of bins." % bitdepth)
/Users/gw18g940/miniconda3/envs/teaching_test/lib/python3.6/site-packages/skimage/measure/_regionprops.py:250: UserWarning: regionprops and image moments (including moments, normalized moments, central moments, and inertia tensor) of 2D images will change from xy coordinates to rc coordinates in version 0.16.
See http://scikit-image.org/docs/0.14.x/release_notes_and_installation.html#deprecations for details on how to avoid this message.
  warn(XY_TO_RC_DEPRECATION_MESSAGE)
/Users/gw18g940/miniconda3/envs/teaching_test/lib/python3.6/site-packages/skimage/measure/_regionprops.py:260: UserWarning: regionprops and image moments (including moments, normalized moments, central moments, and inertia tensor) of 2D images will chang[...]<jupyter_text>We get an error because in that module, we use skimage functions that were not imported **in the module itself**. We have them in the notebook, but they are not accessible from there.<jupyter_code>plt.imshow(nuclei2, cmap = 'gray')
plt.show()<jupyter_output><empty_output>