API Reference
Index
Classes
Object that performs circular coordinates via persistent cohomology of sparse filtrations (Jose Perea 2020). |
|
Object that performs sparse toroidal coordinates via persistent cohomology as in (L. |
|
Object that performs multiscale real projective coordinates via persistent cohomology of sparse filtrations (Jose Perea 2018). |
|
Object that performs multiscale complex projective coordinates via persistent cohomology of sparse filtrations (Perea 2018). |
|
Object that performs multiscale lens coordinates via persistent cohomology of sparse filtrations (Polanco, Perea 2019). |
|
Partitions of unity subordinate to open ball covers using standard bump functions. |
|
Finite samples from topologically nontrivial spaces. |
|
Plotting utilities for DREiMac's output. |
|
Utilities for adding, rotating, and plotting circle-valued maps. |
|
Utilities for manipulating projective space-valued maps. |
|
Utilities for manipulating lens space-valued maps. |
|
Utilities for subsampling from point clouds and distance matrices. |
dreimac.CircularCoords methods
Get circular coordinates. |
dreimac.ToroidalCoords methods
Get toroidal coordinates. |
dreimac.ProjectiveCoords methods
Get real projective coordinates. |
dreimac.ComplexProjectiveCoords methods
Get complex projective coordinates. |
dreimac.LensCoords methods
Get lens coordinates. |
dreimac.PartUnity methods
Linear partition of unity. |
|
Quadratic partition of unity. |
|
Exponential partition of unity. |
dreimac.GeometryExamples methods
Sample a set of line segments, as witnessed by square patches |
|
Sample a set of (smoothened) dots on the plane, as witnessed by square patches |
|
Return points sampled on a 3D torus |
|
Return samples on a klein bottle in 4D |
|
Return samples on a genus two surface in 3D |
|
Samples on a trefoil in 3D |
|
Samples on three concentric noisy circles in 2D. |
|
Samples on a 2-sphere in 3D. |
|
Samples on a circle in 2D. |
|
Distance matrix of a sample of the Moore space M(1,Z/prime). |
dreimac.PlotUtils methods
Plot patches in a best fitting rectangular grid |
|
Plot the boundary of RP2 on the unit circle |
|
Plot the boundary of two sphere hemispheres |
|
Draw a spherical mesh in 3D around the unit sphere. |
|
Make axes of 3D plot have equal scale so that spheres appear as spheres, cubes as cubes, etc. |
dreimac.CircleMapUtils methods
Rotationally offset a circle-valued map. |
|
Given k circle-valued maps on a dataset with n points and an l x k matrix with integer coefficients, return the l linear combinations of the given circle-valued maps induced by the given matrix. |
|
Given a circle map construct an array of the same shape that can be fed as color in a matplotlib scatterplot to simulate the sinebow colormap. |
|
Given points on the circle and a number of levelsets subdivide the circle into the given number of levelsets and return a smoothened membership function to the levelsets. |
dreimac.ProjectiveMapUtils methods
Do a projective stereographic projection |
dreimac.LensMapUtils methods
Project points on the unit sphere in C^2 representing points in the 3-dimensional lens space corresponding to the prime q to points on the unit disk in R^3. |
dreimac.GeometryUtils methods
A Naive O(NM) algorithm to do furthest points sampling, assuming the input is a point cloud specified in Euclidean space. |
|
Get a minmax sample of Euclidean point cloud using an approximation of the geodesic distance. |
Details
- class dreimac.CircularCoords(X, n_landmarks, distance_matrix=False, prime=41, maxdim=1, verbose=False)
Object that performs circular coordinates via persistent cohomology of sparse filtrations (Jose Perea 2020).
Parameters
- X: ndarray(N, d)
A point cloud with N points in d dimensions
- n_landmarks: int
Number of landmarks to use
- distance_matrix: boolean
If true, treat X as a distance matrix instead of a point cloud
- primeint
Field coefficient with which to compute rips on landmarks
- maxdimint
Maximum dimension of homology. Only dimension 1 is needed for circular coordinates, but it may be of interest to see other dimensions (e.g. for a torus)
- verbosebool
Print debug information.
- get_coordinates(perc=0.5, cocycle_idx=0, partunity_fn=<function PartUnity.linear>, standard_range=True, check_cocycle_condition=True)
Get circular coordinates.
Parameters
- percfloat
Percent coverage. Must be between 0 and 1.
- cocycle_idxinteger
Integer representing the index of the persistent cohomology class used to construct the Eilenberg-MacLane coordinate. Persistent cohomology classes are ordered by persistence, from largest to smallest.
- partunity_fn(dist_land_data, r_cover) -> phi
A function from the distances of each landmark to a bump function
- standard_rangebool
Whether to use the parameter perc to choose a filtration parameter that guarantees that the selected cohomology class represents a non-trivial class in the Cech complex.
- check_cocycle_conditionbool
Whether to check, and fix if necessary, that the integer cocycle constructed using finite field coefficients satisfies the cocycle condition.
Returns
- thetasndarray(N)
Circular coordinates
- class dreimac.ToroidalCoords(X, n_landmarks, distance_matrix=False, prime=41, maxdim=1, verbose=False)
Object that performs sparse toroidal coordinates via persistent cohomology as in (L. Scoccola, H. Gakhar, J. Bush, N. Schonsheck, T. Rask, L. Zhou, J. Perea 2022)
Parameters
- X: ndarray(N, d)
A point cloud with N points in d dimensions
- n_landmarks: int
Number of landmarks to use
- distance_matrix: boolean
If true, treat X as a distance matrix instead of a point cloud
- primeint
Field coefficient with which to compute rips on landmarks
- maxdimint
Maximum dimension of homology. Only dimension 1 is needed for circular coordinates, but it may be of interest to see other dimensions (e.g. for a torus)
- get_coordinates(perc=0.5, cocycle_idxs=[0], partunity_fn=<function PartUnity.linear>, standard_range=True, check_cocycle_condition=True)
Get toroidal coordinates.
Parameters
- percfloat
Percent coverage. Must be between 0 and 1.
- cocycle_idxinteger
Integer representing the index of the persistent cohomology class used to construct the Eilenberg-MacLane coordinate. Persistent cohomology classes are ordered by persistence, from largest to smallest.
- partunity_fn(dist_land_data, r_cover) -> phi
A function from the distances of each landmark to a bump function
- standard_rangebool
Whether to use the parameter perc to choose a filtration parameter that guarantees that the selected cohomology class represents a non-trivial class in the Cech complex.
- check_cocycle_conditionbool
Whether to check, and fix if necessary, that the integer cocycle constructed using finite field coefficients satisfies the cocycle condition.
Returns
- thetasndarray(n, N)
List of circular coordinates, with n the length of cocycle_idxs
- class dreimac.ProjectiveCoords(X, n_landmarks, distance_matrix=False, maxdim=1, verbose=False)
Object that performs multiscale real projective coordinates via persistent cohomology of sparse filtrations (Jose Perea 2018).
Parameters
- X: ndarray(N, d)
A point cloud with N points in d dimensions
- n_landmarks: int
Number of landmarks to use
- distance_matrix: boolean
If true, treat X as a distance matrix instead of a point cloud
- maxdimint
Maximum dimension of homology. Only dimension 1 is needed for real projective coordinates, but it may be of interest to see other dimensions (e.g. for a torus)
- partunity_fn: ndarray(n_landmarks, N) -> ndarray(n_landmarks, N)
A partition of unity function
- get_coordinates(perc=0.9, cocycle_idx=0, proj_dim=2, partunity_fn=<function PartUnity.linear>, standard_range=True, projective_dim_red_mode='exponential')
Get real projective coordinates.
Parameters
- percfloat
Percent coverage. Must be between 0 and 1.
- cocycle_idxlist
Integer representing the index of the persistent cohomology class used to construct the Eilenberg-MacLane coordinate. Persistent cohomology classes are ordered by persistence, from largest to smallest.
- proj_diminteger
Dimension down to which to project the data.
- partunity_fn: (dist_land_data, r_cover) -> phi
A function from the distances of each landmark to a bump function
- standard_rangebool
Whether to use the parameter perc to choose a filtration parameter that guarantees that the selected cohomology class represents a class in the Cech complex.
- projective_dim_red_modestring
Either “one-by-one”, “exponential”, or “direct”. How to perform equivariant dimensionality reduction. “exponential” usually works best, being fast without compromising quality.
Returns
- ndarray(N, proj_dim+1)
The projective coordinates
- class dreimac.ComplexProjectiveCoords(X, n_landmarks, distance_matrix=False, prime=41, maxdim=2, verbose=False)
Object that performs multiscale complex projective coordinates via persistent cohomology of sparse filtrations (Perea 2018).
Parameters
- X: ndarray(N, d)
A point cloud with N points in d dimensions
- n_landmarks: int
Number of landmarks to use
- distance_matrix: boolean
If true, treat X as a distance matrix instead of a point cloud
- primeint
Field coefficient with which to compute rips on landmarks
- maxdimint
Maximum dimension of homology. Only dimension 2 is needed for complex projective coordinates, but it may be of interest to see other dimensions
- get_coordinates(perc=0.5, cocycle_idx=0, proj_dim=1, partunity_fn=<function PartUnity.linear>, standard_range=True, check_cocycle_condition=True, projective_dim_red_mode='exponential')
Get complex projective coordinates.
Parameters
- percfloat
Percent coverage. Must be between 0 and 1.
- cocycle_idxinteger
Integer representing the index of the persistent cohomology class used to construct the Eilenberg-MacLane coordinate. Persistent cohomology classes are ordered by persistence, from largest to smallest.
- proj_diminteger
Complex dimension down to which to project the data.
- partunity_fn(dist_land_data, r_cover) -> phi
A function from the distances of each landmark to a bump function
- standard_rangebool
Whether to use the parameter perc to choose a filtration parameter that guarantees that the selected cohomology class represents a non-trivial class in the Cech complex.
- check_cocycle_conditionbool
Whether to check, and fix if necessary, that the integer cocycle constructed using finite field coefficients satisfies the cocycle condition.
- projective_dim_red_modestring
Either “one-by-one”, “exponential”, or “direct”. How to perform equivariant dimensionality reduction. “exponential” usually works best, being fast without compromising quality.
Returns
- ndarray(N, proj_dim+1)
Complex projective coordinates
- class dreimac.LensCoords(X, n_landmarks, distance_matrix=False, prime=3, maxdim=1, verbose=False)
Object that performs multiscale lens coordinates via persistent cohomology of sparse filtrations (Polanco, Perea 2019).
Parameters
- X: ndarray(N, d)
A point cloud with N points in d dimensions
- n_landmarks: int
Number of landmarks to use
- distance_matrix: boolean
If true, treat X as a distance matrix instead of a point cloud
- primeint
Field coefficient with which to compute rips on landmarks
- maxdimint
Maximum dimension of homology. Only dimension 1 is needed for lens coordinates, but it may be of interest to see other dimensions (e.g. for a torus)
- partunity_fn: ndarray(n_landmarks, N) -> ndarray(n_landmarks, N)
A partition of unity function
- get_coordinates(perc=0.9, cocycle_idx=0, complex_dim=2, partunity_fn=<function PartUnity.linear>, standard_range=True, projective_dim_red_mode='exponential')
Get lens coordinates.
Parameters
- percfloat
Percent coverage. Must be between 0 and 1.
- cocycle_idxlist
Add the cocycles together, sorted from most to least persistent.
- complex_diminteger
Describes the dimension of the lens space. The lens space is obtained by taking a quotient of the unit sphere inside the complex plane to the power of complex_dimension.
- partunity_fn: (dist_land_data, r_cover) -> phi
A function from the distances of each landmark to a bump function.
- standard_rangebool
Whether to use the parameter perc to choose a filtration parameter that guarantees that the selected cohomology class represents a class in the Cech complex.
- projective_dim_red_modestring
Either “one-by-one”, “exponential”, or “direct”. How to perform equivariant dimensionality reduction. “exponential” usually works best, being fast without compromising quality.
Returns
- ndarray(N, complex_dim-1)
The lens coordinates
- class dreimac.PartUnity
Partitions of unity subordinate to open ball covers using standard bump functions.
- static exp(ds, r_cover)
Exponential partition of unity.
Parameters
- ds: ndarray(n)
Some subset of distances between landmarks and data points
- r_cover: float
Covering radius
Returns
- varphi: ndarray(n)
The bump function
- class dreimac.GeometryExamples
Finite samples from topologically nontrivial spaces.
- static bullseye()
Samples on three concentric noisy circles in 2D.
Returns
- X: ndarray(n_samples, 2)
2D circles samples
- static genus_two_surface()
Return samples on a genus two surface in 3D
Returns
- X: ndarray(n_samples, 3)
3D surface samples
- static klein_bottle_4d(n_samples, R, r, seed=None)
Return samples on a klein bottle in 4D
Parameters
- n_samplesint
Number of random samples on the projective plane
- R: float
Outer radius
- r: float
Inner radius
- seed: int
Seed to use. If omitted, use the number of samples as a seed
Returns
- X: ndarray(n_samples, 4)
4D klein bottle samples
- static line_patches(dim, n_angles, n_offsets, sigma)
Sample a set of line segments, as witnessed by square patches
Parameters
- dim: int
Patches will be dim x dim
- n_angles: int
Number of angles to sweep between 0 and pi
- n_offsets: int
Number of offsets to sweep from the origin to the edge of the patch
- sigma: float
The blur parameter. Higher sigma is more blur
Returns
- ndarray(n_angles*n_offsets, dim*dim)
An array of all of the patches raveled into dim*dim dimensional Euclidean space
- static moore_space_distance_matrix(rough_n_points=2000, prime=3)
Distance matrix of a sample of the Moore space M(1,Z/prime).
Parameters
- rough_n_samplesint
The function returns a distance matrix with approximately n = rough_n_points * (pi/4) points.
- primeint
The prime associated with the Moore space.
Returns
- ndarray(n, n)
Distance matrix of points on the Moore space M(1,Z/prime).
- static moving_dot(sqrt_num_images, sigma=3, dim=10)
Sample a set of (smoothened) dots on the plane, as witnessed by square patches
Parameters
- sqrt_num_imagesint
The output will consist of sqrt_num_images**2 patches.
- sigmafloat
The blur parameter. Higher sigma is more blur. Default is 3.
- dimint
Patches will be dim x dim. Default is 10.
Returns
- ndarray(sqrt_num_images*sqrt_num_images, dim*dim)
An array of all of the patches raveled into dim*dim dimensional Euclidean space
- static noisy_circle(n_samples, noise_size=0.2, seed=0)
Samples on a circle in 2D.
Parameters
- n_samplesint
Number of points on the circle to sample.
- noise_sizefloat
Maximum perturbation per sample. Must be between 0 and 1.
- seedfloat
Seed for random number generator.
Returns
- X: ndarray(n_samples, 2)
2D circle samples
- static rp2_metric(n_samples, seed=None)
Return a distance matrix of points on the projective plane obtained by identifying antipodal Gaussian random samples of a sphere
Parameters
- n_samplesint
Number of random samples on the projective plane
- seed: int
Seed to use. If omitted, use the number of samples as a seed
Returns
- ndarray(n_samples, 3)
Original points on the sphere
- ndarray(n_samples, n_samples)
Distance matrix of rp2
- static sphere(n_samples)
Samples on a 2-sphere in 3D.
Returns
- X: ndarray(n_samples, 3)
3D sphere samples
- static three_circles()
Samples on two circles of different radii in 2D.
Returns
- X: ndarray(n_samples, 2)
2D circles samples
- static torus_3d(n_samples, R, r, seed=None)
Return points sampled on a 3D torus
Parameters
- n_samplesint
Number of random samples on the torus
- R: float
Outer radius
- r: float
Inner radius
- seed: int
Seed to use. If omitted, use the number of samples as a seed
Returns
- X: ndarray(n_samples, 4)
3D torus samples
- class dreimac.PlotUtils
Plotting utilities for DREiMac’s output.
- static imscatter(X, P, dim, zoom=1, ax=None)
Plot patches in specified locations in R2
Parameters
- Xndarray (N, 2)
The positions of each patch in R2
- Pndarray (N, dim*dim)
An array of all of the patches
- dimint
The dimension of each patch
- axmatplotlib axes, optional
If given, plot on those axes, otherwise plot on current axes by calling gca()
- static plot_2sphere_boundary(ax=None)
Plot the boundary of two sphere hemispheres
Parameters
- axmatplotlib axes, optional
If given, plot on those axes, otherwise plot on current axes by calling gca()
- static plot_3sphere_mesh(n_parallels=10, n_meridians=20, alpha=0.1, ax=None)
Draw a spherical mesh in 3D around the unit sphere.
Parameters
- n_parallelsint
How many parallels to draw.
- n_meridiansint
How many meridians to draw.
- alphafloat
Opacity of the mesh. Must be between 0 and 1.
- axmatplotlib axes, optional
If given, plot on those axes, otherwise plot on current axes by calling gca()
- static plot_patches(P, zoom=1, ax=None)
Plot patches in a best fitting rectangular grid
Parameters
- axmatplotlib axes, optional
If given, plot on those axes, otherwise plot on current axes by calling gca()
- static plot_proj_boundary(ax=None)
Plot the boundary of RP2 on the unit circle
Parameters
- axmatplotlib axes, optional
If given, plot on those axes, otherwise plot on current axes by calling gca()
- static plot_rp2_stereo(S, f, arrowcolor='c', facecolor=(0.15, 0.15, 0.15))
Plot a 2D Stereographic Projection
Parameters
- Sndarray (N, 2)
An Nx2 array of N points to plot on RP2
- fndarray (N) or ndarray (N, 3)
A function with which to color the points, or a list of colors
- static plot_rp3_stereo(ax, S, f, draw_sphere=False)
Plot a 3D Stereographic Projection
Parameters
- axmatplotlib axis
3D subplotting axis
- Sndarray (N, 3)
An Nx3 array of N points to plot on RP3
- fndarray (N) or ndarray (N, 3)
A function with which to color the points, or a list of colors
- draw_sphereboolean
Whether to draw the 2-sphere
- class dreimac.CircleMapUtils
Utilities for adding, rotating, and plotting circle-valued maps.
- static center(circle_map)
Rotationally offset a circle-valued map so that most of the points map to the center of the circle (i.e., pi).
Parameters
- circle_map: ndarray
A numpy array of numbers between 0 and 2pi representing points on the circle.
Returns
- ndarray
A numpy array of numbers between 0 and 2pi representing the rotation of the given points by the given offset.
- static levelset_coloring(circle_map, n_levelsets=4, smoothing=0.25)
Given points on the circle and a number of levelsets subdivide the circle into the given number of levelsets and return a smoothened membership function to the levelsets. This is useful for coloring a dataset X according to a circle-valued map X -> S^1.
Parameters
- circle_map: ndarray
A numpy array of numbers between 0 and 2pi representing points on the circle.
- n_levelset: int, optional, default is 4
Number of levelsets to evenly cover the circle.
- smoothing: float, optional, default is 0.25
How much to smoothen the membership function
Returns
- ndarray
The smoothened membership function of each of the given points in the circle.
- static linear_combination(circle_maps, linear_combination_matrix)
Given k circle-valued maps on a dataset with n points and an l x k matrix with integer coefficients, return the l linear combinations of the given circle-valued maps induced by the given matrix.
Parameters
- circle_maps: list or ndarray(k, n, dtype=float)
A numpy array with rows containing n points in the circle represented as floats between 0 and 2pi.
- linear_combination_matrix: list or ndarray(l, k, dtype=int)
A numpy array encoding l integer linear combinations of the given k circle-valued maps.
Returns
- ndarray(l, n, dtype=float)
A numpy array with rows containing n points in the circle representing the l linear combinations of the given k circle-valued maps.
- static offset(circle_map, offset)
Rotationally offset a circle-valued map.
points on the circle.
- offset: float
A number between 0 and 1 representing a rotational offset.
Returns
- ndarray
A numpy array of numbers between 0 and 2pi representing the rotation of the given points by the given offset.
- static to_sinebow(circle_map)
Given a circle map construct an array of the same shape that can be fed as color in a matplotlib scatterplot to simulate the sinebow colormap.
Parameters
- circle_map: ndarray
A numpy array of numbers between 0 and 2pi representing points on the circle.
Returns
- ndarray
A numpy array of floats to be used as color in a matplotlib scatterplot.
- class dreimac.ProjectiveMapUtils
Utilities for manipulating projective space-valued maps.
- static circle_to_3dnorthpole(x)
Convert a point selected on the circle to a 3D unit vector on the upper hemisphere
Parameters
- x: ndarray(2)
Selected point in the circle
Returns
- x: ndarray(2)
Selected point in the circle (possibly clipped to the circle)
- u: ndarray(3)
Unit vector on the upper hemisphere
- static get_stereo_proj_codim1(pX, u=array([], dtype=float64))
Do a projective stereographic projection
Parameters
- pX: ndarray(N, d)
A collection of N points in d dimensions
- u: ndarray(d)
A unit vector representing the north pole for the stereographic projection
Returns
- S: ndarray(N, d-1)
The stereographically projected coordinates
- static hopf_map(X)
Use the Hopf map to project points on the unit sphere in C^2 representing points in CP^1 to points on the unit sphere in R^3.
Parameters
- Xcomplex ndarray(n,2)
Points on the unit sphere in C^2.
Returns
- real ndarray(n,3)
Points on the unit sphere in R^3.
- static rotmat(a, b=array([], dtype=float64))
Construct a d x d rotation matrix that rotates a vector a so that it coincides with a vector b
Parameters
- andarray (d)
A d-dimensional vector that should be rotated to b
- bndarray(d)
A d-dimensional vector that should end up residing at the north pole (0,0,…,0,1)
- class dreimac.LensMapUtils
Utilities for manipulating lens space-valued maps.
- static lens_3D_to_disk_3D(X, q)
Project points on the unit sphere in C^2 representing points in the 3-dimensional lens space corresponding to the prime q to points on the unit disk in R^3.
Parameters
- Xcomplex ndarray(n,2)
Points on the unit sphere in C^2.
- qint
Prime number such that X represents points in the 3-dimensional lens space corresponding to this prime.
Returns
- real ndarray(n,3)
Points on the unit disk in R^3.
- class dreimac.GeometryUtils
Utilities for subsampling from point clouds and distance matrices.
- static get_csm(X, Y)
Return the Euclidean cross-similarity matrix between the M points in the Mxd matrix X and the N points in the Nxd matrix Y.
Parameters
- Xndarray (M, d)
A matrix holding the coordinates of M points
- Yndarray (N, d)
A matrix holding the coordinates of N points
Returns
- Dndarray (M, N)
An MxN Euclidean cross-similarity matrix
- static get_csm_projarc(X, Y)
Return the projective arc length cross-similarity between two point clouds specified as points on the sphere
Parameters
- Xndarray (M, d)
A matrix holding the coordinates of M points on RP^{d-1}
- Yndarray (N, d)
A matrix holding the coordinates of N points on RP^{d-1}
Returns
- Dndarray (M, N)
An MxN cross-similarity matrix
- static get_greedy_perm_dm(D, M)
A Naive O(NM) algorithm to do furthest points sampling, assuming the input is a N x N distance matrix
Parameters
- Dndarray (N, N)
An N x N distance matrix
- Minteger
Number of landmarks to compute
Return
- result: Dictionary
{‘perm’: An array of indices into X of the greedy permutation ‘lambdas’: Insertion radii of the landmarks ‘DLandmarks’: An MxN array of distances from landmarks to points in the point cloud}
- static get_greedy_perm_pc(X, M, csm_fn=<function GeometryUtils.get_csm>)
A Naive O(NM) algorithm to do furthest points sampling, assuming the input is a point cloud specified in Euclidean space. This saves computation over having compute the full distance matrix if the number of landmarks M << N
Parameters
- Xndarray (N, d)
An Nxd Euclidean point cloud
- Minteger
Number of landmarks to compute
- csm_fn: function X, Y -> D
Cross-similarity function (Euclidean by default)
Return
- result: Dictionary
{‘Y’: An Mxd array of landmarks, ‘perm’: An array of indices into X of the greedy permutation ‘lambdas’: Insertion radii of the landmarks ‘D’: An MxN array of distances from landmarks to points in X}
- static landmark_geodesic_distance(X, n_landmarks, n_neighbors)
Get a minmax sample of Euclidean point cloud using an approximation of the geodesic distance. Return the geodesic distance from the sample (landmarks) to the rest of the points.
Parameters
- Xndarray(n,d)
Point cloud of n points in d-dimensional Euclidean space.
- n_landmarksint
How many landmarks to sample.
- n_neighborsint
How many neighbors to use to compute the nearest neighbors graph used to approximate the geodesic distance.
Return
- dist_landmarks_pointsndarray(n_landmarks,n)
Geodesic distance from landmarks to the rest of the points. The first n_landmarks columns represent the landmarks, and are in the same order as landmarks. The rest of the columns represent the rest of the data points. The permutation that reorders the rows of X to correspond to the columns of dist_landmarks_points is given by perm_all_points.
- perm_all_pointsndarray(n)
Permutation of [0, …, n-1] which makes the rows of X correspond to the columns of dist_landmarks_points.