Increase code coverage of tests in algorithms module

[tathey1]

-test files should go in the tests folder -submodules with a few short python files can have their tests consolidated into a single test file -goal coverage is 90%

[FelixChihTing]

Hi Matt @mfigdore @tathey1 @bvarjavand , I'm testing connected_threshold in adaptive_thresh.py with a simple image as input: img = np.array([[[100, 250], [800, 300]], [[1200, 2000], [3000, 2500]]]) and a seed (1,0,1). I got the following error: TypeError: in method 'BinaryDilate', argument 2 of type 'std::vector< unsigned int,std::allocator< unsigned int > >' Do you know where the problem is?

[mfigdore]

Hi Matt @mfigdore , I'm testing connected_threshold in adaptive_thresh.py with a simple image as input: img = np.array([[[100, 250], [800, 300]], [[1200, 2000], [3000, 2500]]]) and a seed (1,0,1). I got the following error: TypeError: in method 'BinaryDilate', argument 2 of type 'std::vector< unsigned int,std::allocator< unsigned int > >' Do you know where the problem is?

Looks like the error is coming from SimpleITK BinaryDilate at the C++ level, from line 254 of https://github.com/neurodata/brainlit/blob/develop/brainlit/algorithms/generate_fragments/adaptive_thresh.py

Does the tutorial work? https://brainlit.netlify.app/notebooks/algorithms/tutorial_notebook_adaptive_thresh.html

[mfigdore]

Running the code in the function manually with the img and seed you provided, it seems to work for me (although seed might need to be a list [(1,0,1)] -- is that an error in the documentation? (From line 251 being replaced by 252-3 to support multiple seeds) @FelixChihTing

[FelixChihTing]

@mfigdore I got the same error after changing from (1,0,1) to [(1,0,1)].

[mfigdore]

@FelixChihTing Does the tutorial work? Also, I suggest you check with others to see if it's a problem related to your OS or simpleITK version

[FelixChihTing]

@mfigdore the tutorial is out of date, some files are moved, e.g., NeuroglancerSession is now in session.py

[FelixChihTing]

@mfigdore currently only the following block in the tutorial is working: from skimage.data import astronaut # create the viewer and display the image viewer = napari.view_image(astronaut(), rgb=True)

[FelixChihTing]

I believe it's an issue of SimpleITK.

[mfigdore]

I believe it's an issue of SimpleITK.

@FelixChihTing Does the solution proposed of casting to 'int' fix it?

[mfigdore]

@mfigdore currently only the following block in the tutorial is working:

from skimage.data import astronaut

# create the viewer and display the image

viewer = napari.view_image(astronaut(), rgb=True)

@bvarjavand is someone keeping tutorials up-to-date?

[FelixChihTing]

@mfigdore the proposed solution does not fix the problem in my case

[mfigdore]

@FelixChihTing Please see https://github.com/neurodata/brainlit/pull/226 which I believe will address the error.

[FelixChihTing]

@FelixChihTing Please see #226 which I believe will address the error.

@mfigdore Great! that resolve the problem with sitk.BinaryDilate. Other than that, when I run connected_threshold with the following input:

img = np.array([[[100, 250], [800, 300]], [[1200, 2000], [3000, 2500]]])
labels = connected_threshold(img,([1,0,1]))

I got similar error:

TypeError: in method 'ConnectedThreshold', argument 2 of type 'std::vector< std::vector< uint32_t,std::allocator< uint32_t > >,std::allocator< std::vector< uint32_t,std::allocator< uint32_t > > > >'

so in this case, elements in seedList are not recognized as int

[mfigdore]

@FelixChihTing Should be labels = connected_threshold(img, [(1,0,1)] )

[FelixChihTing]

@FelixChihTing Should be labels = connected_threshold(img, [(1,0,1)] )

Ohh, this works!! Thank you!

[FelixChihTing]

Hi @mfigdore @tathey1 , how is the sigma in fast_marching_seg used in calculating the labels? I tried to input a 2D Gaussian image with standard deviation of 3 (left), and the result is shown at the right. I don't understand how the filter set the threshold in this case.

from brainlit.algorithms.generate_fragments.adaptive_thresh import (
    get_seed,
    get_img_T1,
    thres_from_gmm,
    fast_marching_seg,
    level_set_seg,
    connected_threshold,
    confidence_connected_threshold,
    neighborhood_connected_threshold,
    otsu,
    gmm_seg,
)
import SimpleITK as sitk
import numpy as np
import matplotlib.pyplot as plt

Gx = np.array([])
for x in range(0,101):
    Gx = np.insert(Gx,x,np.exp(-((x-50)**2)/(2*(5**2))))

Gxy = np.zeros((101,101))

for x in range(0,101):
    for y in range(0,101):
        Gxy[x,y] = (Gx[x])*(Gx[y])

img = 255*(Gxy/Gxy.max())
lower_threshold = 200
seed = (50,50)
labels = fast_marching_seg(img, seed,sigma=0.5)
plt.subplot(221)
plt.imshow(img)
plt.subplot(222)
plt.imshow(labels)
plt.subplot(223)
plt.plot(Gxy[50,:])
plt.subplot(224)
plt.plot(labels[50,:])
plt.show()```

[mfigdore]

@FelixChihTing Please refer to SimpleITK documentation. https://simpleitk.readthedocs.io/en/master/link_FastMarchingSegmentation_docs.html

A Gaussian with a parameterized sigma is used during the gradient computation to enable the level-set to slow down near edges.

The sigma parameter is used by SimpleITK GradientMagnitudeRecursiveGaussian. See: https://itk.org/ITKExamples/src/Filtering/ImageGradient/ComputeGradientMagnitudeRecursiveGaussian/Documentation.html https://itk.org/Doxygen/html/classitk_1_1GradientMagnitudeRecursiveGaussianImageFilter.html

[FelixChihTing]

@tathey1 When I tried to segment a simple image as below, level_set_seg doesn't work correctly:

img = np.zeros((15,15))
for x in range(0,16):
    for y in range(0,16):
        if x > 2 and x < 5 and y > 2 and y < 5:
            img[x,y] = 200
        if x > 6 and x < 11 and y > 6 and y < 11:
            img[x,y] = 255

seed = (7,7)
labels = level_set_seg(img,seed)
plt.subplot(121)
plt.imshow(img)
plt.title('input')
plt.subplot(122)
plt.imshow(labels)
plt.title('output, seed at (7,7)')
plt.show()

But if I put a Gaussian image at the background, the result will become:

Gx = np.array([])
for x in range(0,16):
    Gx = np.insert(Gx,x,np.exp(-((x-7.5)**2)/(2*(2**2))))

Gxy = np.zeros([16,16])
for x in range(0,16):
    for y in range(0,16):
        Gxy[x,y] = (Gx[x])*(Gx[y])

img = 200*(Gxy/Gxy.max())

for x in range(0,16):
    for y in range(0,16):
        if x > 2 and x < 5 and y > 2 and y < 5:
            img[x,y] = 200
        if x > 6 and x < 11 and y > 6 and y < 11:
            img[x,y] = 255

seed = (7,7)
labels = level_set_seg(img,seed)
plt.subplot(121)
plt.imshow(img)
plt.title('input')
plt.subplot(122)
plt.imshow(labels)
plt.title('output, seed at (7,7)')
plt.show()

[FelixChihTing]

@tathey1 @mfigdore @bvarjavand , could you suggest any document or reading explaining the math and application of fast marching and level set method? thanks

[FelixChihTing]

@tathey1 When I run the notebook, I got the following error ClientError: An error occurred (InvalidAccessKeyId) when calling the GetObject operation: The AWS Access Key Id you provided does not exist in our records.

[FelixChihTing]

Hi @mfigdore, i'm updating the adaptive thresh notebook, may I ask how do you put the seeds to each case to ensure desirable results? In one case, I saw _, seed = adaptive_thresh.get_seed(voxel) to get the seed, but such a seed gives the following error:IndexError: index 60 is outside the extent for dimension 0 with size 52

[mfigdore]

Hi @mfigdore, i'm updating the adaptive thresh notebook, may I ask how do you put the seeds to each case to ensure desirable results? In one case, I saw _, seed = adaptive_thresh.get_seed(voxel) to get the seed, but such a seed gives the following error:IndexError: index 60 is outside the extent for dimension 0 with size 52

The seeds are based on manual annotations. The (x,y,z) coordinates need to be in the correct order so nothing is out of range. @tathey1 or @bvarjavand can you please point Felix to these annotations?