diff options
Diffstat (limited to 'src/cython/cython')
-rw-r--r-- | src/cython/cython/nerve_gic.pyx | 30 | ||||
-rw-r--r--[-rwxr-xr-x] | src/cython/cython/persistence_graphical_tools.py | 393 | ||||
-rw-r--r-- | src/cython/cython/simplex_tree.pyx | 2 | ||||
-rw-r--r-- | src/cython/cython/subsampling.pyx | 8 |
4 files changed, 293 insertions, 140 deletions
diff --git a/src/cython/cython/nerve_gic.pyx b/src/cython/cython/nerve_gic.pyx index 30a14d3b..5f01b379 100644 --- a/src/cython/cython/nerve_gic.pyx +++ b/src/cython/cython/nerve_gic.pyx @@ -38,14 +38,14 @@ cdef extern from "Nerve_gic_interface.h" namespace "Gudhi": double compute_distance_from_confidence_level(double alpha) void compute_distribution(int N) double compute_p_value() - void compute_PD() + vector[pair[double, double]] compute_PD() void find_simplices() void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree) bool read_point_cloud(string off_file_name) double set_automatic_resolution() void set_color_from_coordinate(int k) void set_color_from_file(string color_file_name) - void set_color_from_vector(vector[double] color) + void set_color_from_range(vector[double] color) void set_cover_from_file(string cover_file_name) void set_cover_from_function() void set_cover_from_Euclidean_Voronoi(int m) @@ -67,6 +67,8 @@ cdef extern from "Nerve_gic_interface.h" namespace "Gudhi": void write_info() void plot_DOT() void plot_OFF() + void set_point_cloud_from_range(vector[vector[double]] cloud) + void set_distances_from_range(vector[vector[double]] distance_matrix) # CoverComplex python interface cdef class CoverComplex: @@ -102,6 +104,22 @@ cdef class CoverComplex: """ return self.thisptr != NULL + def set_point_cloud_from_range(self, cloud): + """ Reads and stores the input point cloud from a vector stored in memory. + + :param cloud: Input vector containing the point cloud. + :type cloud: vector[vector[double]] + """ + return self.thisptr.set_point_cloud_from_range(cloud) + + def set_distances_from_range(self, distance_matrix): + """ Reads and stores the input distance matrix from a vector stored in memory. + + :param distance_matrix: Input vector containing the distance matrix. + :type distance_matrix: vector[vector[double]] + """ + return self.thisptr.set_distances_from_range(distance_matrix) + def compute_confidence_level_from_distance(self, distance): """Computes the confidence level of a specific bottleneck distance threshold. @@ -145,7 +163,7 @@ cdef class CoverComplex: def compute_PD(self): """Computes the extended persistence diagram of the complex. """ - self.thisptr.compute_PD() + return self.thisptr.compute_PD() def create_simplex_tree(self): """ @@ -162,7 +180,7 @@ cdef class CoverComplex: self.thisptr.find_simplices() def read_point_cloud(self, off_file): - """Reads and stores the input point cloud. + """Reads and stores the input point cloud from .(n)OFF file. :param off_file: Name of the input .OFF or .nOFF file. :type off_file: string @@ -206,14 +224,14 @@ cdef class CoverComplex: else: print("file " + color_file_name + " not found.") - def set_color_from_vector(self, color): + def set_color_from_range(self, color): """Computes the function used to color the nodes of the simplicial complex from a vector stored in memory. :param color: Input vector of values. :type color: vector[double] """ - self.thisptr.set_color_from_vector(color) + self.thisptr.set_color_from_range(color) def set_cover_from_file(self, cover_file_name): """Creates the cover C from a file containing the cover elements of diff --git a/src/cython/cython/persistence_graphical_tools.py b/src/cython/cython/persistence_graphical_tools.py index 216ab8d6..d7be936f 100755..100644 --- a/src/cython/cython/persistence_graphical_tools.py +++ b/src/cython/cython/persistence_graphical_tools.py @@ -1,8 +1,3 @@ -import matplotlib.pyplot as plt -import matplotlib.patches as mpatches -import numpy as np -import os - """This file is part of the Gudhi Library. The Gudhi library (Geometric Understanding in Higher Dimensions) is a generic C++ library for computational topology. @@ -59,157 +54,295 @@ Only 13 colors for the palette palette = ['#ff0000', '#00ff00', '#0000ff', '#00ffff', '#ff00ff', '#ffff00', '#000000', '#880000', '#008800', '#000088', '#888800', '#880088', '#008888'] - + def plot_persistence_barcode(persistence=[], persistence_file='', alpha=0.6, - max_barcodes=1000, inf_delta=0.1, legend=False): + max_intervals=1000, max_barcodes=1000, + inf_delta=0.1, legend=False): """This function plots the persistence bar code from persistence values list or from a :doc:`persistence file <fileformats>`. - :param persistence: Persistence values list. + :param persistence: Persistence intervals values list grouped by dimension. :type persistence: list of tuples(dimension, tuple(birth, death)). :param persistence_file: A :doc:`persistence file <fileformats>` style name (reset persistence if both are set). :type persistence_file: string - :param alpha: barcode transparency value (0.0 transparent through 1.0 opaque - default is 0.6). + :param alpha: barcode transparency value (0.0 transparent through 1.0 + opaque - default is 0.6). :type alpha: float. - :param max_barcodes: number of maximal barcodes to be displayed. - Set it to 0 to see all, Default value is 1000. - (persistence will be sorted by life time if max_barcodes is set) - :type max_barcodes: int. - :param inf_delta: Infinity is placed at ((max_death - min_birth) x inf_delta). - A reasonable value is between 0.05 and 0.5 - default is 0.1. + :param max_intervals: maximal number of intervals to display. + Selected intervals are those with the longest life time. Set it + to 0 to see all. Default value is 1000. + :type max_intervals: int. + :param inf_delta: Infinity is placed at :code:`((max_death - min_birth) x + inf_delta)` above :code:`max_death` value. A reasonable value is + between 0.05 and 0.5 - default is 0.1. :type inf_delta: float. + :param legend: Display the dimension color legend (default is False). + :type legend: boolean. :returns: A matplotlib object containing horizontal bar plot of persistence (launch `show()` method on it to display it). """ - if persistence_file is not '': - if os.path.isfile(persistence_file): - # Reset persistence - persistence = [] - diag = read_persistence_intervals_grouped_by_dimension(persistence_file=persistence_file) - for key in diag.keys(): - for persistence_interval in diag[key]: - persistence.append((key, persistence_interval)) - else: - print("file " + persistence_file + " not found.") - return None - - if max_barcodes > 0 and max_barcodes < len(persistence): - # Sort by life time, then takes only the max_plots elements - persistence = sorted(persistence, key=lambda life_time: life_time[1][1]-life_time[1][0], reverse=True)[:max_barcodes] - - persistence = sorted(persistence, key=lambda birth: birth[1][0]) - - (min_birth, max_death) = __min_birth_max_death(persistence) - ind = 0 - delta = ((max_death - min_birth) * inf_delta) - # Replace infinity values with max_death + delta for bar code to be more - # readable - infinity = max_death + delta - axis_start = min_birth - delta - # Draw horizontal bars in loop - for interval in reversed(persistence): - if float(interval[1][1]) != float('inf'): - # Finite death case - plt.barh(ind, (interval[1][1] - interval[1][0]), height=0.8, - left = interval[1][0], alpha=alpha, - color = palette[interval[0]], - linewidth=0) - else: - # Infinite death case for diagram to be nicer - plt.barh(ind, (infinity - interval[1][0]), height=0.8, - left = interval[1][0], alpha=alpha, - color = palette[interval[0]], - linewidth=0) - ind = ind + 1 - - if legend: - dimensions = list(set(item[0] for item in persistence)) - plt.legend(handles=[mpatches.Patch(color=palette[dim], - label=str(dim)) for dim in dimensions], - loc='lower right') - plt.title('Persistence barcode') - # Ends plot on infinity value and starts a little bit before min_birth - plt.axis([axis_start, infinity, 0, ind]) - return plt + try: + import matplotlib.pyplot as plt + import matplotlib.patches as mpatches + import numpy as np + import os + + if persistence_file is not '': + if os.path.isfile(persistence_file): + # Reset persistence + persistence = [] + diag = read_persistence_intervals_grouped_by_dimension(persistence_file=persistence_file) + for key in diag.keys(): + for persistence_interval in diag[key]: + persistence.append((key, persistence_interval)) + else: + print("file " + persistence_file + " not found.") + return None + + if max_barcodes is not 1000: + print('Deprecated parameter. It has been replaced by max_intervals') + max_intervals = max_barcodes + + if max_intervals > 0 and max_intervals < len(persistence): + # Sort by life time, then takes only the max_intervals elements + persistence = sorted(persistence, key=lambda life_time: life_time[1][1]-life_time[1][0], reverse=True)[:max_intervals] + + persistence = sorted(persistence, key=lambda birth: birth[1][0]) + + (min_birth, max_death) = __min_birth_max_death(persistence) + ind = 0 + delta = ((max_death - min_birth) * inf_delta) + # Replace infinity values with max_death + delta for bar code to be more + # readable + infinity = max_death + delta + axis_start = min_birth - delta + # Draw horizontal bars in loop + for interval in reversed(persistence): + if float(interval[1][1]) != float('inf'): + # Finite death case + plt.barh(ind, (interval[1][1] - interval[1][0]), height=0.8, + left = interval[1][0], alpha=alpha, + color = palette[interval[0]], + linewidth=0) + else: + # Infinite death case for diagram to be nicer + plt.barh(ind, (infinity - interval[1][0]), height=0.8, + left = interval[1][0], alpha=alpha, + color = palette[interval[0]], + linewidth=0) + ind = ind + 1 + + if legend: + dimensions = list(set(item[0] for item in persistence)) + plt.legend(handles=[mpatches.Patch(color=palette[dim], + label=str(dim)) for dim in dimensions], + loc='lower right') + plt.title('Persistence barcode') + # Ends plot on infinity value and starts a little bit before min_birth + plt.axis([axis_start, infinity, 0, ind]) + return plt + + except ImportError: + print("This function is not available, you may be missing numpy and/or matplotlib.") def plot_persistence_diagram(persistence=[], persistence_file='', alpha=0.6, - band=0., max_plots=1000, inf_delta=0.1, legend=False): - """This function plots the persistence diagram from persistence values list - or from a :doc:`persistence file <fileformats>`. + band=0., max_intervals=1000, max_plots=1000, inf_delta=0.1, legend=False): + """This function plots the persistence diagram from persistence values + list or from a :doc:`persistence file <fileformats>`. - :param persistence: Persistence values list. + :param persistence: Persistence intervals values list grouped by dimension. :type persistence: list of tuples(dimension, tuple(birth, death)). :param persistence_file: A :doc:`persistence file <fileformats>` style name (reset persistence if both are set). :type persistence_file: string - :param alpha: plot transparency value (0.0 transparent through 1.0 opaque - default is 0.6). + :param alpha: plot transparency value (0.0 transparent through 1.0 + opaque - default is 0.6). :type alpha: float. :param band: band (not displayed if :math:`\leq` 0. - default is 0.) :type band: float. - :param max_plots: number of maximal plots to be displayed - Set it to 0 to see all, Default value is 1000. - (persistence will be sorted by life time if max_plots is set) - :type max_plots: int. - :param inf_delta: Infinity is placed at ((max_death - min_birth) x inf_delta). - A reasonable value is between 0.05 and 0.5 - default is 0.1. + :param max_intervals: maximal number of intervals to display. + Selected intervals are those with the longest life time. Set it + to 0 to see all. Default value is 1000. + :type max_intervals: int. + :param inf_delta: Infinity is placed at :code:`((max_death - min_birth) x + inf_delta)` above :code:`max_death` value. A reasonable value is + between 0.05 and 0.5 - default is 0.1. :type inf_delta: float. + :param legend: Display the dimension color legend (default is False). + :type legend: boolean. :returns: A matplotlib object containing diagram plot of persistence (launch `show()` method on it to display it). """ - if persistence_file is not '': - if os.path.isfile(persistence_file): - # Reset persistence - persistence = [] - diag = read_persistence_intervals_grouped_by_dimension(persistence_file=persistence_file) - for key in diag.keys(): - for persistence_interval in diag[key]: - persistence.append((key, persistence_interval)) - else: - print("file " + persistence_file + " not found.") - return None - - if max_plots > 0 and max_plots < len(persistence): - # Sort by life time, then takes only the max_plots elements - persistence = sorted(persistence, key=lambda life_time: life_time[1][1]-life_time[1][0], reverse=True)[:max_plots] - - (min_birth, max_death) = __min_birth_max_death(persistence, band) - ind = 0 - delta = ((max_death - min_birth) * inf_delta) - # Replace infinity values with max_death + delta for diagram to be more - # readable - infinity = max_death + delta - axis_start = min_birth - delta - - # line display of equation : birth = death - x = np.linspace(axis_start, infinity, 1000) - # infinity line and text - plt.plot(x, x, color='k', linewidth=1.0) - plt.plot(x, [infinity] * len(x), linewidth=1.0, color='k', alpha=alpha) - plt.text(axis_start, infinity, r'$\infty$', color='k', alpha=alpha) - # bootstrap band - if band > 0.: - plt.fill_between(x, x, x+band, alpha=alpha, facecolor='red') + try: + import matplotlib.pyplot as plt + import matplotlib.patches as mpatches + import numpy as np + import os - # Draw points in loop - for interval in reversed(persistence): - if float(interval[1][1]) != float('inf'): - # Finite death case - plt.scatter(interval[1][0], interval[1][1], alpha=alpha, - color = palette[interval[0]]) + if persistence_file is not '': + if os.path.isfile(persistence_file): + # Reset persistence + persistence = [] + diag = read_persistence_intervals_grouped_by_dimension(persistence_file=persistence_file) + for key in diag.keys(): + for persistence_interval in diag[key]: + persistence.append((key, persistence_interval)) + else: + print("file " + persistence_file + " not found.") + return None + + if max_plots is not 1000: + print('Deprecated parameter. It has been replaced by max_intervals') + max_intervals = max_plots + + if max_intervals > 0 and max_intervals < len(persistence): + # Sort by life time, then takes only the max_intervals elements + persistence = sorted(persistence, key=lambda life_time: life_time[1][1]-life_time[1][0], reverse=True)[:max_intervals] + + (min_birth, max_death) = __min_birth_max_death(persistence, band) + delta = ((max_death - min_birth) * inf_delta) + # Replace infinity values with max_death + delta for diagram to be more + # readable + infinity = max_death + delta + axis_start = min_birth - delta + + # line display of equation : birth = death + x = np.linspace(axis_start, infinity, 1000) + # infinity line and text + plt.plot(x, x, color='k', linewidth=1.0) + plt.plot(x, [infinity] * len(x), linewidth=1.0, color='k', alpha=alpha) + plt.text(axis_start, infinity, r'$\infty$', color='k', alpha=alpha) + # bootstrap band + if band > 0.: + plt.fill_between(x, x, x+band, alpha=alpha, facecolor='red') + + # Draw points in loop + for interval in reversed(persistence): + if float(interval[1][1]) != float('inf'): + # Finite death case + plt.scatter(interval[1][0], interval[1][1], alpha=alpha, + color = palette[interval[0]]) + else: + # Infinite death case for diagram to be nicer + plt.scatter(interval[1][0], infinity, alpha=alpha, + color = palette[interval[0]]) + + if legend: + dimensions = list(set(item[0] for item in persistence)) + plt.legend(handles=[mpatches.Patch(color=palette[dim], label=str(dim)) for dim in dimensions]) + + plt.title('Persistence diagram') + plt.xlabel('Birth') + plt.ylabel('Death') + # Ends plot on infinity value and starts a little bit before min_birth + plt.axis([axis_start, infinity, axis_start, infinity + delta]) + return plt + + except ImportError: + print("This function is not available, you may be missing numpy and/or matplotlib.") + +def plot_persistence_density(persistence=[], persistence_file='', + nbins=300, bw_method=None, + max_intervals=1000, dimension=None, + cmap=None, legend=False): + """This function plots the persistence density from persistence + values list or from a :doc:`persistence file <fileformats>`. Be + aware that this function does not distinguish the dimension, it is + up to you to select the required one. This function also does not handle + degenerate data set (scipy correlation matrix inversion can fail). + + :param persistence: Persistence intervals values list grouped by dimension. + :type persistence: list of tuples(dimension, tuple(birth, death)). + :param persistence_file: A :doc:`persistence file <fileformats>` + style name (reset persistence if both are set). + :type persistence_file: string + :param nbins: Evaluate a gaussian kde on a regular grid of nbins x + nbins over data extents (default is 300) + :type nbins: int. + :param bw_method: The method used to calculate the estimator + bandwidth. This can be 'scott', 'silverman', a scalar constant + or a callable. If a scalar, this will be used directly as + kde.factor. If a callable, it should take a gaussian_kde + instance as only parameter and return a scalar. If None + (default), 'scott' is used. See + `scipy.stats.gaussian_kde documentation + <http://scipy.github.io/devdocs/generated/scipy.stats.gaussian_kde.html>`_ + for more details. + :type bw_method: str, scalar or callable, optional. + :param max_intervals: maximal number of points used in the density + estimation. + Selected intervals are those with the longest life time. Set it + to 0 to see all. Default value is 1000. + :type max_intervals: int. + :param dimension: the dimension to be selected in the intervals + (default is None to mix all dimensions). + :type dimension: int. + :param cmap: A matplotlib colormap (default is + matplotlib.pyplot.cm.hot_r). + :type cmap: cf. matplotlib colormap. + :param legend: Display the color bar values (default is False). + :type legend: boolean. + :returns: A matplotlib object containing diagram plot of persistence + (launch `show()` method on it to display it). + """ + try: + import matplotlib.pyplot as plt + import numpy as np + from scipy.stats import kde + import os + import math + + if persistence_file is not '': + if os.path.isfile(persistence_file): + # Reset persistence + persistence = [] + diag = read_persistence_intervals_grouped_by_dimension(persistence_file=persistence_file) + for key in diag.keys(): + for persistence_interval in diag[key]: + persistence.append((key, persistence_interval)) + else: + print("file " + persistence_file + " not found.") + return None + + persistence_dim = [] + if dimension is not None: + persistence_dim = [(dim_interval) for dim_interval in persistence if (dim_interval[0] == dimension)] else: - # Infinite death case for diagram to be nicer - plt.scatter(interval[1][0], infinity, alpha=alpha, - color = palette[interval[0]]) - ind = ind + 1 - - if legend: - dimensions = list(set(item[0] for item in persistence)) - plt.legend(handles=[mpatches.Patch(color=palette[dim], label=str(dim)) for dim in dimensions]) - - plt.title('Persistence diagram') - plt.xlabel('Birth') - plt.ylabel('Death') - # Ends plot on infinity value and starts a little bit before min_birth - plt.axis([axis_start, infinity, axis_start, infinity + delta]) - return plt + persistence_dim = persistence + + if max_intervals > 0 and max_intervals < len(persistence_dim): + # Sort by life time, then takes only the max_intervals elements + persistence_dim = sorted(persistence_dim, + key=lambda life_time: life_time[1][1]-life_time[1][0], + reverse=True)[:max_intervals] + + # Set as numpy array birth and death (remove undefined values - inf and NaN) + birth = np.asarray([(interval[1][0]) for interval in persistence_dim if (math.isfinite(interval[1][1]) and math.isfinite(interval[1][0]))]) + death = np.asarray([(interval[1][1]) for interval in persistence_dim if (math.isfinite(interval[1][1]) and math.isfinite(interval[1][0]))]) + + # line display of equation : birth = death + x = np.linspace(death.min(), birth.max(), 1000) + plt.plot(x, x, color='k', linewidth=1.0) + + # Evaluate a gaussian kde on a regular grid of nbins x nbins over data extents + k = kde.gaussian_kde([birth,death], bw_method=bw_method) + xi, yi = np.mgrid[birth.min():birth.max():nbins*1j, death.min():death.max():nbins*1j] + zi = k(np.vstack([xi.flatten(), yi.flatten()])) + + # default cmap value cannot be done at argument definition level as matplotlib is not yet defined. + if cmap is None: + cmap = plt.cm.hot_r + # Make the plot + plt.pcolormesh(xi, yi, zi.reshape(xi.shape), cmap=cmap) + + if legend: + plt.colorbar() + + plt.title('Persistence density') + plt.xlabel('Birth') + plt.ylabel('Death') + return plt + + except ImportError: + print("This function is not available, you may be missing numpy, matplotlib and/or scipy.") diff --git a/src/cython/cython/simplex_tree.pyx b/src/cython/cython/simplex_tree.pyx index e302486b..8397d9d9 100644 --- a/src/cython/cython/simplex_tree.pyx +++ b/src/cython/cython/simplex_tree.pyx @@ -515,7 +515,7 @@ cdef class SimplexTree: :returns: The persistence intervals. :rtype: list of pair of list of int - :note: intervals_in_dim function requires + :note: persistence_pairs function requires :func:`persistence()<gudhi.SimplexTree.persistence>` function to be launched first. """ diff --git a/src/cython/cython/subsampling.pyx b/src/cython/cython/subsampling.pyx index ac09b7a3..e9d61a37 100644 --- a/src/cython/cython/subsampling.pyx +++ b/src/cython/cython/subsampling.pyx @@ -112,7 +112,8 @@ def pick_n_random_points(points=None, off_file='', nb_points=0): return subsampling_n_random_points(points, nb_points) def sparsify_point_set(points=None, off_file='', min_squared_dist=0.0): - """Subsample a point set by picking random vertices. + """Outputs a subset of the input points so that the squared distance + between any two points is greater than or equal to min_squared_dist. :param points: The input point set. :type points: vector[vector[double]]. @@ -122,8 +123,9 @@ def sparsify_point_set(points=None, off_file='', min_squared_dist=0.0): :param off_file: An OFF file style name. :type off_file: string - :param min_squared_dist: Number of points of the subsample. - :type min_squared_dist: unsigned. + :param min_squared_dist: Minimum squared distance separating the output \ + points. + :type min_squared_dist: float. :returns: The subsample point set. :rtype: vector[vector[double]] """ |