Tidy code and add LICENSE and README.md.

.gitignore

@@ -102,3 +102,6 @@ ENV/
 
 # Ignore featurization folder
 /featurization/
+
+# Ignore visualization folder
+/visualization/

README.md

@@ -0,0 +1,14 @@
+Lung Cancer Action Toolkit (LCAT)
+=================================
+The Lung Cancer Action Toolkit (LCAT) provides tools for loading and analyzing chest CT scans.
+
+The Lung Cancer Action Toolkit includes functionality for:
+* Loading LIDC-IDRI data, including
+  * Voxel data, optionally normalizing pixel distances
+  * Nodule data, including radiologist segmentations and radiologist-assigned characteristics
+* Performing segmentation of the lung from chest CTs
+* Performing segmentation of the body from chest CTs
+* Analyzing mouth-to-point distances through the trachea and the bronchioles of the lungs
+* Analyzing nodule size and shape
+
+Detailed documentation for each of these functions is available at TODO.

lcat/analysis/tracheal_distance.py

@@ -26,11 +26,7 @@ def get_tracheal_distances(scan, lung_segmentation):
     phi = np.ma.MaskedArray(seed_boundary, air_boundaries)
 
     # Perform fast marching problem
-    try:
-        distances = skfmm.distance(phi)
-    except:
-        import IPython
-        IPython.embed()
+    distances = skfmm.distance(phi)
 
     # TODO: Return something else?
     return distances

lcat/featurization/__init__.py

@@ -7,7 +7,9 @@
 from . import registry
 
 # Import featurization modules
+from . import body_depth
 from . import center
+from . import characteristics
 from . import region_properties
 from . import tracheal_distance
 

lcat/featurization/body_depth.py

@@ -0,0 +1,51 @@
+"""
+Body depth featurization module.
+"""
+from __future__ import absolute_import
+
+import numpy as np
+import pandas as pd
+import scipy.ndimage
+
+import lcat
+from . import registry
+
+
+@registry.register_featurizer('body_depth')
+def featurize_center(scan):
+    """
+    Featurize the given scan, returning body depth statistics.
+    """
+    # Create data distance placeholder
+    index = pd.Index([], name='nodule_id')
+    data = pd.DataFrame(index=index, columns=['min_body_depth',
+                                              'mean_body_depth',
+                                              'median_body_depth',
+                                              'max_body_depth'])
+
+    # Get the body segmentation for the scan
+    body_shell = lcat.get_body_segmentation(scan)
+
+    # Get body depths
+    body_depths = scipy.ndimage.distance_transform_edt(body_shell)
+
+    # Multiply by unit lengths (assume cubic unit cell)
+    body_depths *= np.mean(scan.unit_cell[0])
+
+    # For each nodule
+    for nodule in scan.nodules:
+        # Create full mask
+        mask = lcat.util.get_full_nodule_mask(nodule, scan.voxels.shape)
+
+        # Select tumor tracheal distances
+        nodule_depths = body_depths[mask]
+
+        # Add attributes to dataframe
+        data.loc[nodule.nodule_id, :] = [
+            np.min(nodule_depths),
+            np.mean(nodule_depths),
+            np.median(nodule_depths),
+            np.max(nodule_depths)
+        ]
+
+    return data

lcat/featurization/characteristics.py

@@ -0,0 +1,42 @@
+"""
+Characteristics featurization module.
+"""
+from __future__ import absolute_import
+
+import numpy as np
+import pandas as pd
+import scipy.ndimage
+
+import lcat
+from . import registry
+
+
+CHARACTERISTICS = [
+	'subtlety',
+	'internalStructure',
+	'calcification',
+	'sphericity',
+	'margin',
+	'lobulation',
+	'spiculation',
+	'texture',
+	'malignancy',
+]
+
+
+@registry.register_featurizer('characteristics')
+def featurize_characteristics(scan):
+    """
+    Featurize the given scan, returning nodule characteristics.
+    """
+    # Create data distance placeholder
+    index = pd.Index([], name='nodule_id')
+    data = pd.DataFrame(index=index, columns=CHARACTERISTICS)
+
+    # For each nodule
+    for nodule in scan.nodules:
+        # Load characteristics
+        data.loc[nodule.nodule_id] = [nodule.characteristics.get(attribute, np.nan)
+                                      for attribute in CHARACTERISTICS]
+
+    return data

lcat/segmentation/__init__.py

@@ -3,4 +3,5 @@
 """
 from __future__ import absolute_import
 
+from .body import get_body_segmentation
 from .lungs import get_lung_segmentation

lcat/segmentation/body.py

@@ -0,0 +1,72 @@
+"""
+Segments a body from a CT scan.
+"""
+import numpy as np
+import skimage
+import skimage.measure
+import skimage.segmentation
+
+import lcat
+
+
+def get_body_segmentation(scan):
+    """
+    Given a `Scan` object representing a chest CT scan, return a binary mask representing the region
+    occupied by the body.
+    """
+    # Threshold the image (threshold is lower limit for lung tissue in HU)
+    foreground = scan.voxels >= -700
+
+    # Identify strongly connected components
+    labels = skimage.measure.label(foreground, connectivity=1)
+
+    # Identify the largest volume
+    body_mask = get_largest_volume(labels)
+
+    # Obtain body envelope
+    envelope_mask = get_body_envelope(body_mask)
+
+    return envelope_mask
+
+
+def get_largest_volume(labels):
+    """
+    Return a binary mask equivalent to the component with the largest volumes in the array `labels`.
+    """
+    return labels == get_top_value(labels[labels != 0])
+
+
+def get_top_value(arr):
+    """
+    Given an ndarray, return the value which occurs most frequently.
+    """
+    # Count values
+    values, counts = np.unique(arr, return_counts=True)
+
+    # Identify top value
+    top_value_index = np.argmax(counts)
+
+    return values[top_value_index]
+
+
+def get_body_envelope(body_mask):
+    """
+    Given a mask representing thresholded lung values, obtain an envelope containing the lung region
+    with no interior holes.
+    """
+    # Invert the mask
+    reversed_mask = np.logical_not(body_mask)
+
+    # Identify connected_components
+    reversed_labels = skimage.measure.label(reversed_mask)
+
+    # Identify inner labels only (on x and y edges, z outer labels can remain)
+    inner_labels = lcat.util.clear_border(reversed_labels, axis=[0, 1])
+
+    # Obtain only the outermost (edge-touching) region
+    outside_mask = np.logical_xor(reversed_labels != 0, inner_labels != 0)
+
+    # Invert the mask again to get the envelope
+    envelope_mask = np.logical_not(outside_mask)
+
+    return envelope_mask

lcat/segmentation/lungs.py

@@ -26,8 +26,8 @@ def get_lung_segmentation(scan):
     threshold = skimage.filters.threshold_otsu(scan.voxels[scan_mask])
 
     # Make sure threshold is within air/lung parameters
-    if threshold <= -1000 or threshold >= -700:
-        threshold = -850
+    if threshold <= -1000 or threshold >= 0:
+        threshold = -500
 
     # Threshold the image
     foreground = scan.voxels >= threshold
@@ -42,14 +42,22 @@ def get_lung_segmentation(scan):
     lung_mask = get_largest_volume(labels)
 
     # Fill edge holes by dilation
-    smoother = skimage.morphology.ball(10, dtype=bool)
-    lung_mask = skimage.morphology.binary_dilation(lung_mask, selem=smoother)
+    # smoother = skimage.morphology.ball(10, dtype=bool)
+    # lung_mask = skimage.morphology.binary_dilation(lung_mask, selem=smoother)
+
+    smoother = skimage.morphology.disk(10, dtype=bool)
+    for z_index in range(lung_mask.shape[-1]):
+        lung_mask[..., z_index] = skimage.morphology.binary_dilation(lung_mask[..., z_index],
+                                                                     selem=smoother)
 
     # Obtain lung envelope
     envelope_mask = get_lung_envelope(lung_mask)
 
     # Erode envelope to revert to proper extent
-    envelope_mask = skimage.morphology.binary_erosion(envelope_mask, selem=smoother)
+    # envelope_mask = skimage.morphology.binary_erosion(envelope_mask, selem=smoother)
+    for z_index in range(lung_mask.shape[-1]):
+        lung_mask[..., z_index] = skimage.morphology.binary_erosion(envelope_mask[..., z_index],
+                                                                    selem=smoother)
 
     return envelope_mask
 
@@ -86,7 +94,7 @@ def get_lung_envelope(lung_mask):
     reversed_labels = skimage.measure.label(reversed_mask)
 
     # Identify inner labels only
-    inner_labels = skimage.segmentation.clear_border(reversed_labels)
+    inner_labels = lcat.util.clear_border(reversed_labels, axis=[0, 1])
 
     # Obtain only the outermost (edge-touching) region
     outside_mask = np.logical_xor(reversed_labels != 0, inner_labels != 0)

lcat/util.py

@@ -65,6 +65,7 @@ def save_slices(voxels, destination_folder, prefix="slice"):
         destination_filepath = os.path.join(destination_folder, filename_template % (z_index))
         scipy.misc.imsave(destination_filepath, voxels[..., z_index])
 
+
 def get_bounding_box(arr):
     """
     Given an array of values, returns an list of tuples, where each tuple represents the extent of
@@ -170,3 +171,24 @@ def clear_border(labels, axis=None, in_place=False):
         labels[labels == value] = 0
 
     return labels
+
+
+def image_from_mask(mask):
+    """
+    Convert a binary mask into a PIL image.
+    """
+    # Import PIL (if necessary)
+    import PIL
+
+    # Create mask image
+    mask_image = PIL.Image.new('1', mask.shape)
+
+    # Reserve memory for pixels
+    pixels = mask_image.load()
+
+    # Store pixel values
+    for i in range(mask_image.size[0]):
+        for j in range(mask_image.size[1]):
+            pixels[i, j] = int(mask[i, j])
+
+    return mask_image