_images/logo.png

Utility Functions#

Utility functions for SRSF Manipulations

moduleauthor:: J. Derek Tucker <jdtuck@sandia.gov>

utility_functions.SqrtMean(gam, parallel=False, cores=-1)[source]#

calculates the srsf of warping functions with corresponding shooting vectors

Parameters:
  • gam – numpy ndarray of shape (M,N) of M warping functions with N samples

  • parallel – run in parallel (default = F)

  • cores – number of cores for parallel (default = -1 (all))

Return type:

2 numpy ndarray and vector

Return mu:

Karcher mean psi function

Return gam_mu:

vector of dim N which is the Karcher mean warping function

Return psi:

numpy ndarray of shape (M,N) of M SRSF of the warping functions

Return vec:

numpy ndarray of shape (M,N) of M shooting vectors

utility_functions.SqrtMeanInverse(gam)[source]#

finds the inverse of the mean of the set of the diffeomorphisms gamma

Parameters:

gam – numpy ndarray of shape (M,N) of N warping functions with M samples

Return type:

vector

Return gamI:

inverse of gam

utility_functions.SqrtMedian(gam)[source]#

calculates the median srsf of warping functions with corresponding shooting vectors

Parameters:

gam – numpy ndarray of shape (M,N) of M warping functions with N samples

Return type:

2 numpy ndarray and vector

Return gam_median:

Karcher median warping function

Return psi_meidan:

vector of dim N which is the Karcher median srsf function

Return psi:

numpy ndarray of shape (M,N) of M SRSF of the warping functions

Return vec:

numpy ndarray of shape (M,N) of M shooting vectors

utility_functions.cumtrapzmid(x, y, c, mid)[source]#

cumulative trapezoidal numerical integration taken from midpoint

Parameters:
  • x – vector of size N describing the time samples

  • y – vector of size N describing the function

  • c – midpointtic

  • mid – midpiont location

Return type:

vector

Return fa:

cumulative integration

utility_functions.diffop(n, binsize=1)[source]#

Creates a second order differential operator

Parameters:
  • n – dimension

  • binsize – dx (default = 1)

Return type:

numpy ndarray

Return m:

matrix describing differential operator

utility_functions.elastic_depth(f, time, method='DP2', lam=0.0, parallel=True)[source]#

calculates the elastic depth between functions in matrix f

Parameters:
  • f – matrix of size MxN (M time points for N functions)

  • time – vector of size M describing the sample points

  • method – method to apply optimization (default=”DP2”) options are “DP”,”DP2”,”RBFGS”,”cRBFGS”

  • lam – controls the elasticity (default = 0.0)

Return type:

scalar

Return amp:

amplitude depth

Return phase:

phase depth

utility_functions.elastic_distance(f1, f2, time, method='DP2', lam=0.0, alpha=None, return_dt_only=True)[source]#

” calculates the distances between function, where f1 is aligned to f2. In other words calculates the elastic distances

Parameters:
  • f1 – vector of size N

  • f2 – vector of size N

  • time – vector of size N describing the sample points

  • method – method to apply optimization (default=”DP2”) options are “DP”,”DP2”,”RBFGS”,”cRBFGS”

  • lam – controls the elasticity (default = 0.0)

  • alpha – makes alpha * dx + (1-alpha) * dy

  • return_dt_only – returns only dt if alpha is set

Return type:

scalar

Return Dy:

amplitude distance

Return Dx:

phase distance

Return Dt:

combined distance

utility_functions.f_K_fold(Nobs, K=5)[source]#

generates sample indices for K-fold cross validation

:param Nobs number of observations :param K number of folds

Return type:

numpy ndarray

Return train:

train indexes (Nobs*(K-1)/K X K)

Return test:

test indexes (Nobs*(1/K) X K)

utility_functions.f_to_srsf(f, time, smooth=False)[source]#

converts f to a square-root slope function (SRSF)

Parameters:
  • f – vector of size N samples

  • time – vector of size N describing the sample points

Return type:

vector

Return q:

srsf of f

utility_functions.geigen(Amat, Bmat, Cmat)[source]#

generalized eigenvalue problem of the form

max tr L’AM / sqrt(tr L’BL tr M’CM) w.r.t. L and M

:param Amat numpy ndarray of shape (M,N) :param Bmat numpy ndarray of shape (M,N) :param Bmat numpy ndarray of shape (M,N)

Return type:

numpy ndarray

Return values:

eigenvalues

Return Lmat:

left eigenvectors

Return Mmat:

right eigenvectors

utility_functions.gradient_spline(time, f, smooth=False)[source]#

This function takes the gradient of f using b-spline smoothing

Parameters:
  • time – vector of size N describing the sample points

  • f – numpy ndarray of shape (M,N) of M functions with N samples

  • smooth – smooth data (default = F)

Return type:

tuple of numpy ndarray

Return f0:

smoothed functions functions

Return g:

first derivative of each function

Return g2:

second derivative of each function

utility_functions.innerprod_q(time, q1, q2)[source]#

calculates the innerproduct between two srsfs

:param time vector descrbing time samples :param q1 vector of srsf 1 :param q2 vector of srsf 2

Return type:

scalar

Return val:

inner product value

utility_functions.invertGamma(gam)[source]#

finds the inverse of the diffeomorphism gamma

Parameters:

gam – vector describing the warping function

Return type:

vector

Return gamI:

inverse of gam

utility_functions.optimum_reparam(q1, time, q2, method='DP2', lam=0.0, penalty='roughness', grid_dim=7)[source]#

calculates the warping to align srsf q2 to q1

Parameters:
  • q1 – vector of size N or array of NxM samples of first SRSF

  • time – vector of size N describing the sample points

  • q2 – vector of size N or array of NxM samples samples of second SRSF

  • method – method to apply optimization (default=”DP2”) options are “DP”,”DP2”,”RBFGS”,”cRBFGS”

  • lam – controls the amount of elasticity (default = 0.0)

  • penalty – penalty type (default=”roughness”) options are “roughness”, “l2gam”, “l2psi”, “geodesic”. Only roughness implemented in all methods. To use others method needs to be “RBFGS” or “cRBFGS”

  • grid_dim – size of the grid, for the DP2 method only (default = 7)

Return type:

vector

Return gam:

describing the warping function used to align q2 with q1

utility_functions.optimum_reparam_pair(q, time, q1, q2, lam=0.0)[source]#

calculates the warping to align srsf pair q1 and q2 to q

Parameters:
  • q – vector of size N or array of NxM samples of first SRSF

  • time – vector of size N describing the sample points

  • q1 – vector of size N or array of NxM samples samples of second SRSF

  • q2 – vector of size N or array of NxM samples samples of second SRSF

  • lam – controls the amount of elasticity (default = 0.0)

Return type:

vector

Return gam:

describing the warping function used to align q2 with q1

utility_functions.outlier_detection(q, time, mq, k=1.5)[source]#

calculates outlier’s using geodesic distances of the SRSFs from the median

Parameters:
  • q – numpy ndarray of N x M of M SRS functions with N samples

  • time – vector of size N describing the sample points

  • mq – median calculated using time_warping.srsf_align()

  • k – cutoff threshold (default = 1.5)

Returns:

q_outlier: outlier functions

utility_functions.randomGamma(gam, num)[source]#

generates random warping functions

Parameters:
  • gam – numpy ndarray of N x M of M of warping functions

  • num – number of random functions

Returns:

rgam: random warping functions

utility_functions.resamplefunction(x, n)[source]#

resample function using n points

Parameters:
  • x – functions

  • n – number of points

Return type:

numpy array

Return xn:

resampled function

utility_functions.rgam(N, sigma, num, mu_gam=None)[source]#

Generates random warping functions

Parameters:
  • N – length of warping function

  • sigma – variance of warping functions

  • num – number of warping functions

:param mu_gam mean warping function (default identity) :return: gam: numpy ndarray of warping functions

utility_functions.smooth_data(f, sparam=1)[source]#

This function smooths a collection of functions using a box filter

Parameters:
  • f – numpy ndarray of shape (M,N) of M functions with N samples

  • sparam – Number of times to run box filter (default = 25)

Return type:

numpy ndarray

Return f:

smoothed functions functions

utility_functions.srsf_to_f(q, time, f0=0.0)[source]#

converts q (srsf) to a function

Parameters:
  • q – vector of size N samples of srsf

  • time – vector of size N describing time sample points

  • f0 – initial value

Return type:

vector

Return f:

function

utility_functions.update_progress(progress)[source]#

This function creates a progress bar

Parameters:

progress – fraction of progress

utility_functions.warp_f_gamma(time, f, gam)[source]#

warps a function f by gam

:param time vector describing time samples :param q vector describing srsf :param gam vector describing warping function

Return type:

numpy ndarray

Return f_temp:

warped srsf

utility_functions.warp_q_gamma(time, q, gam)[source]#

warps a srsf q by gam

:param time vector describing time samples :param q vector describing srsf :param gam vector describing warping function

Return type:

numpy ndarray

Return q_temp:

warped srsf

utility_functions.zero_crossing(Y, q, bt, time, y_max, y_min, gmax, gmin)[source]#

finds zero-crossing of optimal gamma, gam = s*gmax + (1-s)*gmin from elastic regression model

Parameters:
  • Y – response

  • q – predicitve function

  • bt – basis function

  • time – time samples

  • y_max – maximum repsonse for warping function gmax

  • y_min – minimum response for warping function gmin

  • gmax – max warping function

  • gmin – min warping fucntion

Return type:

numpy array

Return gamma:

optimal warping function