mass.util.matrix
Containins basic matrix operations that can be applied to a mass
model.
To assist with various matrix operations using different packages, the
following values can be provided to the array_type
argument to set the
return type of the output. Valid matrix types include:
'dense'
for anumpy.ndarray
'dok'
for ascipy.sparse.dok_matrix
'lil'
for ascipy.sparse.lil_matrix
'DataFrame'
for apandas.DataFrame
'symbolic'
for asympy.MutableDenseMatrix
For all matrix types, species (excluding genes) are the row indicies and reactions are the column indicies.
There are also several methods that are nearly identical to scipy.linalg
methods, with the main exception being that matrix conversions are performed
beforehand to ensure that valid input is passed to the scipy
method.
These methods include:
Module Contents
Functions
|
Create the gradient matrix for a given model. |
|
Create the kappa matrix for a given model. |
|
Create the gamma matrix for a given model. |
|
Get the jacobian matrix for a given model. |
|
Compute an approximate basis for the nullspace of a matrix. |
|
Compute an approximate basis for the left nullspace of a matrix. |
|
Compute an approximate basis for the columnspace of a matrix. |
|
Compute an approximate basis for the rowspace of a matrix. |
|
Estimate the rank (i.e. the dimension of the nullspace) of a matrix. |
|
Get the singular value decomposition of a matrix. |
|
Get the eigenvalues of a matrix. |
|
Convert a matrix to a different type. |
- mass.util.matrix.gradient(model, use_parameter_values=True, use_concentration_values=True, array_type='dense')[source]
Create the gradient matrix for a given model.
- Parameters
model (MassModel) – The
MassModel
to construct the matrix for.use_parameter_values (bool) – Whether to substitute the numerical values for parameters into the matrix. If
True
then numerical values of the kinetic parameters are substituted into the matrix. Otherwise parameters in the matrix are left as symbols. Default isTrue
.use_concentration_values (bool) – Whether to substitute the numerical values for concentrations into the matrix. If
True
then numerical values of the initial conditions are substituted into the matrix. Otherwise species concentrations in the matrix are left as symbols. Default isTrue
.array_type (str) – A string identifiying the desired format for the returned matrix. Default is
'dense'
. See thematrix
module documentation for more information on thearray_type
- Returns
The gradient matrix for the model.
- Return type
matrix of type
array_type
- mass.util.matrix.kappa(model, use_parameter_values=True, use_concentration_values=True, array_type='dense')[source]
Create the kappa matrix for a given model.
Notes
The kappa matrix is the diagnolization of the norms for the rows in the gradient matrix.
- Parameters
model (MassModel) – The
MassModel
to construct the matrix for.use_parameter_values (bool) – Whether to substitute the numerical values for parameters into the matrix. If
True
then numerical values of the kinetic parameters are substituted into the matrix. Otherwise parameters in the matrix are left as symbols. Default isTrue
.use_concentration_values (bool) – Whether to substitute the numerical values for concentrations into the matrix. If
True
then numerical values of the initial conditions are substituted into the matrix. Otherwise species concentrations in the matrix are left as symbols. Default isTrue
.array_type (str) – A string identifiying the desired format for the returned matrix. Default is
'dense'
. See thematrix
module documentation for more information on thearray_type
.
- Returns
The kappa matrix for the model.
- Return type
matrix of type
array_type
- mass.util.matrix.gamma(model, use_parameter_values=True, use_concentration_values=True, array_type='dense')[source]
Create the gamma matrix for a given model.
Notes
The gamma matrix is composed of the 1-norms of the gradient matrix.
- Parameters
model (MassModel) – The
MassModel
to construct the matrix for.use_parameter_values (bool) – Whether to substitute the numerical values for parameters into the matrix. If
True
then numerical values of the kinetic parameters are substituted into the matrix. Otherwise parameters in the matrix are left as symbols. Default isTrue
.use_concentration_values (bool) – Whether to substitute the numerical values for concentrations into the matrix. If
True
then numerical values of the initial conditions are substituted into the matrix. Otherwise species concentrations in the matrix are left as symbols. Default isTrue
.array_type (str) – A string identifiying the desired format for the returned matrix. Default is
'dense'
. See thematrix
module documentation for more information on thearray_type
.
- Returns
The gamma matrix for the model.
- Return type
matrix of type
array_type
- mass.util.matrix.jacobian(model, jacobian_type='species', use_parameter_values=True, use_concentration_values=True, array_type='dense')[source]
Get the jacobian matrix for a given model.
- Parameters
model (MassModel) – The
MassModel
to construct the matrix for.jacobian_type (str) – Either the string
'species'
to obtain the jacobian matrix with respect to species, or the string'reactions'
to obtain the jacobian matrix with respect to the reactions. Default is'reactions'
.use_parameter_values (bool) – Whether to substitute the numerical values for parameters into the matrix. If
True
then numerical values of the kinetic parameters are substituted into the matrix. Otherwise parameters in the matrix are left as symbols. Default isTrue
.use_concentration_values (bool) – Whether to substitute the numerical values for concentrations into the matrix. If
True
then numerical values of the initial conditions are substituted into the matrix. Otherwise species concentrations in the matrix are left as symbols. Default isTrue
.array_type (str) – A string identifiying the desired format for the returned matrix. Default is
'dense'
. See thematrix
module documentation for more information on thearray_type
.
- Returns
The jacobian matrix for the model.
- Return type
matrix of type
array_type
- mass.util.matrix.nullspace(matrix, atol=1e-13, rtol=0, decimal_precision=False)[source]
Compute an approximate basis for the nullspace of a matrix.
The algorithm used by this function is based on singular value decomposition.
Notes
If both
atol
andrtol
are positive, the combined tolerance is the maximum of the two; that is:tol = max(atol, rtol * smax)
Singular values smaller than
tol
are considered to be zero.Similar to the
cobra.util.array.nullspace
function, but includes utlization of thedecimal_precision
in theMassConfiguration
and sets values below the tolerance to 0.Taken from the numpy cookbook and extended.
- Parameters
matrix (array-like) – The matrix to decompose. The matrix should be at most 2-D. A 1-D array with length
k
will be treated as a 2-D with shape(1, k)
.atol (float) – The absolute tolerance for a zero singular value. Singular values smaller than
atol
are considered to be zero.rtol (float) – The relative tolerance. Singular values less than
rtol * smax
are considered to be zero, wheresmax
is the largest singular value.decimal_precision (bool) – Whether to apply the
decimal_precision
set in theMassConfiguration
to the nullspace values before comparing to the tolerance. Default isFalse
.
- Returns
ns – If
matrix
is an array with shape(m, k)
, thenns
will be an array with shape(k, n)
, wheren
is the estimated dimension of the nullspace ofmatrix
. The columns ofns
are a basis for the nullspace; each element in the dot product of the matrix and the nullspace will be approximately 0.- Return type
- mass.util.matrix.left_nullspace(matrix, atol=1e-13, rtol=0, decimal_precision=False)[source]
Compute an approximate basis for the left nullspace of a matrix.
The algorithm used by this function is based on singular value decomposition.
Notes
If both
atol
andrtol
are positive, the combined tolerance is the maximum of the two; that is:tol = max(atol, rtol * smax)
Singular values smaller than
tol
are considered to be zero.- Parameters
matrix (array-like) – The matrix to decompose. The matrix should be at most 2-D. A 1-D array with length
k
will be treated as a 2-D with shape(1, k)
.atol (float) – The absolute tolerance for a zero singular value. Singular values smaller than
atol
are considered to be zero.rtol (float) – The relative tolerance. Singular values less than
rtol * smax
are considered to be zero, wheresmax
is the largest singular value.decimal_precision (bool) – Whether to apply the
decimal_precision
set in theMassConfiguration
to the left nullspace values before comparing to the tolerance. Default isFalse
.
- Returns
lns – If
matrix
is an array with shape(m, k)
, thenlns
will be an array with shape(n, m)
, wheren
is the estimated dimension of the left nullspace ofmatrix
. The rows oflns
are a basis for the left nullspace; each element in the dot product of the matrix and the left nullspace will be approximately 0.- Return type
See also
nullspace()
Base function.
- mass.util.matrix.columnspace(matrix, atol=1e-13, rtol=0, decimal_precision=False)[source]
Compute an approximate basis for the columnspace of a matrix.
This function utilizes the
scipy.linalg.qr()
function to obtain an orthogonal basis for the columnspace of the matrix.Notes
If both
atol
andrtol
are positive, the combined tolerance is the maximum of the two; that is:tol = max(atol, rtol * smax)
Singular values smaller than
tol
are considered to be zero.- Parameters
matrix (array-like) – The matrix to decompose. The matrix should be at most 2-D. A 1-D array with length
k
will be treated as a 2-D with shape(1, k)
.atol (float) – The absolute tolerance for a zero singular value. Singular values smaller than
atol
are considered to be zero.rtol (float) – The relative tolerance. Singular values less than
rtol * smax
are considered to be zero, wheresmax
is the largest singular value.decimal_precision (bool) – Whether to apply the
decimal_precision
set in theMassConfiguration
to the columnspace values before comparing to the tolerance. Default isFalse
.
- Returns
cs – If
matrix
is an array with shape(m, k)
, thencs
will be an array with shape(m, n)
, wheren
is the estimated dimension of the columnspace ofmatrix
. The columns ofcs
are a basis for the columnspace.- Return type
- mass.util.matrix.rowspace(matrix, atol=1e-13, rtol=0, decimal_precision=False)[source]
Compute an approximate basis for the rowspace of a matrix.
This function utilizes the
scipy.linalg.qr()
function to obtain an orthogonal basis for the rowspace of the matrix.Notes
If both
atol
andrtol
are positive, the combined tolerance is the maximum of the two; that is:tol = max(atol, rtol * smax)
Singular values smaller than
tol
are considered to be zero.- Parameters
matrix (array-like) – The matrix to decompose. The matrix should be at most 2-D. A 1-D array with length
k
will be treated as a 2-D with shape(1, k)
.atol (float) – The absolute tolerance for a zero singular value. Singular values smaller than
atol
are considered to be zero.rtol (float) – The relative tolerance. Singular values less than
rtol * smax
are considered to be zero, wheresmax
is the largest singular value.decimal_precision (bool) – Whether to apply the
decimal_precision
set in theMassConfiguration
to the rowspace values before comparing to the tolerance. Default isFalse
.
- Returns
rs – If
matrix
is an array with shape(m, k)
, thenrs
will be an array with shape(n, k)
, wheren
is the estimated dimension of the rowspace ofmatrix
. The columns ofrs
are a basis for the rowspace.- Return type
See also
columnspace()
Base function.
- mass.util.matrix.matrix_rank(matrix, atol=1e-13, rtol=0)[source]
Estimate the rank (i.e. the dimension of the nullspace) of a matrix.
The algorithm used by this function is based on singular value decomposition. Taken from the
scipy
cookbook.Notes
If both
atol
andrtol
are positive, the combined tolerance is the maximum of the two; that is:tol = max(atol, rtol * smax)
Singular values smaller than
tol
are considered to be zero.- Parameters
matrix (array-like) – The matrix to obtain the rank for. The matrix should be at most 2-D. A 1-D array with length
k
will be treated as a 2-D with shape(1, k)
.atol (float) – The absolute tolerance for a zero singular value. Singular values smaller than
atol
are considered to be zero.rtol (float) – The relative tolerance. Singular values less than
rtol * smax
are considered to be zero, wheresmax
is the largest singular value.
- Returns
rank – The estimated rank of the matrix.
- Return type
See also
numpy.linalg.matrix_rank()
mass.util.matrix.matrix_rank()
is nearly identical to this function, but it does not provide the option of the absolute tolerance.
- mass.util.matrix.svd(matrix, **kwargs)[source]
Get the singular value decomposition of a matrix.
kwargs
are passed on toscipy.linalg.svd()
.- Parameters
matrix (array-like) – The matrix to decompose. The matrix should be at most 2-D. A 1-D array with length
k
will be treated as a 2-D with shape(1, k)
.- Returns
U (ndarray) – Unitary matrix having left singular vectors as columns. Of shape
(M, M)
or(M, K)
, depending on full_matrices.s (ndarray) – The singular values, sorted in non-increasing order. Of shape
(K, )
, withK = min(M, N)
.Vh (ndarray) – Unitary matrix having right singular vectors as rows. Of shape
(N, N)
or(K, N)
depending on full_matrices.For
compute_uv=False
, onlys
is returned.
See also
scipy.linalg.svd()
Base function.
- mass.util.matrix.eig(matrix, left=False, right=False, **kwargs)[source]
Get the eigenvalues of a matrix.
kwargs
are passed on toscipy.linalg.eig()
- Parameters
matrix (array-like) – The matrix to decompose. The matrix should be at most 2-D. A 1-D array with length
k
will be treated as a 2-D with shape(1, k)
.left (bool) – Whether to calculate and return left eigenvectors. Default is
False
.right (bool) – Whether to calculate and return right eigenvectors. Default is
True
.
- Returns
w ((M, ) double or complex ndarray) – The eigenvalues, each repeated according to its multiplicity.
vl ((M, M) double or complex ndarray) – The normalized left eigenvector corresponding to the eigenvalue
w[i]
is the column v[:,i]. Only returned ifleft=True
.vr ((M, M) double or complex ndarray) – The normalized right eigenvector corresponding to the eigenvalue
w[i]
is the columnvr[:,i]
. Only returned ifright=True
.
See also
scipy.linalg.eig()
Base function.
- mass.util.matrix.convert_matrix(matrix, array_type, dtype, row_ids=None, col_ids=None)[source]
Convert a matrix to a different type.
- Parameters
matrix (array-like) – The matrix to convert.
array_type (str) – A string identifiying the desired format for the returned matrix. Valid matrix types include
'dense'
,'dok'
,'lil'
,'DataFrame'
, and'symbolic'
See thematrix
module documentation for more information on thearray_type
.dtype (data-type) – The desired array data-type for the matrix.
row_ids (array-like) – The idenfifiers for each row. Only used if type is
'DataFrame'
.col_ids (array-like) – The idenfifiers for each column. Only used if type is
'DataFrame'
.
Warning
This method is NOT the safest way to convert a matrix to another type. To safely convert a matrix into another type, use the
'array_type'
argument in the method that returns the desired matrix.