Package csb :: Package bio :: Package utils
[frames] | no frames]

Package utils

source code

Computational utility functions.

This module defines a number of low-level, numerical, high-performance utility functions like rmsd for example.

Functions
(n,3) numpy.array
average_structure(X)
Calculate an average structure from an ensemble of structures (i.e.
source code
(d,) numpy.array
center_of_mass(x, m=None)
Compute the mean of a set of (optionally weighted) points.
source code
(m,) numpy.array
distance(X, Y)
Distance between X and Y along the last axis.
source code
numpy array
distance_matrix(X, Y=None)
Calculates a matrix of pairwise distances
source code
(m,) numpy.array
distance_sq(X, Y)
Squared distance between X and Y along the last axis.
source code
tuple
fit(X, Y)
Return the translation vector and the rotation matrix minimizing the RMSD between two sets of vectors, i.e.
source code
 
fit_wellordered(X, Y, n_iter=None, n_stdv=2, tol_rmsd=0.5, tol_stdv=0.05, full_output=False)
Matche two arrays onto each other by iteratively throwing out highly deviating entries.
source code
(d,d) numpy.array
inertia_tensor(x, m=None)
Compute the inertia tensor of a set of (optionally weighted) points.
source code
bool
is_mirror_image(X, Y)
Check if two configurations X and Y are mirror images (i.e.
source code
 
probabilistic_fit(X, Y, w=None, niter=10)
Generate a superposition of X,Y where:
source code
(d,) numpy.array
radius_of_gyration(x, m=None)
Compute the radius of gyration of a set of (optionally weighted) points.
source code
float
rmsd(X, Y)
Calculate the root mean squared deviation (RMSD) using Kabsch' formula.
source code
float
rmsd_cur(X, Y)
Calculate the RMSD of two conformations as they are (no fitting is done).
source code
(d,d) numpy.array
second_moments(x, m=None)
Compute the tensor second moments of a set of (optionally weighted) points.
source code
float
tm_score(x, y, L=None, d0=None)
Evaluate the TM-score of two conformations as they are (no fitting is done).
source code
tuple
tm_superimpose(x, y, fit_method=<function fit at 0x3ff3a28>, L=None, d0=None, L_ini_min=4, iL_step=1)
Compute the TM-score of two protein coordinate vector sets.
source code
float
torsion_rmsd(x, y)
Compute the circular RMSD of two phi/psi angle sets.
source code
tuple
wfit(X, Y, w)
Return the translation vector and the rotation matrix minimizing the weighted RMSD between two sets of vectors, i.e.
source code
float
wrmsd(X, Y, w)
Calculate the weighted root mean squared deviation (wRMSD) using Kabsch' formula.
source code
Variables
  __package__ = 'csb.bio.utils'
Function Details

average_structure(X)

source code 

Calculate an average structure from an ensemble of structures (i.e. X is a rank-3 tensor: X[i] is a (N,3) configuration matrix).

Parameters:
  • X (numpy array) - m x n x 3 input vector
Returns: (n,3) numpy.array
average structure

center_of_mass(x, m=None)

source code 

Compute the mean of a set of (optionally weighted) points.

Parameters:
  • x (numpy.array) - array of rank (n,d) where n is the number of points and d the dimension
  • m (numpy.array) - rank (n,) array of masses / weights
Returns: (d,) numpy.array
center of mass

distance(X, Y)

source code 

Distance between X and Y along the last axis.

Parameters:
  • X (numpy array) - m x n input vector
  • Y (numpy array) - m x n input vector
Returns: (m,) numpy.array
scalar or array of length m

distance_matrix(X, Y=None)

source code 

Calculates a matrix of pairwise distances

Parameters:
  • X (numpy array) - m x n input vector
  • Y (numpy array) - k x n input vector or None, which defaults to Y=X
Returns: numpy array
m x k distance matrix

distance_sq(X, Y)

source code 

Squared distance between X and Y along the last axis. For details, see distance.

Returns: (m,) numpy.array
scalar or array of length m

fit(X, Y)

source code 

Return the translation vector and the rotation matrix minimizing the RMSD between two sets of vectors, i.e. if

>>> R,t = fit(X,Y)

then

>>> Y = dot(Y, transpose(R)) + t

will be the fitted configuration.

Parameters:
  • X (numpy array) - 3 x n input vector
  • Y (numpy array) - 3 x n input vector
Returns: tuple
3 x 3 rotation matrix and 3 x 1 translation vector

fit_wellordered(X, Y, n_iter=None, n_stdv=2, tol_rmsd=0.5, tol_stdv=0.05, full_output=False)

source code 

Matche two arrays onto each other by iteratively throwing out highly deviating entries.

(Reference: Nilges et al.: Delineating well-ordered regions in protein structure ensembles).

Parameters:
  • X (numpy array) - 3 x n input vector
  • Y (numpy array) - 3 x n input vector
  • n_stdv - number of standard deviations above which points are considered to be outliers
  • tol_rmsd - tolerance in rmsd
  • tol_stdv - tolerance in standard deviations
  • full_output - also return full history of values calculated by the algorithm

inertia_tensor(x, m=None)

source code 

Compute the inertia tensor of a set of (optionally weighted) points.

Parameters:
  • x (numpy.array) - array of rank (n,d) where n is the number of points and d the dimension
  • m (numpy.array) - rank (n,) array of masses / weights
Returns: (d,d) numpy.array
inertia tensor

is_mirror_image(X, Y)

source code 

Check if two configurations X and Y are mirror images (i.e. their optimal superposition involves a reflection).

Parameters:
  • X (numpy array) - n x 3 input vector
  • Y (numpy array) - n x 3 input vector
Returns: bool

probabilistic_fit(X, Y, w=None, niter=10)

source code 

Generate a superposition of X,Y where:

   R ~ exp(trace(dot(transpose(dot(transpose(X-t), Y)), R)))
   t ~ N(t_opt, 1 / sqrt(N))

radius_of_gyration(x, m=None)

source code 

Compute the radius of gyration of a set of (optionally weighted) points.

Parameters:
  • x (numpy.array) - array of rank (n,d) where n is the number of points and d the dimension
  • m (numpy.array) - rank (n,) array of masses / weights
Returns: (d,) numpy.array
center of mass

rmsd(X, Y)

source code 

Calculate the root mean squared deviation (RMSD) using Kabsch' formula.

Parameters:
  • X (numpy array) - 3 x n input vector
  • Y (numpy array) - 3 x n input vector
Returns: float
rmsd value between the input vectors

rmsd_cur(X, Y)

source code 

Calculate the RMSD of two conformations as they are (no fitting is done). For details, see rmsd.

Returns: float
rmsd value between the input vectors

second_moments(x, m=None)

source code 

Compute the tensor second moments of a set of (optionally weighted) points.

Parameters:
  • x (numpy.array) - array of rank (n,d) where n is the number of points and d the dimension
  • m (numpy.array) - rank (n,) array of masses / weights
Returns: (d,d) numpy.array
second moments

tm_score(x, y, L=None, d0=None)

source code 

Evaluate the TM-score of two conformations as they are (no fitting is done).

Parameters:
  • x (numpy array) - 3 x n input array
  • y (numpy array) - 3 x n input array
  • L (int) - length for normalization (default: len(x))
  • d0 (float) - d0 in Angstroms (default: calculate from L)
Returns: float
computed TM-score

tm_superimpose(x, y, fit_method=<function fit at 0x3ff3a28>, L=None, d0=None, L_ini_min=4, iL_step=1)

source code 

Compute the TM-score of two protein coordinate vector sets. Reference: http://zhanglab.ccmb.med.umich.edu/TM-score

Parameters:
  • x (numpy.array) - 3 x n input vector
  • y (numpy.array) - 3 x n input vector
  • fit_method (function) - a reference to a proper fitting function, e.g. fit or fit_wellordered
  • L (int) - length for normalization (default: len(x))
  • d0 (float) - d0 in Angstroms (default: calculate from L)
  • L_ini_min (int) - minimum length of initial alignment window (increase to speed up but loose precision, a value of 0 disables local alignment initialization)
  • iL_step (int) - initial alignment window shift (increase to speed up but loose precision)
Returns: tuple
rotation matrix, translation vector, TM-score

torsion_rmsd(x, y)

source code 

Compute the circular RMSD of two phi/psi angle sets.

Parameters:
  • x (array) - query phi/psi angles (Nx2 array, in radians)
  • y (array) - subject phi/psi angles (Nx2 array, in radians)
Returns: float

wfit(X, Y, w)

source code 

Return the translation vector and the rotation matrix minimizing the weighted RMSD between two sets of vectors, i.e. if

>>> R,t = wfit(X,Y,w)

then

>>> Y = dot(Y, transpose(R)) + t

will be the fitted configuration.

Parameters:
  • X (numpy array) - 3 x n input vector
  • Y (numpy array) - 3 x n input vector
  • w (numpy array) - input weights
Returns: tuple
3 x 3 rotation matrix and 3 x 1 translation vector

wrmsd(X, Y, w)

source code 

Calculate the weighted root mean squared deviation (wRMSD) using Kabsch' formula.

Parameters:
  • X (numpy array) - 3 x n input vector
  • Y (numpy array) - 3 x n input vector
  • w (numpy array) - input weights
Returns: float
rmsd value between the input vectors