mindboggle.guts package¶

Submodules¶

mindboggle.guts.compute module¶

Compute functions.

Authors:

Arno Klein, 2012-2016 (arno@mindboggle.info) http://binarybottle.com

compute_image_histogram(infile, nbins=100, threshold=0.0)¶

Compute histogram values from nibabel-readable image.

Parameters:	infile (string) – input nibabel-readable image file name nbins (integer) – number of bins threshold (float) – remove values lower than threshold
Returns:	histogram_values – histogram bin values
Return type:	numpy array

Examples

>>> from mindboggle.guts.compute import compute_image_histogram
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> labels_file = fetch_data(urls['freesurfer_labels'], '', '.nii.gz')
>>> nbins = 100
>>> threshold = 0.5
>>> histogram_values = compute_image_histogram(labels_file, nbins,
...                                            threshold)
>>> histogram_values[0:5]
array([102865, 119610,      0,      0,      0])

compute_overlaps(targets, list1, list2, output_file='', save_output=True, verbose=False)¶

Compute overlap for each target between two lists of numbers.

Parameters:

targets (list of integers) – targets to find in lists
list1 (1-D numpy array (or list) of numbers) –
list2 (1-D numpy array (or list) of numbers) –
output_file (string) – (optional) output file name
save_output (bool) – save output file?
verbose (bool) – print statements?

Returns:

dice_overlaps (numpy array) – Dice overlap values
jacc_overlaps (numpy array) – Jaccard overlap values
output_file (string) – output text file name with overlap values

Examples

>>> import numpy as np
>>> from mindboggle.guts.compute import compute_overlaps
>>> from mindboggle.mio.labels import DKTprotocol
>>> list1 = [1021, 1021, 1021, 1021, 1021, 1010, 1010, 1010, 1010, 1010]
>>> list2 = [1003, 1021, 1021, 1021, 1021, 1021, 1003, 1003, 1003, 1003]
>>> dkt = DKTprotocol()
>>> targets = dkt.cerebrum_cortex_DKT31_numbers
>>> output_file = ''
>>> save_output = True
>>> verbose = False
>>> dice_overlaps, jacc_overlaps, output_file = compute_overlaps(targets,
...     list1, list2, output_file, save_output, verbose)
>>> print('{0:.5f}'.format(dice_overlaps[18]))
0.80000
>>> print('{0:.5f}'.format(jacc_overlaps[18]))
0.66667

count_per_label(labels, include_labels=[], exclude_labels=[])¶

Compute the number of times each label occurs.

Parameters:

labels (numpy 1-D array of integers) – labels (e.g., one label per vertex of a mesh)
include_labels (list of integers) – labels to include (if empty, use unique numbers in image volume file)
exclude_labels (list of integers) – labels to be excluded

Returns:

unique_labels (list of integers) – unique label numbers
counts (list of floats) – number of times each label occurs

Examples

>>> from mindboggle.guts.compute import count_per_label
>>> labels = [8,8,8,8,8,10,11,12,10,10,11,11,11,12,12,12,12,13]
>>> include_labels = [9,10,11,12]
>>> exclude_labels = [13]
>>> unique_labels, counts = count_per_label(labels, include_labels,
...                                         exclude_labels)
>>> unique_labels
[9, 10, 11, 12]
>>> counts
[0, 3, 4, 5]

>>> import nibabel as nb
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.mio.labels import DKTprotocol
>>> from mindboggle.guts.compute import count_per_label
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> labels_file = fetch_data(urls['freesurfer_labels'], '', '.nii.gz')
>>> img = nb.load(labels_file)
>>> hdr = img.get_header()
>>> labels = img.get_data().ravel()
>>> dkt = DKTprotocol()
>>> include_labels = dkt.label_numbers
>>> exclude_labels = []
>>> unique_labels, counts = count_per_label(labels,
...     include_labels, exclude_labels)
>>> counts[0:10]
[0, 0, 0, 972, 2414, 2193, 8329, 2941, 1998, 10906]

distcorr(X, Y)¶

Compute the distance correlation function.

Parameters:	X (list or array of numbers) – Y (list or array of numbers) –
Returns:	dcor – distance correlation
Return type:	float

Examples

>>> from mindboggle.guts.compute import distcorr
>>> import numpy as np
>>> a = [1,2,3,4,5]
>>> b = [1,2,9,4,4]
>>> dcor = distcorr(a, b)
>>> np.float("{0:.{1}f}".format(dcor, 5))
0.76268

Copyright (2014-2015) MIT Written by Satrajit S. Ghosh (Apache License v2.0) as part of the mapalign GitHub repository (https://github.com/satra/mapalign)

means_per_label(values, labels, include_labels=[], exclude_labels=[], areas=[])¶

Compute the mean value across vertices per label, optionally taking into account surface area per vertex (UNTESTED).

Formula: average value = sum(a_i * v_i) / total_surface_area, where a_i and v_i are the area and value for each vertex i.

Note ::

This function is different than stats_per_label() in two ways:

It only computes the (weighted) mean and sdev.
It can accept 2-D arrays (such as [x,y,z] coordinates).

Parameters:

values (numpy array of one or more lists of integers or floats) – values to average per label
labels (list or array of integers) – label for each value
include_labels (list of integers) – labels to include
exclude_labels (list of integers) – labels to be excluded
areas (numpy array of floats) – surface areas (if provided, used to normalize means and sdevs)

Returns:

means (list of floats) – mean(s) for each label
sdevs (list of floats) – standard deviation(s) for each label
label_list (list of integers) – unique label numbers
label_areas (list of floats (if normalize_by_area)) – surface area for each labeled set of vertices

Examples

>>> import numpy as np
>>> from mindboggle.mio.vtks import read_scalars, read_vtk
>>> from mindboggle.guts.compute import means_per_label
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> values_file = fetch_data(urls['left_mean_curvature'], '', '.vtk')
>>> labels_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> area_file = fetch_data(urls['left_area'], '', '.vtk')
>>> values, name = read_scalars(values_file, True, True)
>>> labels, name = read_scalars(labels_file)
>>> areas, name = read_scalars(area_file, True, True)
>>> include_labels = []
>>> exclude_labels = [-1]
>>> # Compute mean curvature per label normalized by area:
>>> means, sdevs, label_list, label_areas = means_per_label(values, labels,
...     include_labels, exclude_labels, areas)
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in means[0:5]]
[-1.1793, -1.21405, -2.49318, -3.58116, -3.34987]
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in sdevs[0:5]]
[2.43827, 2.33857, 2.0185, 3.25964, 2.8274]
>>> # Compute mean curvature per label:
>>> areas = []
>>> means, sdevs, label_list, label_areas = means_per_label(values, labels,
...     include_labels, exclude_labels, areas)
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in means[0:5]]
[-0.99077, -0.3005, -1.59342, -2.03939, -2.31815]
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in sdevs[0:5]]
[2.3486, 2.4023, 2.3253, 3.31023, 2.91794]

>>> # FIX: compute mean coordinates per label:
>>> #points, indices, lines, faces, labels, scalar_names, npoints, input_vtk = read_vtk(values_file)
>>> #means, sdevs, label_list, label_areas = means_per_label(points, labels,
>>> #    include_labels, exclude_labels, areas)
>>> #means[0:3]
>>> #sdevs[0:3]

median_abs_dev(X, W=[], precision=1, c=1.0)¶

Compute the (weighted) median absolute deviation.

mad = median(abs(x - median(x))) / c

Parameters:	X (numpy array of floats or integers) – values W (numpy array of floats or integers) – weights precision (integer) – number of decimal places to consider weights c (float) – constant used as divisor for mad computation; c = 0.6745 is used to convert from mad to standard deviation
Returns:	mad – (weighted) median absolute deviation
Return type:	float

Examples

>>> import numpy as np
>>> from mindboggle.guts.compute import median_abs_dev
>>> X = np.array([1,2,4,7,8])
>>> W = np.array([.1,.1,.3,.2,.3])
>>> precision = 1
>>> # [1, 2, 4, 4, 4, 7, 7, 8, 8, 8]
>>> median_abs_dev(X, W, precision)
2.0

pairwise_vector_distances(vectors, save_file=False, normalize=False)¶

Compare every pair of equal-sized vectors.

Parameters:

vectors (array of 1-D lists or arrays of integers or floats) –
save_file (bool) – save file?
normalize (bool) – normalize each element of the vectors?

Returns:

vector_distances (numpy array of integers or floats) – distances between each pair of vectors
outfile (string [optional]) – output filename for pairwise_vector_distances

Examples

>>> import numpy as np
>>> from mindboggle.guts.compute import pairwise_vector_distances
>>> vectors = [[1,2,3],[0,3,5],[0,3.5,5],[1,1,1]]
>>> save_file = False
>>> normalize = False
>>> vector_distances, outfile = pairwise_vector_distances(vectors,
...     save_file, normalize)
>>> print(np.array_str(np.array(vector_distances),
...       precision=5, suppress_small=True))
[[0.      0.8165  0.89753 0.74536]
 [0.      0.      0.16667 1.52753]
 [0.      0.      0.      1.60728]
 [0.      0.      0.      0.     ]]

point_distance(point, points)¶

Compute the Euclidean distance from one point to a second (set) of points.

Parameters:

point (list of three floats) – coordinates for a single point
points (list with one or more lists of three floats) – coordinates for a second point (or multiple points)

Returns:

min_distance (float) – Euclidean distance between two points, or the minimum distance between a point and a set of points
min_index (int) – index of closest of the points (zero if only one)

Examples

>>> from mindboggle.guts.compute import point_distance
>>> point = [1,2,3]
>>> points = [[10,2.0,3], [0,1.5,2]]
>>> point_distance(point, points)
(1.5, 1)

Notes

Future plan is to use scipy.spatial.distance.cdist to compute distances scipy.spatial.distance.cdist is available in scipy v0.12 or later

source_to_target_distances(sourceIDs, targetIDs, points, segmentIDs=[], excludeIDs=[-1])¶

Create a Euclidean distance matrix between source and target points.

Compute the Euclidean distance from each source point to its nearest target point, optionally within each segment.

Example:

Compute fundus-to-feature distances, the minimum distance
from each label boundary vertex (corresponding to a fundus
in the DKT cortical labeling protocol) to all of the
feature vertices in the same fold.

Parameters:

sourceIDs (list of N integers (N is the number of vertices)) – source IDs, where any ID not in excludeIDs is a source point
targetIDs (list of N integers (N is the number of vertices)) – target IDs, where any ID not in excludeIDs is a target point
points (list of N lists of three floats (N is the number of vertices)) – coordinates of all vertices
segmentIDs (list of N integers (N is the number of vertices)) – segment IDs, where each ID not in excludeIDs is considered a different segment (unlike above, where value in sourceIDs or targetIDs doesn’t matter, so long as its not in excludeIDs); source/target distances are computed within each segment
excludeIDs (list of integers) – IDs to exclude

Returns:

distances (numpy array) – distance value for each vertex (default -1)
distance_matrix (numpy array [#points by maximum segment ID + 1]) – distances organized by segments (columns)

stats_per_label(values, labels, include_labels=[], exclude_labels=[], weights=[], precision=1)¶

Compute various statistical measures across vertices per label, optionally using weights (such as surface area per vertex).

Example (area-weighted mean): average value = sum(a_i * v_i) / total_surface_area, where a_i and v_i are the area and value for each vertex i.

Reference:

Weighted skewness and kurtosis unbiased by sample size Lorenzo Rimoldini, arXiv:1304.6564 (2013) http://arxiv.org/abs/1304.6564

Note ::

This function is different than means_per_label() in two ways:

It computes more than simply the (weighted) mean and sdev.
It only accepts 1-D arrays of values.

Parameters:

values (numpy array of individual or lists of integers or floats) – values for all vertices
labels (list or array of integers) – label for each value
include_labels (list of integers) – labels to include
exclude_labels (list of integers) – labels to be excluded
weights (numpy array of floats) – weights to compute weighted statistical measures
precision (integer) – number of decimal places to consider weights

Returns:

medians (list of floats) – median for each label
mads (list of floats) – median absolute deviation for each label
means (list of floats) – mean for each label
sdevs (list of floats) – standard deviation for each label
skews (list of floats) – skew for each label
kurts (list of floats) – kurtosis value for each label
lower_quarts (list of floats) – lower quartile for each label
upper_quarts (list of floats) – upper quartile for each label
label_list (list of integers) – list of unique labels

Examples

>>> import numpy as np
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.guts.compute import stats_per_label
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> values_file = fetch_data(urls['left_mean_curvature'], '', '.vtk')
>>> labels_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> area_file = fetch_data(urls['left_area'], '', '.vtk')
>>> values, name = read_scalars(values_file, True, True)
>>> areas, name = read_scalars(area_file, True, True)
>>> labels, name = read_scalars(labels_file)
>>> include_labels = []
>>> exclude_labels = [-1]
>>> weights = areas
>>> precision = 1
>>> medians, mads, means, sdevs, skews, kurts, lower_quarts, upper_quarts, label_list = stats_per_label(values,
...     labels, include_labels, exclude_labels, weights, precision)
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in medians[0:5]]
[-1.13602, -1.22961, -2.49665, -3.80782, -3.37309]
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in mads[0:5]]
[1.17026, 1.5045, 1.28234, 2.11515, 1.69333]
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in means[0:5]]
[-1.1793, -1.21405, -2.49318, -3.58116, -3.34987]
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in kurts[0:5]]
[2.34118, -0.3969, -0.55787, -0.73993, 0.3807]

sum_per_label(values, labels, include_labels=[], exclude_labels=[])¶

Compute the sum value across vertices per label.

Parameters:

values (numpy array of one or more lists of integers or floats) – values to average per label
labels (list or array of integers) – label for each value
include_labels (list of integers) – labels to include
exclude_labels (list of integers) – labels to be excluded

Returns:

sums (list of floats) – sum for each label
label_list (list of integers) – unique label numbers

Examples

>>> import numpy as np
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.guts.compute import sum_per_label
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> values_file = fetch_data(urls['left_mean_curvature'], '', '.vtk')
>>> labels_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> values, name = read_scalars(values_file, True, True)
>>> labels, name = read_scalars(labels_file)
>>> include_labels = []
>>> exclude_labels = [-1]
>>> # Compute sum area per label:
>>> sums, label_list = sum_per_label(values, labels, include_labels,
...                                  exclude_labels)
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in sums[0:5]]
[-8228.32913, -424.90109, -1865.8959, -8353.33769, -5130.06613]

vector_distance(vector1, vector2, normalize=False)¶

Compute the Euclidean distance between two equal-sized vectors.

Parameters:	vector1 (numpy array of floats) – vector of values vector2 (numpy array of floats) – vector of values normalize (bool) – normalize each element of the vectors?
Returns:	distance – Euclidean distance between two vectors
Return type:	float

Examples

>>> import numpy as np
>>> from mindboggle.guts.compute import vector_distance
>>> vector1 = np.array([1.,2.,3.])
>>> vector2 = np.array([0,1,5])
>>> distance = vector_distance(vector1, vector2)
>>> print('{0:0.5f}'.format(distance))
0.81650

weighted_median(X, W=[], precision=1)¶

Compute a weighted median.

Parameters:	X (numpy array of floats or integers) – values W (numpy array of floats or integers) – weights precision (integer) – number of decimal places to consider weights
Returns:	wmedian – weighted median
Return type:	float

Examples

>>> import numpy as np
>>> from mindboggle.guts.compute import weighted_median
>>> X = np.array([1,2,4,7,8])
>>> W = np.array([.1,.1,.3,.2,.3])
>>> precision = 1
>>> # [1, 2, 4, 4, 4, 7, 7, 8, 8, 8]
>>> weighted_median(X, W, precision)
5.5

weighted_to_repeated_values(X, W=[], precision=1)¶

Create a list of repeated values from weighted values.

This is useful for computing weighted statistics (ex: weighted median).

Adapted to allow for fractional weights from

http://stackoverflow.com/questions/966896/: code-golf-shortest-code-to-find-a-weighted-median

Parameters:	X (numpy array of floats or integers) – values W (numpy array of floats or integers) – weights precision (integer) – number of decimal places to consider weights
Returns:	repeat_values – repeated values according to their weight
Return type:	numpy array of floats or integers

Examples

>>> import numpy as np
>>> from mindboggle.guts.compute import weighted_to_repeated_values
>>> X = np.array([1,2,4,7,8])
>>> W = np.array([.1,.1,.3,.2,.3])
>>> precision = 1
>>> weighted_to_repeated_values(X, W, precision)
[1, 2, 4, 4, 4, 7, 7, 8, 8, 8]

mindboggle.guts.graph module¶

mindboggle.guts.kernels module¶

Kernels.

Authors:

Eliezer Stavsky, 2012 (eli.stavsky@gmail.com)
Arno Klein, 2016 (arno@mindboggle.info)

inverse_distance(x1, x2, epsilon)¶

This function constructs weighted edges of a graph, where the weight is the inverse of the distance between two nodes.

Parameters:	x1 (Nx1 numpy array) – x2 (Nx1 numpy array) – epsilon (float) –
Returns:	d
Return type:	float

Examples

>>> import numpy as np
>>> from mindboggle.guts.kernels import inverse_distance
>>> x1 = np.array([0.1,0.2,0.4,0])
>>> x2 = np.array([0.1,0.3,0.5,0])
>>> epsilon = 0.05
>>> d = inverse_distance(x1, x2, epsilon)
>>> print('{0:0.5f}'.format(d))
5.22408

rbf_kernel(x1, x2, sigma)¶

Compute normalized and unnormalized graph Laplacians.

Parameters:	x1 (Nx1 numpy array) – x2 (Nx1 numpy array) – sigma (float) –
Returns:	rbf
Return type:	float

Examples

>>> import numpy as np
>>> from mindboggle.guts.kernels import rbf_kernel
>>> x1 = np.array([0.1,0.2,0.4,0])
>>> x2 = np.array([0.1,0.3,0.5,0])
>>> sigma = 0.5
>>> rbf = rbf_kernel(x1, x2, sigma)
>>> print('{0:0.5f}'.format(rbf))
0.96079

mindboggle.guts.mesh module¶

Operations on surface mesh vertices.

Authors:

Forrest Bao, 2012 (forrest.bao@gmail.com)
Arno Klein, 2012-2016 (arno@mindboggle.info) http://binarybottle.com

area_of_faces(points, faces)¶

Compute the areas of all triangles on the mesh.

Parameters:	points (list of lists of 3 floats) – x,y,z coordinates for each vertex of the structure faces (list of lists of 3 integers) – 3 indices to vertices that form a triangle on the mesh
Returns:	area – area[i] is the area of the i-th triangle
Return type:	1-D numpy array

Examples

>>> import numpy as np
>>> from mindboggle.guts.mesh import area_of_faces
>>> from mindboggle.mio.vtks import read_vtk
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> input_vtk = fetch_data(urls['left_area'], '', '.vtk')
>>> points, f1, f2, faces, f3, f4, f5, f6 = read_vtk(input_vtk)
>>> area = area_of_faces(points, faces)
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in area[0:5]]
[0.21703, 0.27139, 0.29033, 0.1717, 0.36011]

decimate(points, faces, reduction=0.75, smooth_steps=25, scalars=[], save_vtk=False, output_vtk='')¶

Decimate vtk triangular mesh with vtk.vtkDecimatePro.

Parameters:

points (list of lists of floats) – each element is a list of 3-D coordinates of a vertex on a surface mesh
faces (list of lists of integers) – each element is list of 3 indices of vertices that form a face on a surface mesh
reduction (float) – fraction of mesh faces to remove
smooth_steps (integer) – number of smoothing steps
scalars (list of integers or floats) – optional scalars for output VTK file
save_vtk (bool) – output decimated vtk file?
output_vtk (string) – output decimated vtk file name

Returns:

points (list of lists of floats) – decimated points
faces (list of lists of integers) – decimated faces
scalars (list of integers or floats) – scalars for output VTK file
output_vtk (string) – output decimated vtk file

Examples

>>> # Example: Twins-2-1 left postcentral pial surface, 0.75 decimation:
>>> from mindboggle.guts.mesh import decimate
>>> from mindboggle.mio.vtks import read_vtk
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> input_vtk = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> points, f1, f2, faces, scalars, f3, f4, f5 = read_vtk(input_vtk)
>>> reduction = 0.5
>>> smooth_steps = 25
>>> save_vtk = True
>>> output_vtk = 'decimate.vtk'
>>> points2, faces2, scalars, output_vtk = decimate(points, faces,
...     reduction, smooth_steps, scalars, save_vtk, output_vtk)
>>> (len(points), len(points2))
(145069, 72535)
>>> (len(faces), len(faces2))
(290134, 145066)

View decimated surface (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> plot_surfaces('decimate.vtk') # doctest: +SKIP

decimate_file(input_vtk, reduction=0.5, smooth_steps=100, save_vtk=True, output_vtk='')¶

Decimate vtk triangular mesh file with vtk.vtkDecimatePro.

Parameters:	input_vtk (string) – input vtk file with triangular surface mesh reduction (float) – fraction of mesh faces to remove do_smooth (bool) – smooth after decimation? save_vtk (bool) – output decimated vtk file? output_vtk (string) – output decimated vtk file name
Returns:	output_vtk – output decimated vtk file
Return type:	string

Examples

>>> from mindboggle.guts.mesh import decimate_file
>>> from mindboggle.mio.vtks import read_vtk
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> input_vtk = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> save_vtk = True
>>> output_vtk = 'decimate.vtk'
>>> reduction = 0.5
>>> smooth_steps = 25
>>> output_vtk = decimate_file(input_vtk, reduction, smooth_steps,
...     save_vtk, output_vtk)
>>> f1, f2, f3, faces1, f4, f5, npoints1, f6 = read_vtk(input_vtk)
>>> f1, f2, f3, faces2, f4, f5, npoints2, f6 = read_vtk('decimate.vtk')
>>> (npoints1, npoints2)
(145069, 72535)
>>> (len(faces1), len(faces2))
(290134, 145066)

View decimated surface (skip test):

>>> from mindboggle.mio.plots import plot_surfaces
>>> plot_surfaces('decimate.vtk') # doctest: +SKIP

dilate(indices, nedges, neighbor_lists)¶

Dilate region on a surface mesh.

Parameters:	indices (list of integers) – indices of vertices to dilate nedges (integer) – number of edges to dilate across neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex
Returns:	dilated_indices – indices of original vertices with dilated vertices
Return type:	list of integers

Examples

>>> import numpy as np
>>> from mindboggle.guts.mesh import dilate, find_neighbors_from_file
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> vtk_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> neighbor_lists = find_neighbors_from_file(vtk_file)
>>> nedges = 3
>>> # Select a single fold:
>>> folds, name = read_scalars(folds_file, True, True)
>>> fold_number = 4
>>> indices = [i for i,x in enumerate(folds) if x == fold_number]
>>> dilated_indices = dilate(indices, nedges, neighbor_lists)
>>> (len(indices), len(dilated_indices))
(1151, 1545)
>>> dilated_indices[0:10]
[50317, 50324, 50325, 50326, 50327, 50332, 50333, 50334, 50339, 50340]

Write results to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars
>>> IDs = -1 * np.ones(len(folds)) # doctest: +SKIP
>>> IDs[dilated_indices] = 2 # doctest: +SKIP
>>> IDs[indices] = 1 # doctest: +SKIP
>>> rewrite_scalars(vtk_file, 'dilate.vtk', IDs, 'dilated_fold', IDs) # doctest: +SKIP
>>> plot_surfaces('dilate.vtk') # doctest: +SKIP

erode(indices, nedges, neighbor_lists)¶

Erode region on a surface mesh.

Parameters:	indices (list of integers) – indices of vertices to erode nedges (integer) – number of edges to erode across neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex
Returns:	eroded_indices – indices of original vertices without eroded vertices
Return type:	list of integers

Examples

>>> import numpy as np
>>> from mindboggle.guts.mesh import erode, find_neighbors_from_file
>>> from mindboggle.mio.vtks import read_scalars, rewrite_scalars
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> vtk_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> neighbor_lists = find_neighbors_from_file(vtk_file)
>>> nedges = 3
>>> # Select a single fold:
>>> folds, name = read_scalars(folds_file, True, True)
>>> fold_number = 4
>>> indices = [i for i,x in enumerate(folds) if x == fold_number]
>>> eroded_indices = erode(indices, nedges, neighbor_lists)
>>> (len(indices), len(eroded_indices))
(1151, 809)

Write results to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> IDs = -1 * np.ones(len(folds)) # doctest: +SKIP
>>> IDs[indices] = 1 # doctest: +SKIP
>>> IDs[eroded_indices] = 2 # doctest: +SKIP
>>> rewrite_scalars(vtk_file, 'erode.vtk', IDs, 'eroded_fold', IDs) # doctest: +SKIP
>>> plot_surfaces('erode.vtk') # doctest: +SKIP

extract_edge(indices, neighbor_lists)¶

Erode region on a surface mesh to extract the region’s edge.

Parameters:	indices (list of integers) – indices of vertices to erode neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex
Returns:	edge_indices – indices of eroded vertices
Return type:	list of integers

Examples

>>> import numpy as np
>>> from mindboggle.guts.mesh import extract_edge
>>> from mindboggle.guts.mesh import find_neighbors_from_file
>>> from mindboggle.mio.vtks import read_scalars, rewrite_scalars
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> vtk_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> neighbor_lists = find_neighbors_from_file(vtk_file)
>>> # Select a single fold:
>>> folds, name = read_scalars(folds_file, True, True)
>>> fold_number = 4
>>> indices = [i for i,x in enumerate(folds) if x == fold_number]
>>> edge_indices = extract_edge(indices, neighbor_lists)
>>> (len(indices), len(edge_indices))
(1151, 111)

Write results to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> IDs = -1 * np.ones(len(folds)) # doctest: +SKIP
>>> IDs[indices] = 1 # doctest: +SKIP
>>> IDs[edge_indices] = 2 # doctest: +SKIP
>>> rewrite_scalars(vtk_file, 'extract_edge.vtk', IDs, 'edges_of_fold', IDs) # doctest: +SKIP
>>> plot_surfaces('extract_edge.vtk') # doctest: +SKIP

find_adjacent_faces(faces)¶

For each face in a list of faces, find adjacent faces.

Parameters:	faces (list of lists of three integers) – the integers for each face are indices to vertices, starting from zero (-1 indicates no result for a given face or vertex)
Returns:	adjacent_faces – list 1 indexes three faces adjacent to the three face’s edges; list 2 indexes three vertices opposite the adjacent faces: adjacent_faces[i][0] = [face0, face1, face2], neighbors of face i (face0 is the neighbor of face i facing vertex0) adjacent_faces[i][1] = [vertex0, vertex1, vertex2] for face i (vertex0 is the vertex of face0 not in face i)
Return type:	list of pairs of lists of three integers

Examples

>>> # Simple example:
>>> from mindboggle.guts.mesh import find_adjacent_faces
>>> faces = [[0,1,2],[0,2,3],[0,3,4],[0,1,4],[4,3,1]]
>>> adjacent_faces = find_adjacent_faces(faces)
>>> adjacent_faces[0:2]
[[[-1, 1, 3], [-1, 3, 4]], [[-1, 2, 0], [-1, 4, 1]]]

find_complete_faces(indices, faces)¶

Given a set of vertices, find the ones that make complete faces.

Parameters:	indices (list of integers) – indices to connected vertices faces (list of lists of three integers) – the integers for each face are indices to vertices, starting from zero
Returns:	indices_complete – indices to vertices making up complete faces
Return type:	list of integers

Examples

>>> from mindboggle.guts.mesh import find_complete_faces
>>> faces = [[0,2,3], [2,3,7], [4,7,8], [3,2,5]]
>>> indices = [3,7,2,5,9,4]
>>> find_complete_faces(indices, faces)
[2, 3, 7, 5]

find_edges(faces)¶

Find all edges on a mesh

Parameters:	faces (list of lists of three integers) – the integers for each face are indices to vertices, starting from zero
Returns:	edges – each element is a 2-tuple of vertex ids representing an edge
Return type:	list of lists of integers

Examples

>>> # Simple example:
>>> from mindboggle.guts.mesh import find_edges
>>> faces=[[0,1,2], [0,1,4], [1,2,3], [0,2,5]]
>>> find_edges(faces)
[[0, 1], [1, 2], [0, 2], [1, 4], [0, 4], [2, 3], [1, 3], [2, 5], [0, 5]]

find_endpoints(indices, neighbor_lists)¶

Extract endpoints from connected set of vertices.

Parameters:	indices (list of integers) – indices to connected vertices neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex
Returns:	indices_endpoints – indices to endpoints of connected vertices
Return type:	list of integers

Examples

>>> # Find endpoints of fundus in a fold:
>>> from mindboggle.guts.mesh import find_endpoints
>>> from mindboggle.guts.mesh import find_neighbors_from_file
>>> from mindboggle.mio.fetch_data import prep_tests
>>> from mindboggle.mio.vtks import read_scalars
>>> urls, fetch_data = prep_tests()
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> fundus_file = fetch_data(urls['left_fundus_per_fold'], '', '.vtk')
>>> folds, name = read_scalars(folds_file, True, True)
>>> fundi, name = read_scalars(fundus_file, True, True)
>>> background_value = -1
>>> # Limit number of folds to speed up the test:
>>> limit_folds = True
>>> if limit_folds:
...     fold_numbers = [2]
...     i0 = [i for i,x in enumerate(folds) if x not in fold_numbers]
...     folds[i0] = background_value
...     fundi[i0] = background_value
...     indices = [i for i,x in enumerate(fundi) if x != background_value]
>>> neighbor_lists = find_neighbors_from_file(fundus_file)
>>> indices_endpoints = find_endpoints(indices, neighbor_lists)
>>> indices_endpoints[0:5]
[32782, 35142, 45244, 49010, 63051]

View endpoints (skip test):

>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> fundi[indices_endpoints] = 50 # doctest: +SKIP
>>> rewrite_scalars(fundus_file, 'find_endpoints.vtk', fundi,
...     'endpoints', folds, background_value) # doctest: +SKIP
>>> plot_surfaces('find_endpoints.vtk') # doctest: +SKIP

find_faces_at_edges(faces)¶

For each edge on the mesh, find the two faces that share the edge.

Parameters:	faces (list of lists of three integers) – the integers for each face are indices to vertices, starting from zero
Returns:	faces_at_edges – keys are tuples of two vertex IDs and values are 2-tuples of face IDs
Return type:	dictionary

Examples

>>> # Simple example:
>>> from mindboggle.guts.mesh import find_faces_at_edges
>>> faces=[[0,1,2], [0,1,4], [1,2,3], [0,2,5]]
>>> faces_at_edges = find_faces_at_edges(faces)
>>> faces_at_edges[(0,2)]
[0, 3]
>>> faces_at_edges[(2,1)]
[0, 2]

Notes ::: The faces are assumed to be triangular.

find_faces_at_vertices(faces, npoints)¶

For each vertex, find all faces containing this vertex. Note: faces do not have to be triangles.

Parameters:	faces (list of lists of three integers) – the integers for each face are indices to vertices, starting from zero npoints (integer) – number of vertices on the mesh
Returns:	faces_at_vertices – faces_at_vertices[i] is a list of faces that contain the i-th vertex
Return type:	list of lists of integers

Examples

>>> # Simple example:
>>> from mindboggle.guts.mesh import find_faces_at_vertices
>>> faces = [[0,1,2],[0,2,3],[0,3,4],[0,1,4],[4,3,1]]
>>> npoints = 5
>>> find_faces_at_vertices(faces, npoints)
[[0, 1, 2, 3], [0, 3, 4], [0, 1], [1, 2, 4], [2, 3, 4]]

find_faces_with_vertex(index, faces)¶

For a given vertex, find all faces containing this vertex. Note: faces do not have to be triangles.

Parameters:	index (integer) – index to a vertex faces (list of lists of three integers) – the integers for each face are indices to vertices, starting from zero
Returns:	faces_with_vertex – the integers for each face are indices to vertices, starting from zero
Return type:	list of lists of three integers

Examples

>>> # Simple example:
>>> from mindboggle.guts.mesh import find_faces_with_vertex
>>> faces = [[0,1,2],[0,2,3],[0,3,4],[0,1,4],[4,3,1]]
>>> index = 3
>>> find_faces_with_vertex(index, faces)
[[0, 2, 3], [0, 3, 4], [4, 3, 1]]

find_neighborhood(neighbor_lists, indices, nedges=1)¶

Find neighbors in the neighborhood of given surface mesh vertices.

For indices to surface mesh vertices, find unique indices for vertices in the neighborhood of the vertices.

Parameters:	neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex indices (list of integers) – indices of surface vertices nedges (integer) – number of edges to propagate from indices
Returns:	neighborhood – indices to vertices in neighborhood
Return type:	list of integers

Examples

>>> from mindboggle.guts.mesh import find_neighborhood
>>> neighbor_lists = [[0,1],[0,2],[1,4,5],[2],[],[0,1,4,5]]
>>> indices = [1,3,4]
>>> neighborhood = find_neighborhood(neighbor_lists, indices, 2)
>>> neighborhood
[0, 2, 5]

find_neighbors(faces, npoints)¶

Generate the list of unique, sorted indices of neighboring vertices for all vertices in the faces of a triangular mesh.

Parameters:	faces (list of lists of three integers) – the integers for each face are indices to vertices, starting from zero npoints (integer) – number of vertices on the mesh
Returns:	neighbor_lists – each list contains indices to neighboring vertices for each vertex
Return type:	list of lists of integers

Examples

>>> # Simple example:
>>> from mindboggle.guts.mesh import find_neighbors
>>> faces = [[0,1,2],[0,2,3],[0,3,4],[0,1,4],[4,3,1]]
>>> npoints = 5
>>> neighbor_lists = find_neighbors(faces, npoints)
>>> neighbor_lists
[[1, 2, 3, 4], [0, 2, 4, 3], [0, 1, 3], [0, 2, 4, 1], [0, 3, 1]]

Real example:

>>> import numpy as np
>>> from mindboggle.guts.mesh import find_neighbors
>>> from mindboggle.mio.vtks import read_faces_points
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> vtk_file = fetch_data(urls['left_mean_curvature'], '', '.vtk')
>>> faces, points, npoints = read_faces_points(vtk_file)
>>> neighbor_lists = find_neighbors(faces, npoints)
>>> neighbor_lists[0:3]
[[1, 4, 48, 49], [0, 4, 5, 49, 2], [1, 5, 6, 49, 50, 54]]

Write results to vtk file and view (skip test):

>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> index = 100 # doctest: +SKIP
>>> IDs = -1 * np.ones(len(neighbor_lists)) # doctest: +SKIP
>>> IDs[index] = 1 # doctest: +SKIP
>>> IDs[neighbor_lists[index]] = 2 # doctest: +SKIP
>>> rewrite_scalars(vtk_file, 'find_neighbors.vtk', IDs, 'neighbors', IDs) # doctest: +SKIP
>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> plot_surfaces('find_neighbors.vtk') # doctest: +SKIP

find_neighbors_from_file(input_vtk)¶

Generate the list of unique, sorted indices of neighboring vertices for all vertices in the faces of a triangular mesh in a VTK file.

Parameters:	input_vtk (string) – name of input VTK file containing surface mesh
Returns:	neighbor_lists – each list contains indices to neighboring vertices for each vertex
Return type:	list of lists of integers

Examples

>>> import numpy as np
>>> from mindboggle.guts.mesh import find_neighbors_from_file
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> vtk_file = fetch_data(urls['left_mean_curvature'], '', '.vtk')
>>> neighbor_lists = find_neighbors_from_file(vtk_file)
>>> neighbor_lists[0:3]
[[1, 4, 48, 49], [0, 4, 5, 49, 2], [1, 5, 6, 49, 50, 54]]

Write results to vtk file and view (skip test):

>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> index = 100 # doctest: +SKIP
>>> IDs = -1 * np.ones(len(neighbor_lists)) # doctest: +SKIP
>>> IDs[index] = 1 # doctest: +SKIP
>>> IDs[neighbor_lists[index]] = 2 # doctest: +SKIP
>>> rewrite_scalars(vtk_file, 'find_neighbors_from_file.vtk', IDs,
...                 'neighbors', IDs) # doctest: +SKIP
>>> plot_surfaces('find_neighbors_from_file.vtk') # doctest: +SKIP

find_neighbors_vertex(faces, index)¶

Find neighbors to a surface mesh vertex.

For a set of surface mesh faces and the index of a surface vertex, find unique indices for neighboring vertices.

Parameters:	faces (list of lists of three integers) – the integers for each face are indices to vertices, starting from zero index (int) – index of surface vertex
Returns:	neighbor_lists – each list contains indices to neighboring vertices for each vertex
Return type:	list of lists of integers

Examples

>>> from mindboggle.guts.mesh import find_neighbors_vertex
>>> faces = [[0,1,2],[0,2,3],[0,3,4],[0,1,4]]
>>> index = 1
>>> neighbor_lists = find_neighbors_vertex(faces, index)
>>> neighbor_lists
[0, 2, 4]

keep_faces(faces, indices)¶

Remove surface mesh faces whose three vertices are not all in “indices”.

Parameters:	faces (list of lists of three integers) – the integers for each face are indices to vertices, starting from zero indices (list of integers) – indices to vertices of the surface mesh that are to be retained
Returns:	faces – reduced number of faces
Return type:	list of lists of three integers

Examples

>>> from mindboggle.guts.mesh import keep_faces
>>> faces = [[1,2,3], [2,3,7], [4,7,8], [3,2,5]]
>>> indices = [0,1,2,3,4,5]
>>> keep_faces(faces, indices)
[[1, 2, 3], [3, 2, 5]]

reindex_faces_0to1(faces)¶

Convert 0-indices (Python) to 1-indices (Matlab) for all face indices.

Parameters:	faces (list of lists of integers) – each sublist contains 3 0-indices of vertices that form a face on a surface mesh
Returns:	faces – each sublist contains 3 1-indices of vertices that form a face on a surface mesh
Return type:	list of lists of integers

Examples

>>> from mindboggle.guts.mesh import reindex_faces_0to1
>>> faces = [[0,2,3], [2,3,7], [4,7,8], [3,2,5]]
>>> reindex_faces_0to1(faces)
[[1, 3, 4], [3, 4, 8], [5, 8, 9], [4, 3, 6]]

reindex_faces_points(faces, points=[])¶

Renumber indices in faces and remove points (coordinates) not in faces.

Parameters:

faces (list of lists of integers) – each sublist contains 3 indices of vertices that form a face on a surface mesh
points (list of lists of floats (optional)) – each sublist contains 3-D coordinates of a vertex on a surface mesh

Returns:

new_faces (list of lists of integers) – each sublist contains 3 (renumbered) indices of vertices that form a face on a surface mesh
new_points (list of lists of floats) – each (new) sublist contains 3-D coordinates of a vertex on a surface mesh
original_indices (list integers) – list of indices to original points

Examples

>>> from mindboggle.guts.mesh import reindex_faces_points
>>> # Reindex faces:
>>> faces = [[8,2,3], [2,3,7], [4,7,8], [3,2,5]]
>>> new_faces, new_points, original_indices = reindex_faces_points(faces,
...     points=[])
>>> new_faces
[[5, 0, 1], [0, 1, 4], [2, 4, 5], [1, 0, 3]]

Reindex faces of a limited number of folds of the brain:

>>> import numpy as np
>>> from mindboggle.guts.mesh import keep_faces
>>> from mindboggle.mio.vtks import read_faces_points
>>> from mindboggle.mio.vtks import read_scalars, rewrite_scalars
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> folds, name = read_scalars(folds_file, True, True)
>>> fold_numbers = [4]
>>> indices = [i for i,x in enumerate(folds) if x in fold_numbers]
>>> i0 = [i for i,x in enumerate(folds) if x not in fold_numbers]
>>> background_value = -1
>>> folds[i0] = background_value
>>> faces, points, npoints = read_faces_points(folds_file)
>>> faces = keep_faces(faces, indices)
>>> faces[0:3]
[[51535, 50324, 51529], [50317, 50325, 50326], [50324, 50332, 50333]]
>>> new_faces, new_points, original_indices = reindex_faces_points(faces,
...     points)
>>> new_faces[0:3]
[[277, 690, 276], [689, 691, 692], [690, 698, 699]]
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in points[0]]
[-13.7924, -76.0973, -2.57594]
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in new_points[0]]
[-13.7802, -12.3814, 57.4042]

View reindexed fold on surface (skip test):

>>> from mindboggle.mio.plots import plot_surfaces
>>> plot_surfaces('reindex_faces_points.vtk') # doctest: +SKIP

remove_neighbor_lists(neighbor_lists, indices)¶

Remove all but a given set of indices from surface mesh neighbor lists.

Note :: SLOW!

Parameters:	neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex indices (list of integers) – indices to vertices of the surface mesh
Returns:	neighbor_lists – each list has indices to remaining neighboring vertices for each vertex
Return type:	list of lists of integers

Examples

>>> from mindboggle.guts.mesh import remove_neighbor_lists
>>> neighbor_lists = [[1,2,3], [2,3,7], [12,43], [4,7,8], [3,2,5]]
>>> indices = [0,1,2,3,4,5]
>>> remove_neighbor_lists(neighbor_lists, indices)
[[1, 2, 3], [2, 3], [], [4], [2, 3, 5]]

rescale_by_label(input_vtk, labels_or_file, save_file=False, output_filestring='rescaled_scalars', background_value=-1, verbose=False)¶

Rescale scalars for each label (such as depth values within each fold).

Default is to normalize the scalar values of a VTK file by a percentile value in each vertex’s surface mesh for each label.

Parameters:

input_vtk (string) – name of VTK file with a scalar value for each vertex
labels_or_file (list or string) – label number for each vertex or name of VTK file with index scalars
save_file (bool) – save output VTK file?
output_filestring (string (if save_file)) – name of output file
background_value (integer or float) – background value
verbose (bool) – print statements?

Returns:

rescaled_scalars (list of floats) – scalar values rescaled for each label, for label numbers not equal to -1
rescaled_scalars_file (string (if save_file)) – name of output VTK file with rescaled scalar values for each label

Examples

>>> # Rescale depths by neighborhood within each label:
>>> import numpy as np
>>> from mindboggle.guts.mesh import rescale_by_label
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.mio.plots import plot_surfaces
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> input_vtk = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> labels_or_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> save_file = True
>>> output_filestring = 'rescale_by_label'
>>> background_value = -1
>>> verbose = False
>>> rescaled, rescaled_label_file = rescale_by_label(input_vtk,
...     labels_or_file, save_file, output_filestring, background_value, verbose)
>>> scalars1, name = read_scalars(input_vtk)
>>> print('{0:0.5f}, {1:0.5f}'.format(max(scalars1), max(rescaled)))
34.95560, 1.00000
>>> print('{0:0.5f}, {1:0.5f}'.format(np.mean(scalars1), np.mean(rescaled)))
7.43822, 0.30677

View rescaled scalar values on surface (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> plot_surfaces(rescaled_label_file) # doctest: +SKIP

rescale_by_neighborhood(input_vtk, indices=[], nedges=10, p=99, set_max_to_1=True, save_file=False, output_filestring='rescaled_scalars', background_value=-1)¶

Rescale the scalar values of a VTK file by a percentile value in each vertex’s surface mesh neighborhood.

Parameters:

input_vtk (string) – name of VTK file with a scalar value for each vertex
indices (list of integers (optional)) – indices of scalars to normalize
nedges (integer) – number or edges from vertex, defining the size of its neighborhood
p (float in range of [0,100]) – percentile used to normalize each scalar
set_max_to_1 (bool) – set all rescaled values greater than 1 to 1.0?
save_file (bool) – save output VTK file?
output_filestring (string (if save_file)) – name of output file
background_value (integer) – background value

Returns:

rescaled_scalars (list of floats) – rescaled scalar values
rescaled_scalars_file (string (if save_file)) – name of output VTK file with rescaled scalar values

Examples

>>> import numpy as np
>>> from mindboggle.guts.mesh import rescale_by_neighborhood
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.mio.plots import plot_surfaces
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> input_vtk = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> indices = []
>>> nedges = 10
>>> p = 99
>>> set_max_to_1 = True
>>> save_file = True
>>> output_filestring = 'rescale_by_neighborhood'
>>> background_value = -1
>>> rescaled, rescaled_file = rescale_by_neighborhood(input_vtk,
...     indices, nedges, p, set_max_to_1, save_file, output_filestring,
...     background_value)
>>> scalars1, name = read_scalars(input_vtk)
>>> print('{0:0.5f}, {1:0.5f}'.format(max(scalars1), max(rescaled)))
34.95560, 1.00000
>>> print('{0:0.5f}, {1:0.5f}'.format(np.mean(scalars1), np.mean(rescaled)))
7.43822, 0.44950

View rescaled scalar values on surface (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> plot_surfaces(rescaled_file) # doctest: +SKIP

topo_test(index, values, neighbor_lists)¶

Test to see if vertex is a “simple point”.

A simple point is a vertex that when added to or removed from an object (e.g., a curve) on a surface mesh does not alter the object’s topology.

“Simple” is not to be mistaken with the following usage: “A vertex is usually assigned one of five possible classifications: simple, complex, boundary, interior edge, or corner vertex. A simple vertex is surrounded by a closed fan of triangles”.

Parameters:

index (integer) – index of vertex
values (numpy array of integers or floats) – values for all vertices
neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex

Returns:

sp (bool) – simple point or not?
n_inside (integer) – number of neighboring vertices with a value greater than threshold

Examples

>>> # Square with a center vertex:
>>> # indices [[0,1,2],[3,4,6],[7,8,9]] = 0 and indices [2,4,6] = 1:
>>> import numpy as np
>>> from mindboggle.guts.mesh import topo_test
>>> values = np.array([0,0,1,0,1,0,1,0,0])
>>> neighbor_lists = [[1,3],[0,2,3,4],[1,4,5],
...                   [0,1,4,6],[1,2,3,5,6,7],[2,4,7,8],
...                   [3,4,7],[4,5,6,8],[5,7]]
>>> sps = []
>>> for index in range(9):
...     sp, n_inside = topo_test(index, values, neighbor_lists)
...     sps.append(sp)
>>> sps
[False, True, True, True, False, True, True, True, False]

mindboggle.guts.paths module¶

Curves, tracks, skeletons connecting surface mesh vertices.

Authors:

Arno Klein, 2012-2016 (arno@mindboggle.info) http://binarybottle.com
Yrjo Hame, 2012 (yrjo.hame@gmail.com)

connect_points_erosion(S, neighbor_lists, outer_anchors, inner_anchors=[], values=[], erode_ratio=0.1, erode_min_size=10, save_steps=[], save_vtk='', background_value=-1, verbose=False)¶

Connect mesh vertices with a skeleton of 1-vertex-thick curves by erosion.

This algorithm iteratively removes simple topological points and endpoints, optionally in order of lowest to highest values.

Parameters:	S (numpy array of integers) – values for all vertices (disregard background values) outer_anchors (list of integers) – indices of vertices to connect inner_anchors (list of integers) – more vertices to connect; they are removed if they result in endpoints neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex values (numpy array of floats) – values for S elements, to optionally remove points in order of lowest to highest values erode_ratio (float) – fraction of indices to test for removal at each iteration (if values) erode_min_size (integer) – minimum number of vertices when considering erode_ratio save_steps (list of integers (optional)) – iterations at which to save incremental VTK file save_vtk (string) – name of VTK file to transfer incremental values (if save_steps) background_value (integer or float) – background value verbose (bool) – print statements?
Returns:	skeleton – indices to vertices of skeleton
Return type:	list of integers

Examples

>>> # Extract a skeleton to connect endpoints in a fold:
>>> import numpy as np
>>> from mindboggle.guts.paths import connect_points_erosion
>>> from mindboggle.guts.paths import find_outer_endpoints
>>> from mindboggle.mio.vtks import read_scalars, read_vtk
>>> from mindboggle.guts.compute import median_abs_dev
>>> from mindboggle.guts.paths import find_max_values
>>> from mindboggle.guts.mesh import find_neighbors_from_file
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> curv_file = fetch_data(urls['left_mean_curvature'], '', '.vtk')
>>> depth_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> points, f1,f2,f3, curvs, f4,f5,f6 = read_vtk(curv_file, True,True)
>>> depths, name = read_scalars(depth_file, True, True)
>>> folds, name = read_scalars(folds_file, True, True)
>>> values = depths * curvs
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in values[0:5]]
[-0.11778, -0.35642, -0.80759, -0.25654, -0.04411]
>>> neighbor_lists = find_neighbors_from_file(curv_file)
>>> background_value = -1
>>> # Limit number of folds to speed up the test:
>>> limit_folds = True
>>> if limit_folds:
...     fold_numbers = [4] #[4, 6]
...     indices = [i for i,x in enumerate(folds) if x in fold_numbers]
...     i0 = [i for i,x in enumerate(folds) if x not in fold_numbers]
...     folds[i0] = background_value
... else:
...     indices = range(len(values))
>>> # Outer anchors:
>>> min_separation = 10
>>> verbose = False
>>> outer_anchors, tracks = find_outer_endpoints(indices, neighbor_lists,
...                             values, depths, min_separation,
...                             background_value, verbose)
>>> outer_anchors[0:10]
[50324, 66986, 75661]
>>> # Inner anchors:
>>> values0 = [x for x in values if x > 0]
>>> thr = np.median(values0) + 2 * median_abs_dev(values0)
>>> inner_anchors = find_max_values(points, values, min_separation, thr)
>>> inner_anchors[0:10]
[61455, 41761, 67978, 72621, 78546, 40675, 73745, 98736, 125536, 119813]
>>> erode_ratio = 0.10
>>> erode_min_size = 10
>>> save_steps = [] #list(range(0,500,50))
>>> save_vtk = depth_file
>>> S = np.copy(folds)
>>> skeleton = connect_points_erosion(S, neighbor_lists,
...     outer_anchors, inner_anchors, values, erode_ratio, erode_min_size,
...     save_steps, save_vtk, background_value, verbose)
>>> skeleton[0:10]
[50324, 50333, 50339, 51552, 51560, 52707, 52716, 52724, 52725, 53893]

Write out vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> folds[skeleton] = 10 # doctest: +SKIP
>>> folds[outer_anchors] = 15 # doctest: +SKIP
>>> rewrite_scalars(depth_file, 'connect_points_erosion.vtk',
...                 folds, 'skeleton', folds, background_value) # doctest: +SKIP
>>> plot_surfaces('connect_points_erosion.vtk') # doctest: +SKIP

connect_points_hmmf(indices_points, indices, L, neighbor_lists, wN_max=1.0, do_erode=True, background_value=-1, verbose=False)¶

Connect mesh vertices with a skeleton of 1-vertex-thick curves using HMMF.

The goal of this algorithm is to assign each vertex a locally optimal Hidden Markov Measure Field (HMMF) value and to connect vertices according to a cost function that penalizes vertices that do not have high likelihood values and have HMMF values different than their neighbors.

We initialize the HMMF values with likelihood values normalized to the interval (0.5, 1.0] (to guarantee correct topology) and take those values that are greater than the likelihood threshold (1 for each anchor point).

We iteratively update each HMMF value if it is near the likelihood threshold such that a H_step makes it cross the threshold, and the vertex is a “simple point” (its addition/removal alters topology).

Parameters for computing the cost and cost gradients:

wL: weight influence of likelihood on the cost function wN: weight influence of neighbors on the cost function H_step: the amount that the HMMF values are incremented

Parameters:	indices_points (list of integers) – indices of vertices to connect (should contain > 1) indices (list of integers) – indices of vertices through which to connect points L (numpy array of floats) – likelihood values for all vertices in mesh neighbor_lists (list of lists of integers) – indices to neighboring vertices for each vertex in mesh wN_max (float) – maximum neighborhood weight (trust prior more for smoother fundi) do_erode (bool) – erode to create skeleton? background_value (integer) – background value verbose (bool) – print statements?
Returns:	skeleton – indices to vertices connecting the points
Return type:	list of integers

Examples

>>> # Connect vertices according to (usually likelihood) values in a fold:
>>> import numpy as np
>>> from mindboggle.guts.paths import connect_points_hmmf
>>> from mindboggle.guts.paths import find_outer_endpoints
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.guts.compute import median_abs_dev
>>> from mindboggle.guts.mesh import find_neighbors_from_file
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> url1 = urls['left_mean_curvature']
>>> url2 = urls['left_travel_depth']
>>> url3 = urls['left_folds']
>>> curv_file = fetch_data(url1, '', '.vtk')
>>> depth_file = fetch_data(url2, '', '.vtk')
>>> folds_file = fetch_data(url3, '', '.vtk')
>>> curvs, name = read_scalars(curv_file, True, True)
>>> depths, name = read_scalars(depth_file, True, True)
>>> folds, name = read_scalars(folds_file, True, True)
>>> L = curvs * depths
>>> print(np.array_str(L[0:5], precision=5, suppress_small=True))
[-0.11778 -0.35642 -0.80759 -0.25654 -0.04411]
>>> neighbor_lists = find_neighbors_from_file(curv_file)
>>> background_value = -1
>>> # Limit number of folds to speed up the test:
>>> limit_folds = True
>>> if limit_folds:
...     fold_numbers = [4] #[4, 6]
...     indices = [i for i,x in enumerate(folds) if x in fold_numbers]
...     i0 = [i for i,x in enumerate(folds) if x not in fold_numbers]
...     folds[i0] = background_value
... else:
...     indices = range(len(L))
>>> # Outer anchors:
>>> min_separation = 10
>>> verbose = False
>>> indices_points, tracks = find_outer_endpoints(indices, neighbor_lists,
...                             L, depths, min_separation,
...                             background_value, verbose)
>>> wN_max = 2.0
>>> do_erode = True
>>> skeleton = connect_points_hmmf(indices_points, indices, L,
...     neighbor_lists, wN_max, do_erode, background_value, verbose)
>>> skeleton[0:10]
[50324, 50333, 50339, 51552, 51560, 51561, 52708, 52716, 53891, 53892]

Write out vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> folds_copy = np.copy(folds) # doctest: +SKIP
>>> folds[skeleton] = 100 # doctest: +SKIP
>>> folds[indices_points] = 120 # doctest: +SKIP
>>> rewrite_scalars(depth_file, 'connect_points_hmmf.vtk',
...                 folds, 'skeleton', folds, -1) # doctest: +SKIP
>>> plot_surfaces('connect_points_hmmf.vtk') # doctest: +SKIP

find_max_values(points, values, min_separation=10, thr=0.5)¶

Find points with maximal values that are not too close together.

Steps

1. Sort values and find values above the threshold.

2. Initialize special points with the maximum value,
   remove this value, and loop through the remaining high values.

3. If there are no nearby special points,
   assign the maximum value vertex as a special point.

Parameters:	points (numpy array of floats) – coordinates for all vertices values (list (or array) of integers) – values of some kind to maximize over for all vertices min_separation (integer) – minimum number of edges between maximum value vertices thr (float) – value threshold in [0,1]
Returns:	highest – subset of surface mesh vertex indices with the highest values
Return type:	list of integers

Examples

>>> # Extract a skeleton to connect endpoints in a fold:
>>> import numpy as np
>>> from mindboggle.guts.paths import find_max_values
>>> from mindboggle.mio.vtks import read_scalars, read_vtk
>>> from mindboggle.guts.compute import median_abs_dev
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> curv_file = fetch_data(urls['left_mean_curvature'], '', '.vtk')
>>> depth_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> points, f1,f2,f3, curvs, f4,f5,f6 = read_vtk(curv_file, True,True)
>>> depths, name = read_scalars(depth_file, True, True)
>>> values = depths * curvs
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in values[0:5]]
[-0.11778, -0.35642, -0.80759, -0.25654, -0.04411]
>>> min_separation = 10
>>> values0 = [x for x in values if x > 0]
>>> thr = np.median(values0) + 2 * median_abs_dev(values0)
>>> inner_anchors = find_max_values(points, values, min_separation, thr)
>>> inner_anchors[0:10]
[61455, 41761, 67978, 72621, 78546, 40675, 73745, 98736, 125536, 119813]

View anchors in surface fold (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> values[inner_anchors] = np.max(values) + 10 # doctest: +SKIP
>>> rewrite_scalars(depth_file, 'find_max_values.vtk',
...                 values, 'find_max_values', [], -1) # doctest: +SKIP
>>> plot_surfaces('find_max_values.vtk') # doctest: +SKIP

find_outer_endpoints(indices, neighbor_lists, values, values_seeding, min_separation=10, background_value=-1, verbose=False)¶

Find vertices on the boundary of a surface mesh region that are the endpoints to multiple, high-value tracks from the region’s center.

This algorithm propagates multiple tracks from seed vertices at a given depth within a region of a surface mesh to the boundary of the region (via the track_segments() function). The tracks terminate at boundary vertices that can serve as endpoints of fundus curves running along the depths of a fold.

Parameters:

indices (list of integers) – indices of the vertices to segment (such as a fold in a surface mesh)
neighbor_lists (list of lists of integers) – indices to neighboring vertices for each vertex
values (numpy array of floats) – values for all vertices (e.g., fundus likelihood values)
values_seeding (numpy array of floats) – values for all vertices to threshold for initial seeds
min_separation (integer) – minimum number of edges between anchor vertices
background_value (integer or float) – background value
verbose (bool) – print statements?

Returns:

endpoints (list of integers) – indices of surface mesh vertices that are endpoints
endtracks (list of lists of integers) – indices to track vertices

Examples

>>> # Extract a skeleton to connect endpoints in a fold:
>>> import numpy as np
>>> from mindboggle.guts.paths import find_outer_endpoints
>>> from mindboggle.mio.vtks import read_scalars, read_vtk
>>> from mindboggle.guts.compute import median_abs_dev
>>> from mindboggle.guts.mesh import find_neighbors_from_file
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> curv_file = fetch_data(urls['left_mean_curvature'], '', '.vtk')
>>> depth_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> curvs, name = read_scalars(curv_file, True, True)
>>> depths, name = read_scalars(depth_file, True, True)
>>> folds, name = read_scalars(folds_file, True, True)
>>> values = curvs * depths
>>> print(np.array_str(values[0:5], precision=5, suppress_small=True))
[-0.11778 -0.35642 -0.80759 -0.25654 -0.04411]
>>> neighbor_lists = find_neighbors_from_file(curv_file)
>>> background_value = -1
>>> # Limit number of folds to speed up the test:
>>> limit_folds = True
>>> if limit_folds:
...     fold_numbers = [4] #[4, 6]
...     indices = [i for i,x in enumerate(folds) if x in fold_numbers]
...     i0 = [i for i,x in enumerate(folds) if x not in fold_numbers]
...     folds[i0] = background_value
... else:
...     indices = range(len(folds))
>>> # Outer anchors:
>>> min_separation = 10
>>> verbose = False
>>> outer_anchors, tracks = find_outer_endpoints(indices, neighbor_lists,
...                             values, depths, min_separation,
...                             background_value, verbose)
>>> outer_anchors[0:10]
[50324, 66986, 75661]

View anchors in fold on surface (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> folds[outer_anchors] = 100 # doctest: +SKIP
>>> rewrite_scalars(depth_file, 'find_outer_endpoints.vtk',
...                 folds, 'tracks_endpoints_on_folds', folds) # doctest: +SKIP
>>> plot_surfaces('find_outer_endpoints.vtk') # doctest: +SKIP

smooth_skeletons(skeletons, bounds, vtk_file, likelihoods, wN_max=1.0, do_erode=True, save_file=False, output_file='', background_value=-1, verbose=False)¶

Smooth skeleton by dilation followed by connect_points_hmmf().

Steps ::

Segment skeleton into separate sets of connected vertices.
For each skeleton segment, extract endpoints.
Dilate skeleton segment.
Connect endpoints through dilated segment by connect_points_hmmf().
Store smoothed output from #4.

Parameters:

skeletons (list of integers) – skeleton number for each vertex
bounds (list of integers) – region number for each vertex; constrains smoothed skeletons
vtk_file (string) – file from which to extract neighboring vertices for each vertex
likelihoods (list of integers) – fundus likelihood value for each vertex
wN_max (float) – maximum neighborhood weight (trust prior more for smoother skeletons)
do_erode (bool) – erode skeleton?
save_file (bool) – save output VTK file?
output_file (string) – output VTK file
background_value (integer or float) – background value
verbose (bool) – print statements?

Returns:

skeletons (list of integers) – skeleton numbers for all vertices
n_skeletons (integer) – number of skeletons
skeletons_file (string (if save_file)) – name of output VTK file with skeleton numbers

Examples

>>> # Smooth skeleton to extract fundus from one or more folds:
>>> import numpy as np
>>> from mindboggle.mio.plots import plot_surfaces
>>> from mindboggle.guts.paths import smooth_skeletons
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.guts.compute import median_abs_dev
>>> from mindboggle.guts.mesh import find_neighbors_from_file
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> curv_file = fetch_data(urls['left_mean_curvature'], '', '.vtk')
>>> depth_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> fundus_file = fetch_data(urls['left_fundus_per_fold'], '', '.vtk')
>>> curvs, name = read_scalars(curv_file, True, True)
>>> depths, name = read_scalars(depth_file, True, True)
>>> vtk_file = curv_file
>>> likelihoods = depths * curvs
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in likelihoods[0:5]]
[-0.11778, -0.35642, -0.80759, -0.25654, -0.04411]
>>> bounds, name = read_scalars(folds_file, True, True)
>>> skeletons, name = read_scalars(fundus_file, True, True)
>>> background_value = -1
>>> # Limit number of folds to speed up the test:
>>> limit_folds = True
>>> if limit_folds:
...     fold_numbers = [7] #[4, 6]
...     i0 = [i for i,x in enumerate(bounds) if x not in fold_numbers]
...     bounds[i0] = background_value
...     skeletons[i0] = background_value
>>> wN_max = 1.0
>>> do_erode = True
>>> save_file = True
>>> output_file = 'smooth_skeletons.vtk'
>>> verbose = False
>>> smoothed_skeletons, n_skeletons, skel_file = smooth_skeletons(skeletons,
...     bounds, vtk_file, likelihoods, wN_max, do_erode, save_file,
...     output_file, background_value, verbose)
>>> np.where(np.array(smoothed_skeletons)!=-1)[0][0:8]
array([112572, 113435, 113454, 113469, 114294, 114295, 114296, 114312])

Write out vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> iskels = [i for i,x in enumerate(smoothed_skeletons)
...           if x != background_value] # doctest: +SKIP
>>> bounds[iskels] = 100 # doctest: +SKIP
>>> rewrite_scalars(depth_file, 'smooth_skeletons_no_background.vtk',
...                 bounds, 'skeleton', bounds, -1) # doctest: +SKIP
>>> plot_surfaces('smooth_skeletons.vtk') # doctest: +SKIP

track_segments(seed, segments, neighbor_lists, values, sink, background_value=-1)¶

Build a track from a seed vertex through concentric segments of a mesh.

This function builds a track from an initial seed vertex through concentric segments along high-value vertices of a surface mesh optionally terminating at any of a set of sink vertices.

Parameters:	seed (integer) – index to initial seed vertex from which to grow a track segments (list of lists of integers) – indices to vertices for each concentric segment neighbor_lists (list of lists of integers) – indices to neighboring vertices for each vertex values (numpy array of floats) – values for all vertices that help to guide a track sink (list of integers) – indices for vertices that end a track background_value (integer) – background value
Returns:	track – indices of ordered vertices for a single track
Return type:	list of integers

Examples

>>> # Track from deepest point in a fold to its boundary:
>>> import numpy as np
>>> from mindboggle.guts.paths import track_segments
>>> from mindboggle.guts.segment import extract_borders
>>> from mindboggle.guts.segment import segment_rings
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.guts.mesh import find_neighbors_from_file
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> values_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> folds, name = read_scalars(folds_file, True, True)
>>> background_value = -1
>>> # Limit number of folds to speed up the test:
>>> limit_folds = True
>>> if limit_folds:
...     fold_numbers = [4] #[4, 6]
...     indices = [i for i,x in enumerate(folds) if x in fold_numbers]
...     i0 = [i for i,x in enumerate(folds) if x not in fold_numbers]
...     folds[i0] = background_value
... else:
...     indices = range(len(folds))
>>> neighbor_lists = find_neighbors_from_file(values_file)
>>> values, name = read_scalars(values_file, True, True)
>>> seeds = [x for x in indices if values[x] > np.median(values[indices])]
>>> segments = segment_rings(indices, seeds, neighbor_lists, 1, background_value)
>>> #seed = indices[np.argmax(values[indices])]
>>> seed = segments[0][np.argmax(values[segments[0]])]
>>> seed
55181
>>> # Extract boundary:
>>> D = np.ones(len(values))
>>> D[indices] = 2
>>> borders, f1,f2 = extract_borders(list(range(len(values))), D, neighbor_lists)
>>> sink = borders
>>> sink[0:10]
[49083, 49084, 50316, 50317, 50318, 50319, 50324, 50325, 50326, 50327]
>>> track = track_segments(seed, segments, neighbor_lists, values, sink,
...                        background_value)
>>> track[0:10]
[55181, 53892, 52725, 52717, 52708, 51561, 51553, 50340, 50333, 50324]

View track in fold on surface (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> folds[track] = 10 # doctest: +SKIP
>>> folds[seed] = 15 # doctest: +SKIP
>>> rewrite_scalars(values_file, 'track_segments.vtk', folds,
...                 'track', folds) # doctest: +SKIP
>>> plot_surfaces('track_segments.vtk') # doctest: +SKIP

mindboggle.guts.rebound module¶

mindboggle.guts.relabel module¶

Relabel surface or volume or annot files.

Authors:

Arno Klein, 2012-2016 (arno@mindboggle.info) http://binarybottle.com

keep_volume_labels(input_file, labels_to_keep, output_file='', second_file='')¶

Keep only given labels in an image volume (or use to mask second volume).

Parameters:	input_file (string) – labeled nibabel-readable (e.g., nifti) file labels_to_keep (list of integers) – labels to keep output_file (string) – output file name second_file (string) – second nibabel-readable file (keep/erase voxels in this file instead)
Returns:	output_file – output file name
Return type:	string

Examples

>>> # Remove right hemisphere labels
>>> import os
>>> from mindboggle.guts.relabel import keep_volume_labels
>>> from mindboggle.mio.labels import DKTprotocol
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> input_file = fetch_data(urls['freesurfer_labels'], '', '.nii.gz')
>>> second_file = ''
>>> labels_to_keep = list(range(1000, 1036))
>>> output_file = 'keep_volume_labels.nii.gz'
>>> output_file = keep_volume_labels(input_file, labels_to_keep,
...                                  output_file, second_file)

View nifti file (skip test):

>>> from mindboggle.mio.plots import plot_volumes
>>> plot_volumes(output_file) # doctest: +SKIP

overwrite_volume_labels(source, target, output_file='', ignore_labels=[0], erase_labels=True, background_value=-1)¶

For every label in a source image, optionally erase all voxels in the target image with this label (if erase_labels is True), and for every voxel in the source image with this label, assign the label to the corresponding voxel in the target image. The source and target images must have the same volume dimensions.

Parameters:	source (string) – labeled nibabel-readable (e.g., nifti) file target (string) – labeled nibabel-readable (e.g., nifti) file output_file (string) – labeled nibabel-readable (e.g., nifti) file ignore_labels (list) – list of source labels to ignore erase_labels (bool) – erase target labels (that are in source) before overwriting? background_value (integer) – background value (if erase_labels==True)
Returns:	output_file – labeled nibabel-readable (e.g., nifti) file
Return type:	string

Examples

>>> # Overwrite DKT25 with DKT31 labels
>>> import os
>>> from mindboggle.guts.relabel import overwrite_volume_labels
>>> from mindboggle.mio.labels import DKTprotocol
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> source = fetch_data(urls['freesurfer_labels'], '', '.nii.gz')
>>> target = fetch_data(urls['ants_labels'], '', '.nii.gz')
>>> output_file = 'overwrite_volume_labels.nii.gz'
>>> ignore_labels = [0]
>>> erase_labels = False
>>> background_value = -1
>>> output_file = overwrite_volume_labels(source, target, output_file,
...     ignore_labels, erase_labels, background_value)

View nifti file (skip test):

>>> from mindboggle.mio.plots import plot_volumes
>>> plot_volumes(output_file) # doctest: +SKIP

relabel_surface(vtk_file, hemi='', old_labels=[], new_labels=[], erase_remaining=True, erase_labels=[], erase_value=-1, output_file='')¶

Relabel surface in a VTK file.

Parameters:	vtk_file (string) – input labeled VTK file hemi (string) – hemisphere (‘lh’ or ‘rh’ or ‘’) if set, add 1000 to left and 2000 to right hemisphere labels; old_labels (list of integers) – old labels (empty list if labels drawn from vtk scalars); may be used in conjunction with hemi new_labels (list of integers) – new labels (empty list if labels drawn from vtk scalars); may be used in conjunction with hemi erase_remaining (bool) – set all values not in old_labels to erase_value? erase_labels (list of integers) – values to erase (set to erase_value) erase_value (integer) – set vertices with labels in erase_labels to this value output_file (string) – new vtk file name
Returns:	output_file – new vtk file name
Return type:	string

Examples

>>> import numpy as np
>>> from mindboggle.guts.relabel import relabel_surface
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> vtk_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> hemi = 'lh'
>>> old_labels = [1003,1009,1030]
>>> new_labels = [0,500,1000]
>>> erase_remaining = True
>>> erase_labels = [0]
>>> erase_value = -1
>>> output_file = 'relabel_surface.vtk'
>>> output_file = relabel_surface(vtk_file, hemi, old_labels, new_labels,
...     erase_remaining, erase_labels, erase_value, output_file)
>>> labels, name = read_scalars(output_file, True, True)
>>> np.unique(labels)
array([  -1, 1000, 1500, 2000])

View relabeled surface file (skip test):

>>> from mindboggle.mio.plots import plot_surfaces
>>> plot_surfaces(output_file) # doctest: +SKIP

relabel_volume(input_file, old_labels, new_labels, output_file='')¶

Relabel volume labels.

Parameters:	input_file (string) – labeled nibabel-readable (e.g., nifti) file old_labels (list of integers) – old labels new_labels (list of integers) – new labels output_file (string) – output file name
Returns:	output_file – output file name
Return type:	string

Examples

>>> # Convert DKT31 to DKT25 labels
>>> import os
>>> from mindboggle.guts.relabel import relabel_volume
>>> from mindboggle.mio.labels import DKTprotocol
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> input_file = fetch_data(urls['freesurfer_labels'], '', '.nii.gz')
>>> dkt = DKTprotocol()
>>> old_labels = dkt.cerebrum_cortex_numbers + dkt.cerebrum_noncortex_numbers
>>> ctx = [5000 for x in dkt.cerebrum_cortex_numbers]
>>> nonctx = [6000 for x in dkt.cerebrum_noncortex_numbers]
>>> new_labels = ctx + nonctx
>>> output_file = 'relabel_volume.nii.gz'
>>> output_file = relabel_volume(input_file, old_labels, new_labels,
...                              output_file)

View nifti file (skip test):

>>> from mindboggle.mio.plots import plot_volumes # doctest: +SKIP
>>> plot_volumes(output_file) # doctest: +SKIP

remove_volume_labels(input_file, labels_to_remove, output_file='', second_file='')¶

Remove labels from an image volume (or corresponding voxels in a 2nd volume).

Parameters:	input_file (string) – labeled nibabel-readable (e.g., nifti) file labels_to_remove (list of integers) – labels to remove output_file (string) – output file name second_file (string) – second nibabel-readable file (erase voxels in this file instead)
Returns:	output_file – output file name
Return type:	string

Examples

>>> # Remove subcortical labels
>>> import os
>>> from mindboggle.guts.relabel import remove_volume_labels
>>> from mindboggle.mio.labels import DKTprotocol
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> input_file = fetch_data(urls['freesurfer_labels'], '', '.nii.gz')
>>> os.rename(input_file, input_file + '.nii.gz')
>>> input_file = input_file + '.nii.gz'
>>> second_file = ''
>>> labels_to_remove = list(range(1,300)) # Remove noncortex (+aseg) labels
>>> labels_to_remove.extend([1000,1001,2000,2001])
>>> labels_to_remove.extend(list(range(2000,2036))) # Remove right cortex labels
>>> output_file = 'remove_volume_labels.nii.gz'
>>> output_file = remove_volume_labels(input_file, labels_to_remove,
...                                    output_file, second_file)

View nifti file (skip test):

>>> from mindboggle.mio.plots import plot_volumes # doctest: +SKIP
>>> plot_volumes(output_file) # doctest: +SKIP

mindboggle.guts.segment module¶

Functions for segmenting a surface mesh.

Authors:

Yrjo Hame, 2012 (yrjo.hame@gmail.com)
Arno Klein, 2012-2016 (arno@mindboggle.info) http://binarybottle.com

combine_2labels_in_2volumes(file1, file2, label1=3, label2=2, output_file='')¶

Combine the voxels of one label from two files and overlay them on the second label combined from the two files.

An example application is to combine two segmentation volumes, such as from FreeSurfer and ANTs, to obtain a single cortex=2 and noncortex=3 segmentation file, by taking the union of cortex voxels from the segmentations, the union of the noncortex voxels from the segmentations, and overwriting intersecting cortex and noncortex voxels with noncortex (3) labels.

ANTs tends to include more cortical gray matter at the periphery of the brain than Freesurfer, and FreeSurfer tends to include more white matter that extends deep into gyral folds than ANTs, so this function attempts to remedy their differences by overlaying ANTs cortical gray with FreeSurfer white matter.

Example preprocessing steps

1. Run Freesurfer and antsCorticalThickness.sh on T1-weighted image.
2. Convert FreeSurfer volume labels (e.g., wmparc.mgz or aparc+aseg.mgz)
   to cortex (2) and noncortex (3) segments using relabel_volume()
   function [refer to labels.rst or FreeSurferColorLUT labels file].
3. Convert ANTs Atropos-segmented volume (tmpBrainSegmentation.nii.gz)
   to cortex and noncortex segments, by converting 1-labels to 0 and
   4-labels to 3 with the relabel_volume() function
   (the latter is to include deep-gray matter with noncortical tissues).

Parameters:	file1 (string) – name of nibabel-readable file with labels in {label1, label2} file2 (string) – name of nibabel-readable file with labels in {label1, label2} label1 (integer) – source label number label2 (integer) – target label number (to be overwritten by label1 where they intersect) output_file (string) – output file name
Returns:	output_file – name of combined segmented nibabel-readable file with labels in {label1, label2}
Return type:	string

Examples

>>> # Example following documentation above:
>>> import os
>>> from mindboggle.guts.segment import combine_2labels_in_2volumes
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> file1 = fetch_data(urls['freesurfer_segmentation'], '', '.nii.gz')
>>> file2 = fetch_data(urls['ants_segmentation'], '', '.nii.gz')
>>> label1 = 3
>>> label2 = 2
>>> output_file = 'combine_2labels_in_2volumes.nii.gz'
>>> output_file = combine_2labels_in_2volumes(file1, file2, label1,
...                                           label2, output_file)

View nifti file (skip test):

>>> from mindboggle.mio.plots import plot_volumes # doctest: +SKIP
>>> plot_volumes(output_file) # doctest: +SKIP

extract_borders(indices, labels, neighbor_lists, ignore_values=[], return_label_pairs=False)¶

Detect surface label borders in a collection of vertices such as a region.

Label borders are the set of all vertices whose neighbors do not share the same label.

Parameters:

indices (list of integers) – indices to (a subset of) vertices
labels (numpy array of integers) – label numbers for all vertices
neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex
ignore_values (list of integers) – integers to ignore (e.g., background)
return_label_pairs (bool) – return label pairs?

Returns:

border_indices (list of integers) – indices to label boundary vertices
border_label_tuples (list of lists of sorted pairs of integers) – sorted label pairs
unique_border_label_tuples (list of sorted pairs of integers) – unique, sorted label pairs

Examples

>>> # Small example:
>>> from mindboggle.guts.segment import extract_borders
>>> indices = [0,1,2,4,5,8,9]
>>> labels = [10, 20, 30, 40, 50, 60, 70, 80, 90, 100, -1, -1]
>>> neighbor_lists = [[1,2,3], [1,2], [2,3], [2], [4,7], [3,2,3],
...                   [2,3,7,8], [2,6,7], [3,4,8], [7], [7,8], [3,2,3]]
>>> border_indices, border_label_tuples, ubLT = extract_borders(indices,
...     labels, neighbor_lists, [], True)
>>> border_indices
[0, 1, 2, 4, 5, 8]
>>> border_label_tuples
[[20, 30, 40], [20, 30], [30, 40], [50, 80], [30, 40], [40, 50, 90]]
>>> ubLT
[[20, 30, 40], [20, 30], [30, 40], [50, 80], [40, 50, 90]]

Real example – extract sulcus label boundaries:

>>> import numpy as np
>>> from mindboggle.guts.mesh import find_neighbors
>>> from mindboggle.guts.segment import extract_borders
>>> from mindboggle.mio.vtks import read_vtk
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> label_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> f1,f2,f3, faces, labels, f4, npoints, f5 = read_vtk(label_file,
...                                                     True, True)
>>> neighbor_lists = find_neighbors(faces, npoints)
>>> ignore_values = []
>>> return_label_pairs = True
>>> indices_borders, label_pairs, f1 = extract_borders(list(range(npoints)),
...     labels, neighbor_lists, ignore_values, return_label_pairs)
>>> indices_borders[0:10]
[115, 116, 120, 121, 125, 126, 130, 131, 281, 286]
>>> label_pairs[0:5]
[[1005, 1011], [1005, 1011], [1005, 1011], [1005, 1011], [1005, 1011]]

Write borders on surfaces to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars
>>> IDs = -1 * np.ones(npoints) # doctest: +SKIP
>>> IDs[indices_borders] = 1 # doctest: +SKIP
>>> rewrite_scalars(label_file, 'extract_borders.vtk',
...                 IDs, 'borders') # doctest: +SKIP
>>> plot_surfaces('extract_borders.vtk') # doctest: +SKIP

Write just the borders to vtk file and view (skip test):

>>> rewrite_scalars(label_file, 'extract_borders_no_background.vtk',
...                 IDs, 'borders', IDs) # doctest: +SKIP
>>> plot_surfaces('extract_borders_no_background.vtk') # doctest: +SKIP

extract_borders_2nd_surface(labels_file, values_file='', output_file='', background_value=-1)¶

Extract borders (between labels) on a surface. Option: extract border values on a second surface.

Parameters:

labels_file (string) – file name for surface mesh with labels
values_file (string) – file name for surface mesh with values to extract along borders
output_file (string) – full path to output file
background_value (integer or float) – background value

Returns:

border_file (string) – full path to output file
border_values (numpy array) – values for all vertices (background for vertices off label borders)
indices_borders (list of integers) – indices to borders vertices

Examples

>>> # Extract depth values along label borders in sulci (mask):
>>> import numpy as np
>>> from mindboggle.guts.segment import extract_borders_2nd_surface
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> label_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> values_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> background_value = -1
>>> output_file = 'extract_borders_2nd_surface.vtk'
>>> border_file, values, I = extract_borders_2nd_surface(label_file,
...     values_file, output_file, background_value)
>>> [np.float("{0:.{1}f}".format(x, 5)) for x in np.unique(values)[0:8]]
[-1.0, 0.0, 0.00012, 0.00023, 0.00032, 0.00044, 0.00047, 0.00051]
>>> I[0:10]
[115, 116, 120, 121, 125, 126, 130, 131, 281, 286]

Write depth values on label borders to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> plot_surfaces(border_file) # doctest: +SKIP

propagate(points, faces, region, seeds, labels, max_iters=500, tol=0.001, sigma=10, background_value=-1, verbose=False)¶

Propagate labels to segment a surface into contiguous regions, starting from seed vertices.

Imports : mindboggle.guts.rebound

Parameters:	points (array (or list) of lists of three integers) – coordinates for all vertices faces (list of lists of three integers) – indices to three vertices per face (indices start from zero) region (list (or array) of integers) – values > background_value: inclusion in a region for all vertices seeds (numpy array of integers) – seed numbers for all vertices labels (numpy array of integers) – label numbers for all vertices max_iters (integer) – maximum number of iterations to run graph-based learning algorithm tol (float) – threshold to assess convergence of the algorithm sigma (float) – gaussian kernel parameter background_value (integer) – background value verbose (bool) – print statements?
Returns:	segments – region numbers for all vertices
Return type:	numpy array of integers

Examples

>>> # Propagate labels between label boundary segments in a single fold:
>>> import numpy as np
>>> import mindboggle.guts.rebound as rb
>>> from mindboggle.guts.mesh import find_neighbors
>>> from mindboggle.guts.segment import extract_borders
>>> from mindboggle.guts.segment import propagate
>>> from mindboggle.mio.vtks import read_scalars, read_vtk
>>> from mindboggle.mio.labels import DKTprotocol
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> label_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> dkt = DKTprotocol()
>>> folds, name = read_scalars(folds_file, True, True)
>>> points, f1,f2, faces, labels, f3, npoints, f4 = read_vtk(label_file,
...     True, True)
>>> neighbor_lists = find_neighbors(faces, npoints)
>>> background_value = -1
>>> # Limit number of folds to speed up the test:
>>> limit_folds = True
>>> if limit_folds:
...     fold_numbers = [4, 7] #[4, 6]
...     indices_fold = [i for i,x in enumerate(folds) if x in fold_numbers]
...     i0 = [i for i,x in enumerate(folds) if x not in fold_numbers]
...     folds[i0] = background_value
... else:
...     indices_fold = range(len(values))
>>> # Extract the boundary for this fold:
>>> indices_borders, label_pairs, foo = extract_borders(indices_fold,
...     labels, neighbor_lists, [], True)
>>> # Select boundary segments in the sulcus labeling protocol:
>>> seeds = background_value * np.ones(npoints)
>>> for ilist,label_pair_list in enumerate(dkt.sulcus_label_pair_lists):
...     I = [x for i,x in enumerate(indices_borders)
...          if np.sort(label_pairs[i]).tolist() in label_pair_list]
...     seeds[I] = ilist
>>> verbose = False
>>> region = folds
>>> max_iters = 500
>>> tol = 0.001
>>> sigma = 10
>>> segments = propagate(points, faces, region, seeds, labels,
...                      max_iters, tol, sigma, background_value, verbose)
>>> np.unique(segments)[0:10]
array([-1.,  3., 12., 22.])
>>> len_segments = [len(np.where(segments == x)[0])
...                 for x in np.unique(segments) if x != background_value]
>>> len_segments[0:10]
[1152, 388, 116]

Write results to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> rewrite_scalars(label_file, 'propagate.vtk',
...                 segments, 'segments', folds, background_value) # doctest: +SKIP
>>> plot_surfaces('propagate.vtk') # doctest: +SKIP

segment_by_filling_borders(regions, neighbor_lists, background_value=-1, verbose=False)¶

Fill borders (contours) on a surface mesh to segment vertices into contiguous regions.

Steps ::

Extract region borders (assumed to be closed contours)
Segment borders into separate, contiguous borders
For each boundary
1. Find the neighbors to either side of the boundary
2. Segment the neighbors into exterior and interior sets of neighbors
3. Find the interior (smaller) sets of neighbors
4. Fill the contours formed by the interior neighbors

Parameters:	regions (numpy array of integers) – region numbers for all vertices neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex background_value (integer or float) – background value verbose (bool) – print statements?
Returns:	segments – region numbers for all vertices
Return type:	numpy array of integers

Examples

>>> # Segment folds by extracting their borders and filling them in separately:
>>> import numpy as np
>>> from mindboggle.guts.segment import segment_by_filling_borders
>>> from mindboggle.guts.mesh import find_neighbors
>>> from mindboggle.mio.vtks import read_vtk
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> depth_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> f1,f2,f3, faces, depths, f4, npoints, input_vtk = read_vtk(depth_file,
...                                                            True, True)
>>> background_value = -1
>>> regions = background_value * np.ones(npoints)
>>> regions[depths > 0.50] = 1
>>> neighbor_lists = find_neighbors(faces, npoints)
>>> verbose = False
>>> segments = segment_by_filling_borders(regions, neighbor_lists,
...                                       background_value, verbose)
>>> len(np.unique(segments))
19
>>> len_segments = []
>>> for useg in np.unique(segments):
...     len_segments.append(len(np.where(segments == useg)[0]))
>>> len_segments[0:10]
[19446, 8619, 13846, 23, 244, 101687, 16, 792, 210, 76]

Write results to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> rewrite_scalars(depth_file, 'segment_by_filling_borders.vtk',
...                 segments, 'segments', [], -1) # doctest: +SKIP
>>> plot_surfaces('segment_by_filling_borders.vtk') # doctest: +SKIP

segment_by_region(data, regions=[], surface_file='', save_file=False, output_file='', background_value=-1, verbose=False)¶

Segment data on one surface by regions on a corresponding surface.

For example use sulcus regions to segment fundus data.

Parameters:

data (list of integers) – numbers for all vertices (including background_value)
regions (numpy array or list of integers) – number for each vertex, used to filter and label data
surface_file (string (if save_file)) – VTK file
save_file (bool) – save output VTK file?
output_file (string) – output VTK file
background_value (integer or float) – background value
verbose (bool) – print statements?

Returns:

segment_per_region (list of integers) – region numbers for all vertices (or background_value)
n_segments (integer) – number of segments
segment_per_region_file (string (if save_file)) – output VTK file with region numbers for data vertices

Examples

>>> # Segment fundus data by sulcus regions:
>>> import numpy as np
>>> from mindboggle.guts.segment import segment_by_region
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> fundus_file = fetch_data(urls['left_fundus_per_fold'], '', '.vtk')
>>> data, name = read_scalars(fundus_file, True, True)
>>> surface_file = fetch_data(urls['left_sulci'], '', '.vtk')
>>> regions, name = read_scalars(surface_file, True, True)
>>> save_file = True
>>> output_file = 'segment_by_region.vtk'
>>> background_value = -1
>>> verbose = False
>>> o1, o2, segment_per_region_file = segment_by_region(data, regions,
...     surface_file, save_file, output_file, background_value, verbose)
>>> segment_numbers = [x for x in np.unique(o1) if x != background_value]
>>> lens = []
>>> for segment_number in segment_numbers:
...     lens.append(len([x for x in o1 if x == segment_number]))
>>> lens[0:10]
[333, 229, 419, 251, 281, 380, 162, 137, 256, 29]

View result (skip test):

>>> from mindboggle.mio.plots import plot_surfaces
>>> plot_surfaces(segment_per_region_file) # doctest: +SKIP

segment_regions(vertices_to_segment, neighbor_lists, min_region_size=1, seed_lists=[], keep_seeding=False, spread_within_labels=False, labels=[], label_lists=[], values=[], max_steps='', background_value=-1, verbose=False)¶

Segment vertices of surface into contiguous regions by seed growing, starting from zero or more lists of seed vertices.

Parameters:	vertices_to_segment (list of integers) – indices to mesh vertices to be segmented neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex min_region_size (integer) – minimum size of segmented set of vertices seed_lists (list of lists, or empty list) – each list contains indices to seed vertices to segment vertices_to_segment keep_seeding (bool) – grow from new seeds even after all seed lists have fully grown spread_within_labels (bool) – grow seeds only by vertices with labels in the seed labels? labels (list of integers (required only if spread_within_labels)) – label numbers for all vertices label_lists (list of lists of integers (required only if spread_within_labels)) – List of unique labels for each seed list to grow into (If empty, set to unique labels for each seed list) values (list of floats (default empty)) – values for all vertices for use in preferentially directed segmentation (segment in direction of lower values) max_steps (integer (or empty string for infinity)) – maximum number of segmentation steps to take for each seed list background_value (integer or float) – background value verbose (bool) – print statements?
Returns:	segments – region numbers for all vertices
Return type:	numpy array of integers

Examples

>>> # Segment deep regions with or without seeds:
>>> import numpy as np
>>> from mindboggle.guts.segment import segment_regions
>>> from mindboggle.mio.vtks import read_vtk
>>> from mindboggle.guts.mesh import find_neighbors
>>> from mindboggle.mio.fetch_data import prep_tests
>>> background_value = -1
>>> urls, fetch_data = prep_tests()
>>> depth_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> f1,f2,f3, faces, depths, f4, npoints, t5 = read_vtk(depth_file,
...                                                     True, True)
>>> vertices_to_segment = np.where(depths > 0.50)[0].tolist()  # (sped up)
>>> neighbor_lists = find_neighbors(faces, npoints)

Example 1: without seed lists

>>> segments = segment_regions(vertices_to_segment, neighbor_lists)
>>> len(np.unique(segments))
92
>>> len_segments = [len(np.where(segments == x)[0])
...                 for x in np.unique(segments) if x != background_value]
>>> len_segments[0:10]
[110928, 4, 1399, 1274, 5, 139, 255, 12, 5, 1686]

Write results to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> rewrite_scalars(depth_file, 'segment_regions_no_seeds.vtk', segments,
...                 'segments', [], -1) # doctest: +SKIP
>>> plot_surfaces('segment_regions_no_seeds.vtk') # doctest: +SKIP

Example 2: with seed lists

>>> from mindboggle.guts.segment import extract_borders
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.mio.labels import DKTprotocol
>>> dkt = DKTprotocol()
>>> label_lists = [np.unique(np.ravel(x)) for x in dkt.sulcus_label_pair_lists]
>>> labels_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> labels, name = read_scalars(labels_file)
>>> indices_borders, label_pairs, foo = extract_borders(vertices_to_segment,
...     labels, neighbor_lists, ignore_values=[], return_label_pairs=True)
>>> seed_lists = []
>>> for label_pair_list in dkt.sulcus_label_pair_lists:
...     seed_lists.append([x for i,x in enumerate(indices_borders)
...                        if np.sort(label_pairs[i]).tolist()
...                        in label_pair_list])
>>> keep_seeding = True
>>> spread_within_labels = True
>>> values = []
>>> max_steps = ''
>>> background_value = -1
>>> verbose = False
>>> segments = segment_regions(vertices_to_segment,
...     neighbor_lists, 1, seed_lists, keep_seeding, spread_within_labels,
...     labels, label_lists, values, max_steps, background_value, verbose)
>>> len(np.unique(segments))
122
>>> np.unique(segments)[0:10]
array([-1.,  1.,  3.,  4.,  5.,  6.,  7.,  8.,  9., 11.])
>>> len_segments = [len(np.where(segments == x)[0])
...                 for x in np.unique(segments) if x != background_value]
>>> len_segments[0:10]
[6962, 8033, 5965, 5598, 7412, 3636, 3070, 5244, 3972, 6144]

Write results to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> rewrite_scalars(depth_file, 'segment_regions.vtk', segments,
...     'segments_from_seeds', [], -1) # doctest: +SKIP
>>> plot_surfaces('segment_regions.vtk') # doctest: +SKIP

segment_rings(region, seeds, neighbor_lists, step=1, background_value=-1, verbose=False)¶

Iteratively segment a region of surface mesh as concentric segments.

Parameters:	region (list of integers) – indices of region vertices to segment (such as a fold) seeds (list of integers) – indices of seed vertices neighbor_lists (list of lists of integers) – indices to neighboring vertices for each vertex step (integer) – number of segmentation steps before assessing segments background_value (integer) – background value verbose (bool) – print statements?
Returns:	segments – indices to vertices for each concentric segment
Return type:	list of lists of integers

Examples

>>> import numpy as np
>>> from mindboggle.mio.vtks import read_scalars
>>> from mindboggle.guts.mesh import find_neighbors_from_file
>>> from mindboggle.guts.segment import extract_borders
>>> from mindboggle.guts.segment import segment_rings
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> vtk_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> folds_file = fetch_data(urls['left_folds'], '', '.vtk')
>>> values, name = read_scalars(vtk_file, True, True)
>>> neighbor_lists = find_neighbors_from_file(vtk_file)
>>> background_value = -1
>>> fold, name = read_scalars(folds_file)
>>> indices = [i for i,x in enumerate(fold) if x != background_value]
>>> # Initialize seeds with the boundary of thresholded indices:
>>> use_threshold = True
>>> if use_threshold:
...     # Threshold at the median depth or within maximum values in boundary:
...     threshold = np.median(values[indices]) #+ np.std(values[indices])
...     indices_high = [x for x in indices if values[x] >= threshold]
...     # Make sure threshold is within the maximum values of the boundary:
...     B = np.ones(len(values))
...     B[indices] = 2
...     borders, foo1, foo2 = extract_borders(list(range(len(B))), B, neighbor_lists)
...     borders = [x for x in borders if values[x] != background_value]
...     if list(frozenset(indices_high).intersection(borders)):
...         threshold = np.max(values[borders]) + np.std(values[borders])
...         indices_high = [x for x in indices if values[x] >= threshold]
...     # Extract threshold boundary vertices as seeds:
...     B = background_value * np.ones(len(values))
...     B[indices_high] = 2
...     seeds, foo1, foo2 = extract_borders(list(range(len(values))), B, neighbor_lists)
... # Or initialize P with the maximum value point:
... else:
...     seeds = [indices[np.argmax(values[indices])]]
...     indices_high = []
>>> indices = list(frozenset(indices).difference(indices_high))
>>> indices = list(frozenset(indices).difference(seeds))
>>> step = 1
>>> verbose = False
>>> segments = segment_rings(indices, seeds, neighbor_lists, step,
...                          background_value, verbose)
>>> len(segments)
56
>>> [len(x) for x in segments][0:10]
[5540, 5849, 6138, 5997, 4883, 3021, 1809, 1165, 842, 661]
>>> segments[0][0:10]
[65539, 65540, 98308, 98316, 131112, 131121, 131122, 131171, 131175, 131185]

Write results to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import read_scalars, rewrite_scalars # doctest: +SKIP
>>> S = background_value * np.ones(len(values)) # doctest: +SKIP
>>> for i, segment in enumerate(segments): S[segment] = i # doctest: +SKIP
>>> rewrite_scalars(vtk_file, 'segment_rings.vtk', S, 'segment_rings',
...                 [], -1) # doctest: +SKIP
>>> plot_surfaces('segment_rings.vtk') # doctest: +SKIP

select_largest(points, faces, exclude_labels=[-1], areas=None, reindex=True, background_value=-1, verbose=False)¶

Select the largest segment (connected set of indices) in a surface mesh.

In case a surface patch is fragmented, we select the largest fragment, remove extraneous triangular faces, and reindex indices.

Parameters:

points (list of lists of 3 floats) – x,y,z coordinates for each vertex of the structure
faces (list of lists of 3 integers) – 3 indices to vertices that form a triangle on the mesh
exclude_labels (list of integers) – background values to exclude
areas (numpy array or list of floats (or None)) – surface area scalar values for all vertices
reindex (bool) – reindex indices in faces?
background_value (integer or float) – background value
verbose (bool) – print statements?

Returns:

points (list of lists of 3 floats) – x,y,z coordinates for each vertex of the structure
faces (list of lists of 3 integers) – 3 indices to vertices that form a triangle on the mesh

Examples

>>> # Two surface patches:
>>> import numpy as np
>>> from mindboggle.mio.vtks import read_scalars, read_vtk
>>> from mindboggle.guts.mesh import keep_faces
>>> from mindboggle.guts.segment import select_largest
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> label_file = fetch_data(urls['left_freesurfer_labels'], '', '.vtk')
>>> area_file = fetch_data(urls['left_area'], '', '.vtk')
>>> exclude_labels = [-1]
>>> points, indices, f1, faces, labels, f2,f3,f4 = read_vtk(label_file,
...                                                         True, True)
>>> I28 = [i for i,x in enumerate(labels) if x==1028] # superior frontal
>>> I20 = [i for i,x in enumerate(labels) if x==1020] # pars triangularis
>>> I28.extend(I20)
>>> faces = keep_faces(faces, I28)
>>> areas, u1 = read_scalars(area_file, True, True)
>>> reindex = True
>>> background_value = -1
>>> verbose = False
>>> points2, faces2 = select_largest(points, faces, exclude_labels, areas,
...                                  reindex, background_value, verbose)
>>> points2[0]
[-1.4894376993179321, 53.260337829589844, 56.523414611816406]
>>> points2[1]
[-2.2537832260131836, 53.045711517333984, 56.23670959472656]
>>> points2[2]
[-13.091879844665527, 56.41604232788086, 59.330955505371094]
>>> faces2[0]
[7704, 7306, 7703]
>>> faces2[1]
[7694, 7704, 7705]
>>> faces2[2]
[7703, 8119, 7704]

Write two surfaces to vtk file and view (skip test):

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import write_vtk # doctest: +SKIP
>>> scalars = np.zeros(np.shape(labels)) # doctest: +SKIP
>>> scalars[I28] = 1 # doctest: +SKIP
>>> segments_file = 'select_largest_two_labels.vtk' # doctest: +SKIP
>>> write_vtk(vtk_file, points, indices, [], faces) # doctest: +SKIP
>>> plot_surfaces(vtk_file) # doctest: +SKIP

Write larger surface to vtk file and view (skip test):

>>> vtk_file = 'select_largest.vtk' # doctest: +SKIP
>>> write_vtk(vtk_file, points2, list(range(len(points2))), [], faces2) # doctest: +SKIP
>>> plot_surfaces(vtk_file) # doctest: +SKIP

split_brain(image_file, label_file, left_labels, right_labels)¶

Split a brain using left/right labels.

Parameters:

image_file (string) – nibabel-readable image volume
label_file (string) – nibabel-readable labeled brain image volume
left_labels (list of integers) – left label numbers
right_labels (list of integers) – right label numbers

Returns:

left_brain (string) – name of nibabel-readable image volume with left half of brain image
right_brain (string) – name of nibabel-readable image volume with right half of brain image

Examples

>>> from mindboggle.guts.segment import split_brain
>>> from mindboggle.mio.labels import DKTprotocol
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> image_file = fetch_data(urls['freesurfer_segmentation'], '', '.nii.gz')
>>> label_file = fetch_data(urls['freesurfer_labels'], '', '.nii.gz')
>>> dkt = DKTprotocol()
>>> left_labels = dkt.left_cerebrum_numbers
>>> right_labels = dkt.right_cerebrum_numbers
>>> left_brain, right_brain = split_brain(image_file, label_file,
...                                       left_labels, right_labels)

Write results to nifti file and view (skip test):

>>> from mindboggle.mio.plots import plot_volumes # doctest: +SKIP
>>> plot_volumes([left_brain, right_brain]) # doctest: +SKIP

watershed(depths, points, indices, neighbor_lists, min_size=1, depth_factor=0.25, depth_ratio=0.1, tolerance=0.01, regrow=True, background_value=-1, verbose=False)¶

Segment surface vertices into contiguous regions by a watershed algorithm.

Segment surface mesh into contiguous “watershed basins” by seed growing from an iterative selection of the deepest vertices.

Steps

1. Grow segments from an iterative selection of the deepest seeds.
2. Regrow segments from the resulting seeds, until each seed's
    segment touches a boundary.
3. Use the segment() function to fill in the rest.
4. Merge segments if their seeds are too close to each other
    or their depths are very different.

Note

Despite the above precautions, the order of seed selection in segment()
could possibly influence the resulting borders between adjoining
segments (vs. propagate(), which is slower and insensitive to depth,
but is not biased by seed order).

In the example below, white spots indicate incomplete segmentation.

Parameters:

depths (numpy array of floats) – depth values for all vertices
points (list of lists of floats) – each element is a list of 3-D coordinates of a vertex on a surface mesh
indices (list of integers) – indices to mesh vertices to be segmented
min_size (index) – the minimum number of vertices in a basin
neighbor_lists (list of lists of integers) – each list contains indices to neighboring vertices for each vertex
depth_factor (float) – factor to determine whether to merge two neighboring watershed catchment basins – they are merged if the Euclidean distance between their basin seeds is less than this fraction of the maximum Euclidean distance between points having minimum and maximum depths
depth_ratio (float) – the minimum fraction of depth for a neighboring shallower watershed catchment basin (otherwise merged with the deeper basin)
tolerance (float) – tolerance for detecting differences in depth between vertices
regrow (bool) – regrow segments from watershed seeds?
background_value (integer or float) – background value
verbose (bool) – print statements?

Returns:

segments (list of integers) – region numbers for all vertices
seed_indices (list of integers) – list of indices to seed vertices

Examples

>>> # Perform watershed segmentation on the deeper portions of a surface:
>>> import numpy as np
>>> from mindboggle.guts.mesh import find_neighbors
>>> from mindboggle.guts.segment import watershed
>>> from mindboggle.mio.vtks import read_vtk
>>> from mindboggle.mio.fetch_data import prep_tests
>>> urls, fetch_data = prep_tests()
>>> depth_file = fetch_data(urls['left_travel_depth'], '', '.vtk')
>>> points, indices, f1, faces, depths, f2, npoints, f3 = read_vtk(depth_file,
...     return_first=True, return_array=True)
>>> indices = np.where(depths > 0.01)[0]  # high to speed up
>>> neighbor_lists = find_neighbors(faces, npoints)
>>> min_size = 50
>>> depth_factor = 0.25
>>> depth_ratio = 0.1
>>> tolerance = 0.01
>>> regrow = True
>>> background_value = -1
>>> verbose = False
>>> segments, seed_indices = watershed(depths, points,
...     indices, neighbor_lists, min_size, depth_factor, depth_ratio,
...     tolerance, regrow, background_value, verbose)
>>> len(np.unique(segments))
202
>>> len_segments = []
>>> for useg in np.unique(segments):
...     len_segments.append(len(np.where(segments == useg)[0]))
>>> len_segments[0:10]
[2976, 4092, 597, 1338, 1419, 1200, 1641, 220, 1423, 182]

Write watershed segments and seeds to vtk file and view (skip test). Note: white spots indicate incomplete segmentation:

>>> from mindboggle.mio.plots import plot_surfaces # doctest: +SKIP
>>> from mindboggle.mio.vtks import rewrite_scalars # doctest: +SKIP
>>> segments_array = np.array(segments)
>>> segments_array[seed_indices] = 1.5 * np.max(segments_array)
>>> rewrite_scalars(depth_file, 'watershed.vtk',
...                 segments_array, 'segments', [], -1) # doctest: +SKIP
>>> plot_surfaces('watershed.vtk') # doctest: +SKIP

mindboggle.guts.utilities module¶

Utility functions.

Authors:

Arno Klein, 2012-2016 (arno@mindboggle.info) http://binarybottle.com

execute(cmd, type='os')¶

Execute command by either subprocess.call or os.system.

Parameters:	cmd (sequence (string also permitted if type=='os')) – command with arguments type (string) – how to execute {os, subprocess}

Examples

>>> from mindboggle.guts.utilities import execute # doctest: +SKIP
>>> cmd = ['date', '-r', '0'] # doctest: +SKIP
>>> type = 'subprocess' # doctest: +SKIP
>>> execute(cmd, type) # doctest: +SKIP
Wed Dec 31 19:00:00 EST 1969
>>> type = 'os' # doctest: +SKIP
>>> execute(cmd, type) # doctest: +SKIP
Wed Dec 31 19:00:00 EST 1969
>>> cmd = 'date -r 0' # doctest: +SKIP
>>> execute(cmd) # doctest: +SKIP
Wed Dec 31 19:00:00 EST 1969

list_strings(string1='', string2='', string3='', string4='')¶

Put strings in a list.

Parameters:	string1 (string) – string2 (string) – string3 (string) – string4 (string) –
Returns:	string_list
Return type:	list of strings

Examples

>>> from mindboggle.guts.utilities import list_strings
>>> string1 = 'a b c'
>>> string2 = 'd e f'
>>> string3 = ''
>>> string4 = 'j k l'
>>> string_list = list_strings(string1, string2, string3, string4)
>>> string_list
['a b c', 'd e f', 'j k l']

mindboggle.guts package¶

Submodules¶

mindboggle.guts.compute module¶

mindboggle.guts.graph module¶

mindboggle.guts.kernels module¶

mindboggle.guts.mesh module¶

mindboggle.guts.paths module¶

mindboggle.guts.rebound module¶

mindboggle.guts.relabel module¶

mindboggle.guts.segment module¶

mindboggle.guts.utilities module¶

Module contents¶

Mindboggle

Navigation

Related Topics