[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]
Linear algebra functions. More...
Classes | |
class | LeastAngleRegressionOptions |
Pass options to leastAngleRegression(). More... | |
class | Matrix |
Functions | |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | abs (MultiArrayView< 2, T, C > const &v) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | acos (MultiArrayView< 2, T, C > const &v) |
template<class T , class C1 , class C2 , class C3 > | |
void | add (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class T , class C > | |
int | argMax (MultiArrayView< 2, T, C > const &a) |
Find the index of the maximum element in a matrix. | |
template<class T , class C , class UnaryFunctor > | |
int | argMaxIf (MultiArrayView< 2, T, C > const &a, UnaryFunctor condition) |
Find the index of the maximum element in a matrix subject to a condition. | |
template<class T , class C > | |
int | argMin (MultiArrayView< 2, T, C > const &a) |
Find the index of the minimum element in a matrix. | |
template<class T , class C , class UnaryFunctor > | |
int | argMinIf (MultiArrayView< 2, T, C > const &a, UnaryFunctor condition) |
Find the index of the minimum element in a matrix subject to a condition. | |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | asin (MultiArrayView< 2, T, C > const &v) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | atan (MultiArrayView< 2, T, C > const &v) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | ceil (MultiArrayView< 2, T, C > const &v) |
template<class T , class C1 , class C2 > | |
bool | choleskyDecomposition (MultiArrayView< 2, T, C1 > const &A, MultiArrayView< 2, T, C2 > &L) |
template<class T , class C1 , class C2 , class C3 > | |
void | choleskySolve (MultiArrayView< 2, T, C1 > const &L, MultiArrayView< 2, T, C2 > const &b, MultiArrayView< 2, T, C3 > &x) |
template<class T , class C > | |
MultiArrayIndex | columnCount (const MultiArrayView< 2, T, C > &x) |
template<... > | |
void | columnStatistics (...) |
template<class T , class C > | |
MultiArrayView< 2, T, C > | columnVector (MultiArrayView< 2, T, C > const &m, MultiArrayIndex d) |
template<class T , class C > | |
MultiArrayView< 2, T, C > | columnVector (MultiArrayView< 2, T, C > const &m, MultiArrayShape< 2 >::type first, int end) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | cos (MultiArrayView< 2, T, C > const &v) |
template<class T , class C > | |
TemporaryMatrix< T > | covarianceMatrixOfColumns (MultiArrayView< 2, T, C > const &features) |
Compute the covariance matrix between the columns of a matrix features. | |
template<class T1 , class C1 , class T2 , class C2 > | |
void | covarianceMatrixOfColumns (MultiArrayView< 2, T1, C1 > const &features, MultiArrayView< 2, T2, C2 > &covariance) |
Compute the covariance matrix between the columns of a matrix features. | |
template<class T , class C > | |
TemporaryMatrix< T > | covarianceMatrixOfRows (MultiArrayView< 2, T, C > const &features) |
Compute the covariance matrix between the rows of a matrix features. | |
template<class T1 , class C1 , class T2 , class C2 > | |
void | covarianceMatrixOfRows (MultiArrayView< 2, T1, C1 > const &features, MultiArrayView< 2, T2, C2 > &covariance) |
Compute the covariance matrix between the rows of a matrix features. | |
template<class T , class C1 , class C2 , class C3 > | |
void | cross (const MultiArrayView< 1, T, C1 > &x, const MultiArrayView< 1, T, C2 > &y, MultiArrayView< 1, T, C3 > &r) |
template<class T , class C1 , class C2 > | |
TemporaryMatrix< T > | cross (const MultiArrayView< 2, T, C1 > &x, const MultiArrayView< 2, T, C2 > &y) |
template<class T , class C1 , class C2 , class C3 > | |
void | cross (const MultiArrayView< 2, T, C1 > &x, const MultiArrayView< 2, T, C2 > &y, MultiArrayView< 2, T, C3 > &r) |
template<class T , class C1 > | |
NumericTraits< T >::Promote | determinant (MultiArrayView< 2, T, C1 > const &a, std::string method="default") |
template<class T , class C > | |
TemporaryMatrix< T > | diagonalMatrix (MultiArrayView< 2, T, C > const &v) |
template<class T , class C1 , class C2 > | |
void | diagonalMatrix (MultiArrayView< 2, T, C1 > const &v, MultiArrayView< 2, T, C2 > &r) |
template<class T , class C1 , class C2 > | |
NormTraits< T >::SquaredNormType | dot (const MultiArrayView< 1, T, C1 > &x, const MultiArrayView< 1, T, C2 > &y) |
template<class T , class C1 , class C2 > | |
NormTraits< T >::SquaredNormType | dot (const MultiArrayView< 2, T, C1 > &x, const MultiArrayView< 2, T, C2 > &y) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | exp (MultiArrayView< 2, T, C > const &v) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | floor (MultiArrayView< 2, T, C > const &v) |
template<class T > | |
TemporaryMatrix< T > | identityMatrix (MultiArrayIndex size) |
template<class T , class C > | |
void | identityMatrix (MultiArrayView< 2, T, C > &r) |
template<class T , class C > | |
TemporaryMatrix< T > | inverse (const MultiArrayView< 2, T, C > &v) |
template<class T , class C1 , class C2 > | |
bool | inverse (const MultiArrayView< 2, T, C1 > &v, MultiArrayView< 2, T, C2 > &res) |
template<class T , class C > | |
bool | isSymmetric (const MultiArrayView< 2, T, C > &v) |
template<class T , class C1 , class C2 > | |
TemporaryMatrix< T > | joinHorizontally (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T , class C1 , class C2 > | |
TemporaryMatrix< T > | joinVertically (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<... > | |
unsigned int | leastAngleRegression (...) |
template<class T , class C1 , class C2 , class C3 > | |
bool | leastSquares (MultiArrayView< 2, T, C1 > const &A, MultiArrayView< 2, T, C2 > const &b, MultiArrayView< 2, T, C3 > &x, std::string method="QR") |
template<... > | |
bool | linearSolve (...) |
template<class T , class C1 , class C2 , class C3 > | |
bool | linearSolveLowerTriangular (const MultiArrayView< 2, T, C1 > &l, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > x) |
template<class T , class C1 , class C2 , class C3 > | |
bool | linearSolveUpperTriangular (const MultiArrayView< 2, T, C1 > &r, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > x) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | log (MultiArrayView< 2, T, C > const &v) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | log10 (MultiArrayView< 2, T, C > const &v) |
template<class T , class C1 > | |
T | logDeterminant (MultiArrayView< 2, T, C1 > const &a) |
template<class T , class C1 , class C2 > | |
TemporaryMatrix< T > | mmul (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T , class C1 , class C2 , class C3 > | |
void | mmul (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<... > | |
unsigned int | nonnegativeLeastSquares (...) |
template<class T , class C1 , class C2 , class C3 > | |
bool | nonsymmetricEigensystem (MultiArrayView< 2, T, C1 > const &a, MultiArrayView< 2, std::complex< T >, C2 > &ew, MultiArrayView< 2, T, C3 > &ev) |
template<class T , class ALLOC > | |
Matrix< T, ALLLOC >::NormType | norm (const Matrix< T, ALLLOC > &a) |
template<class T > | |
TemporaryMatrix< T > | ones (MultiArrayIndex rows, MultiArrayIndex cols) |
template<class T , class A , int N, class DATA , class DERIVED > | |
TinyVector< T, N > | operator* (const Matrix< T, A > &a, const TinyVectorBase< T, N, DATA, DERIVED > &b) |
template<class T , class C , class U > | |
TemporaryMatrix< T > | operator* (const MultiArrayView< 2, T, C > &a, PointWise< U > b) |
template<class T , class C > | |
TemporaryMatrix< T > | operator* (const MultiArrayView< 2, T, C > &a, T b) |
template<class T , class C1 , class C2 > | |
TemporaryMatrix< T > | operator* (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T , int N, class DATA , class DERIVED , class A > | |
TinyVector< T, N > | operator* (const TinyVectorBase< T, N, DATA, DERIVED > &a, const Matrix< T, A > &b) |
template<class T , class C > | |
TemporaryMatrix< T > | operator* (T a, const MultiArrayView< 2, T, C > &b) |
template<class T , class C > | |
TemporaryMatrix< T > | operator+ (const MultiArrayView< 2, T, C > &a, T b) |
template<class T , class C1 , class C2 > | |
TemporaryMatrix< T > | operator+ (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T , class C > | |
TemporaryMatrix< T > | operator+ (T a, const MultiArrayView< 2, T, C > &b) |
template<class T , class C > | |
TemporaryMatrix< T > | operator- (const MultiArrayView< 2, T, C > &a) |
template<class T , class C > | |
TemporaryMatrix< T > | operator- (const MultiArrayView< 2, T, C > &a, T b) |
template<class T , class C1 , class C2 > | |
TemporaryMatrix< T > | operator- (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T , class C > | |
TemporaryMatrix< T > | operator- (T a, const MultiArrayView< 2, T, C > &b) |
template<class T , class C , class U > | |
TemporaryMatrix< T > | operator/ (const MultiArrayView< 2, T, C > &a, PointWise< U > b) |
template<class T , class C > | |
TemporaryMatrix< T > | operator/ (const MultiArrayView< 2, T, C > &a, T b) |
template<class T , class C > | |
TemporaryMatrix< T > | operator/ (T a, const MultiArrayView< 2, T, C > &b) |
template<class T , class C > | |
TemporaryMatrix< T > | outer (const MultiArrayView< 2, T, C > &x) |
template<class T , class C1 , class C2 > | |
TemporaryMatrix< T > | outer (const MultiArrayView< 2, T, C1 > &x, const MultiArrayView< 2, T, C2 > &y) |
template<class T , class C1 , class C2 , class C3 > | |
void | outer (const MultiArrayView< 2, T, C1 > &x, const MultiArrayView< 2, T, C2 > &y, MultiArrayView< 2, T, C3 > &r) |
template<class T , int N> | |
TemporaryMatrix< T > | outer (const TinyVector< T, N > &x) |
template<class T , class C1 , class C2 > | |
TemporaryMatrix< T > | pdiv (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T , class C1 , class C2 , class C3 > | |
void | pdiv (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class T , class C1 , class C2 > | |
TemporaryMatrix< T > | pmul (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b) |
template<class T , class C1 , class C2 , class C3 > | |
void | pmul (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class POLYNOMIAL , class VECTOR > | |
bool | polynomialRealRootsEigenvalueMethod (POLYNOMIAL const &p, VECTOR &roots, bool) |
template<class POLYNOMIAL , class VECTOR > | |
bool | polynomialRootsEigenvalueMethod (POLYNOMIAL const &poly, VECTOR &roots, bool polishRoots) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | pow (MultiArrayView< 2, T, C > const &v, T exponent) |
template<... > | |
void | prepareColumns (...) |
Standardize the columns of a matrix according to given DataPreparationGoals . | |
template<... > | |
void | prepareRows (...) |
Standardize the rows of a matrix according to given DataPreparationGoals . | |
template<class T , class C1 , class C2 , class C3 > | |
bool | qrDecomposition (MultiArrayView< 2, T, C1 > const &a, MultiArrayView< 2, T, C2 > &q, MultiArrayView< 2, T, C3 > &r, double epsilon=0.0) |
template<class T , class C > | |
TemporaryMatrix< T > | repeatMatrix (MultiArrayView< 2, T, C > const &v, unsigned int verticalCount, unsigned int horizontalCount) |
template<class T , class C1 , class C2 > | |
void | repeatMatrix (MultiArrayView< 2, T, C1 > const &v, MultiArrayView< 2, T, C2 > &r, unsigned int verticalCount, unsigned int horizontalCount) |
template<class T , class C1 , class C2 , class C3 > | |
bool | reverseElimination (const MultiArrayView< 2, T, C1 > &r, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > x) |
template<class T , class C1 , class C2 , class C3 > | |
bool | ridgeRegression (MultiArrayView< 2, T, C1 > const &A, MultiArrayView< 2, T, C2 > const &b, MultiArrayView< 2, T, C3 > &x, double lambda) |
template<class T , class C1 , class C2 , class C3 , class Array > | |
bool | ridgeRegressionSeries (MultiArrayView< 2, T, C1 > const &A, MultiArrayView< 2, T, C2 > const &b, MultiArrayView< 2, T, C3 > &x, Array const &lambda) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | round (MultiArrayView< 2, T, C > const &v) |
template<class T , class C > | |
MultiArrayIndex | rowCount (const MultiArrayView< 2, T, C > &x) |
template<... > | |
void | rowStatistics (...) |
template<class T , class C > | |
MultiArrayView< 2, T, C > | rowVector (MultiArrayView< 2, T, C > const &m, MultiArrayIndex d) |
template<class T , class C > | |
MultiArrayView< 2, T, C > | rowVector (MultiArrayView< 2, T, C > const &m, MultiArrayShape< 2 >::type first, MultiArrayIndex end) |
template<class T , class C1 , class C2 > | |
void | sdiv (const MultiArrayView< 2, T, C1 > &a, T b, MultiArrayView< 2, T, C2 > &r) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | sign (MultiArrayView< 2, T, C > const &v) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | sin (MultiArrayView< 2, T, C > const &v) |
template<class T , class C1 , class C2 , class C3 , class C4 > | |
unsigned int | singularValueDecomposition (MultiArrayView< 2, T, C1 > const &A, MultiArrayView< 2, T, C2 > &U, MultiArrayView< 2, T, C3 > &S, MultiArrayView< 2, T, C4 > &V) |
template<class T , class C1 , class C2 > | |
void | smul (const MultiArrayView< 2, T, C1 > &a, T b, MultiArrayView< 2, T, C2 > &r) |
template<class T , class C2 , class C3 > | |
void | smul (T a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | sq (MultiArrayView< 2, T, C > const &v) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | sqrt (MultiArrayView< 2, T, C > const &v) |
template<class T , class ALLOC > | |
Matrix< T, ALLLOC >::SquaredNormType | squaredNorm (const Matrix< T, ALLLOC > &a) |
template<class T , class C1 , class C2 , class C3 > | |
void | sub (const MultiArrayView< 2, T, C1 > &a, const MultiArrayView< 2, T, C2 > &b, MultiArrayView< 2, T, C3 > &r) |
template<class T , class C > | |
MultiArrayView< 2, T, C > | subVector (MultiArrayView< 2, T, C > const &m, int first, int end) |
template<class T , class C1 , class C2 , class C3 > | |
bool | symmetricEigensystem (MultiArrayView< 2, T, C1 > const &a, MultiArrayView< 2, T, C2 > &ew, MultiArrayView< 2, T, C3 > &ev) |
template<class T , class C1 , class C2 , class C3 > | |
bool | symmetricEigensystemNoniterative (MultiArrayView< 2, T, C1 > const &a, MultiArrayView< 2, T, C2 > &ew, MultiArrayView< 2, T, C3 > &ev) |
template<class T , class C > | |
linalg::TemporaryMatrix< T > | tan (MultiArrayView< 2, T, C > const &v) |
template<class T , class C > | |
NumericTraits< T >::Promote | trace (MultiArrayView< 2, T, C > const &m) |
template<class T , class C1 , class C2 > | |
void | transpose (const MultiArrayView< 2, T, C1 > &v, MultiArrayView< 2, T, C2 > &r) |
template<class T , class C > | |
MultiArrayView< 2, T, StridedArrayTag > | transpose (MultiArrayView< 2, T, C > const &v) |
template<class T , class C1 , class C2 , class C3 , class C4 > | |
bool | weightedLeastSquares (MultiArrayView< 2, T, C1 > const &A, MultiArrayView< 2, T, C2 > const &b, MultiArrayView< 2, T, C3 > const &weights, MultiArrayView< 2, T, C4 > &x, std::string method="QR") |
template<class T , class C1 , class C2 , class C3 , class C4 > | |
bool | weightedRidgeRegression (MultiArrayView< 2, T, C1 > const &A, MultiArrayView< 2, T, C2 > const &b, MultiArrayView< 2, T, C3 > const &weights, MultiArrayView< 2, T, C4 > &x, double lambda) |
Linear algebra functions.
Namespace vigra/linalg
hold VIGRA's linear algebra functionality. But most of its contents is exported into namespace vigra
via using
directives.
bool symmetricEigensystem | ( | MultiArrayView< 2, T, C1 > const & | a, |
MultiArrayView< 2, T, C2 > & | ew, | ||
MultiArrayView< 2, T, C3 > & | ev | ||
) |
Compute the eigensystem of a symmetric matrix.
a is a real symmetric matrix, ew is a single-column matrix holding the eigenvalues, and ev is a matrix of the same size as a whose columns are the corresponding eigenvectors. Eigenvalues will be sorted from largest to smallest. The algorithm returns false
when it doesn't converge. It can be applied in-place, i.e. &a == &ev
is allowed. The code of this function was adapted from JAMA.
#include <vigra/eigensystem.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool symmetricEigensystemNoniterative | ( | MultiArrayView< 2, T, C1 > const & | a, |
MultiArrayView< 2, T, C2 > & | ew, | ||
MultiArrayView< 2, T, C3 > & | ev | ||
) |
Fast computation of the eigensystem of a 2x2 or 3x3 symmetric matrix.
The function works like symmetricEigensystem(), but uses fast analytic formula to avoid iterative computations.
#include <vigra/eigensystem.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool nonsymmetricEigensystem | ( | MultiArrayView< 2, T, C1 > const & | a, |
MultiArrayView< 2, std::complex< T >, C2 > & | ew, | ||
MultiArrayView< 2, T, C3 > & | ev | ||
) |
Compute the eigensystem of a square, but not necessarily symmetric matrix.
a is a real square matrix, ew is a single-column matrix holding the possibly complex eigenvalues, and ev is a matrix of the same size as a whose columns are the corresponding eigenvectors. Eigenvalues will be sorted from largest to smallest magnitude. The algorithm returns false
when it doesn't converge. It can be applied in-place, i.e. &a == &ev
is allowed. The code of this function was adapted from JAMA.
#include <vigra/eigensystem.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool polynomialRootsEigenvalueMethod | ( | POLYNOMIAL const & | poly, |
VECTOR & | roots, | ||
bool | polishRoots | ||
) |
Compute the roots of a polynomial using the eigenvalue method.
poly is a real polynomial (compatible to vigra::PolynomialView), and roots a complex valued vector (compatible to std::vector
with a value_type
compatible to the type POLYNOMIAL::Complex
) to which the roots are appended. The function calls nonsymmetricEigensystem() with the standard companion matrix yielding the roots as eigenvalues. It returns false
if it fails to converge.
#include <vigra/eigensystem.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool polynomialRealRootsEigenvalueMethod | ( | POLYNOMIAL const & | p, |
VECTOR & | roots, | ||
bool | |||
) |
Compute the real roots of a real polynomial using the eigenvalue method.
poly is a real polynomial (compatible to vigra::PolynomialView), and roots a real valued vector (compatible to std::vector
with a value_type
compatible to the type POLYNOMIAL::Real
) to which the roots are appended. The function calls polynomialRootsEigenvalueMethod() and throws away all complex roots. It returns false
if it fails to converge. The parameter polishRoots
is ignored (it is only here for syntax compatibility with polynomialRealRoots()).
#include <vigra/eigensystem.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool inverse | ( | const MultiArrayView< 2, T, C1 > & | v, |
MultiArrayView< 2, T, C2 > & | res | ||
) |
Create the inverse or pseudo-inverse of matrix v.
If the matrix v is square, res must have the same shape and will contain the inverse of v. If v is rectangular, res must have the transposed shape of v. The inverse is then computed in the least-squares sense, i.e. res will be the pseudo-inverse (Moore-Penrose inverse). The function returns true
upon success, and false
if v is not invertible (has not full rank). The inverse is computed by means of QR decomposition. This function can be applied in-place.
#include <vigra/linear_solve.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
Create the inverse or pseudo-inverse of matrix v.
The result is returned as a temporary matrix. If the matrix v is square, the result will have the same shape and contains the inverse of v. If v is rectangular, the result will have the transposed shape of v. The inverse is then computed in the least-squares sense, i.e. res will be the pseudo-inverse (Moore-Penrose inverse). The inverse is computed by means of QR decomposition. If v is not invertible, vigra::PreconditionViolation
exception is thrown. Usage:
#include <vigra/linear_solve.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
NumericTraits< T >::Promote determinant | ( | MultiArrayView< 2, T, C1 > const & | a, |
std::string | method = "default" |
||
) |
Compute the determinant of a square matrix.
method must be one of the following:
ContractViolation
exception is thrown. #include <vigra/linear_solve.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
Compute the logarithm of the determinant of a symmetric positive definite matrix.
This is useful to avoid multiplication of very large numbers in big matrices. It is implemented by means of Cholesky decomposition.
#include <vigra/linear_solve.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool choleskyDecomposition | ( | MultiArrayView< 2, T, C1 > const & | A, |
MultiArrayView< 2, T, C2 > & | L | ||
) |
Cholesky decomposition.
A must be a symmetric positive definite matrix, and L will be a lower triangular matrix, such that (up to round-off errors):
This implementation cannot be applied in-place, i.e. &L == &A
is an error. If A is not symmetric, a ContractViolation
exception is thrown. If it is not positive definite, the function returns false
.
#include <vigra/linear_solve.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool qrDecomposition | ( | MultiArrayView< 2, T, C1 > const & | a, |
MultiArrayView< 2, T, C2 > & | q, | ||
MultiArrayView< 2, T, C3 > & | r, | ||
double | epsilon = 0.0 |
||
) |
QR decomposition.
a contains the original matrix, results are returned in q and r, where q is a orthogonal matrix, and r is an upper triangular matrix, such that (up to round-off errors):
If a doesn't have full rank, the function returns false
. The decomposition is computed by householder transformations. It can be applied in-place, i.e. &a == &q
or &a == &r
are allowed.
#include <vigra/linear_solve.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool reverseElimination | ( | const MultiArrayView< 2, T, C1 > & | r, |
const MultiArrayView< 2, T, C2 > & | b, | ||
MultiArrayView< 2, T, C3 > | x | ||
) |
Deprecated, use linearSolveUpperTriangular().
bool linearSolveUpperTriangular | ( | const MultiArrayView< 2, T, C1 > & | r, |
const MultiArrayView< 2, T, C2 > & | b, | ||
MultiArrayView< 2, T, C3 > | x | ||
) |
Solve a linear system with upper-triangular coefficient matrix.
The square matrix r must be an upper-triangular coefficient matrix as can, for example, be obtained by means of QR decomposition. If r doesn't have full rank the function fails and returns false
, otherwise it returns true
. The lower triangular part of matrix r will not be touched, so it doesn't need to contain zeros.
The column vectors of matrix b are the right-hand sides of the equation (several equations with the same coefficients can thus be solved in one go). The result is returned int x, whose columns contain the solutions for the corresponding columns of b. This implementation can be applied in-place, i.e. &b == &x
is allowed. The following size requirements apply:
#include <vigra/linear_solve.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool linearSolveLowerTriangular | ( | const MultiArrayView< 2, T, C1 > & | l, |
const MultiArrayView< 2, T, C2 > & | b, | ||
MultiArrayView< 2, T, C3 > | x | ||
) |
Solve a linear system with lower-triangular coefficient matrix.
The square matrix l must be a lower-triangular coefficient matrix. If l doesn't have full rank the function fails and returns false
, otherwise it returns true
. The upper triangular part of matrix l will not be touched, so it doesn't need to contain zeros.
The column vectors of matrix b are the right-hand sides of the equation (several equations with the same coefficients can thus be solved in one go). The result is returned in x, whose columns contain the solutions for the corresponding columns of b. This implementation can be applied in-place, i.e. &b == &x
is allowed. The following size requirements apply:
#include <vigra/linear_solve.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
void choleskySolve | ( | MultiArrayView< 2, T, C1 > const & | L, |
MultiArrayView< 2, T, C2 > const & | b, | ||
MultiArrayView< 2, T, C3 > & | x | ||
) |
Solve a linear system when the Cholesky decomposition of the left hand side is given.
The square matrix L must be a lower-triangular matrix resulting from Cholesky decomposition of some positive definite coefficient matrix.
The column vectors of matrix b are the right-hand sides of the equation (several equations with the same matrix L can thus be solved in one go). The result is returned in x, whose columns contain the solutions for the corresponding columns of b. This implementation can be applied in-place, i.e. &b == &x
is allowed. The following size requirements apply:
#include <vigra/linear_solve.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool linearSolve | ( | ... | ) |
Solve a linear system.
Declarations:
A is the coefficient matrix, and the column vectors in b are the right-hand sides of the equation (so, several equations with the same coefficients can be solved in one go). The result is returned in res, whose columns contain the solutions for the corresponding columns of b. The number of columns of A must equal the number of rows of both b and res, and the number of columns of b and res must match. If right-hand-side and result are specified as TinyVector, the number of columns of A must equal N.
method must be one of the following:
Compute the solution by means of Cholesky decomposition. The coefficient matrix A must by symmetric positive definite. If this is not the case, the function returns false
.
(default) Compute the solution by means of QR decomposition. The coefficient matrix A can be square or rectangular. In the latter case, it must have more rows than columns, and the solution will be computed in the least squares sense. If A doesn't have full rank, the function returns false
.
Compute the solution by means of singular value decomposition. The coefficient matrix A can be square or rectangular. In the latter case, it must have more rows than columns, and the solution will be computed in the least squares sense. If A doesn't have full rank, the function returns false
.
A'*A*x = A'*b
. This only makes sense when the equation is to be solved in the least squares sense, i.e. when A is a rectangular matrix with more rows than columns. If A doesn't have full column rank, the function returns false
. This function can be applied in-place, i.e. &b == &res
or &A == &res
are allowed (provided they have the required shapes).
The following size requirements apply:
#include <vigra/linear_solve.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
MultiArrayIndex rowCount | ( | const MultiArrayView< 2, T, C > & | x | ) |
Number of rows of a matrix represented as a MultiArrayView<2, ...>
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
MultiArrayIndex columnCount | ( | const MultiArrayView< 2, T, C > & | x | ) |
Number of columns of a matrix represented as a MultiArrayView<2, ...>
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
MultiArrayView< 2, T, C > rowVector | ( | MultiArrayView< 2, T, C > const & | m, |
MultiArrayIndex | d | ||
) |
Create a row vector view for row d of the matrix m
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
MultiArrayView< 2, T, C > columnVector | ( | MultiArrayView< 2, T, C > const & | m, |
MultiArrayIndex | d | ||
) |
Create a column vector view for column d of the matrix m
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
void transpose | ( | const MultiArrayView< 2, T, C1 > & | v, |
MultiArrayView< 2, T, C2 > & | r | ||
) |
transpose matrix v. The result is written into r which must have the correct (i.e. transposed) shape.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
Check whether matrix m is symmetric.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
MultiArrayView< 2, T, C > rowVector | ( | MultiArrayView< 2, T, C > const & | m, |
MultiArrayShape< 2 >::type | first, | ||
MultiArrayIndex | end | ||
) |
Create a row vector view of the matrix m starting at element first and ranging to column end (non-inclusive).
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
MultiArrayView< 2, T, C > columnVector | ( | MultiArrayView< 2, T, C > const & | m, |
MultiArrayShape< 2 >::type | first, | ||
int | end | ||
) |
Create a column vector view of the matrix m starting at element first and ranging to row end (non-inclusive).
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
MultiArrayView< 2, T, C > subVector | ( | MultiArrayView< 2, T, C > const & | m, |
int | first, | ||
int | end | ||
) |
Create a sub vector view of the vector m starting at element first and ranging to row end (non-inclusive).
Note: This function may only be called when either rowCount(m) == 1
or columnCount(m) == 1
, i.e. when m really represents a vector. Otherwise, a PreconditionViolation exception is raised.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
NumericTraits< T >::Promote trace | ( | MultiArrayView< 2, T, C > const & | m | ) |
Compute the trace of a square matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
calculate the squared Frobenius norm of a matrix. Equal to the sum of squares of the matrix elements.
#include <vigra/matrix.hxx> Namespace: vigra
calculate the Frobenius norm of a matrix. Equal to the root of the sum of squares of the matrix elements.
#include <vigra/matrix.hxx> Namespace: vigra
initialize the given square matrix as an identity matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
TemporaryMatrix< T > identityMatrix | ( | MultiArrayIndex | size | ) |
create an identity matrix of the given size. Usage:
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
TemporaryMatrix< T > ones | ( | MultiArrayIndex | rows, |
MultiArrayIndex | cols | ||
) |
create matrix of ones of the given size. Usage:
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
void diagonalMatrix | ( | MultiArrayView< 2, T, C1 > const & | v, |
MultiArrayView< 2, T, C2 > & | r | ||
) |
make a diagonal matrix from a vector. The vector is given as matrix v, which must either have a single row or column. The result is written into the square matrix r.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
TemporaryMatrix< T > diagonalMatrix | ( | MultiArrayView< 2, T, C > const & | v | ) |
create a diagonal matrix from a vector. The vector is given as matrix v, which must either have a single row or column. The result is returned as a temporary matrix. Usage:
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
MultiArrayView< 2, T, StridedArrayTag > transpose | ( | MultiArrayView< 2, T, C > const & | v | ) |
create the transpose of matrix v. This does not copy any data, but only creates a transposed view to the original matrix. A copy is only made when the transposed view is assigned to another matrix. Usage:
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
TemporaryMatrix< T > joinVertically | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b | ||
) |
Create new matrix by concatenating two matrices a and b vertically, i.e. on top of each other. The two matrices must have the same number of columns. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > joinHorizontally | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b | ||
) |
Create new matrix by concatenating two matrices a and b horizontally, i.e. side by side. The two matrices must have the same number of rows. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
void repeatMatrix | ( | MultiArrayView< 2, T, C1 > const & | v, |
MultiArrayView< 2, T, C2 > & | r, | ||
unsigned int | verticalCount, | ||
unsigned int | horizontalCount | ||
) |
Initialize a matrix with repeated copies of a given matrix.
Matrix r will consist of verticalCount downward repetitions of v, and horizontalCount side-by-side repetitions. When v has size m
by n
, r must have size (m*verticalCount)
by (n*horizontalCount)
.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > repeatMatrix | ( | MultiArrayView< 2, T, C > const & | v, |
unsigned int | verticalCount, | ||
unsigned int | horizontalCount | ||
) |
Create a new matrix by repeating a given matrix.
The resulting matrix r will consist of verticalCount downward repetitions of v, and horizontalCount side-by-side repetitions, i.e. it will be of size (m*verticalCount)
by (n*horizontalCount)
when v has size m
by n
. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
void add | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b, | ||
MultiArrayView< 2, T, C3 > & | r | ||
) |
add matrices a and b. The result is written into r. All three matrices must have the same shape.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator+ | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b | ||
) |
add matrices a and b. The two matrices must have the same shape. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator+ | ( | const MultiArrayView< 2, T, C > & | a, |
T | b | ||
) |
add scalar b to every element of the matrix a. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator+ | ( | T | a, |
const MultiArrayView< 2, T, C > & | b | ||
) |
add scalar a to every element of the matrix b. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
void sub | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b, | ||
MultiArrayView< 2, T, C3 > & | r | ||
) |
subtract matrix b from a. The result is written into r. All three matrices must have the same shape.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator- | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b | ||
) |
subtract matrix b from a. The two matrices must have the same shape. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
negate matrix a. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator- | ( | const MultiArrayView< 2, T, C > & | a, |
T | b | ||
) |
subtract scalar b from every element of the matrix a. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator- | ( | T | a, |
const MultiArrayView< 2, T, C > & | b | ||
) |
subtract every element of the matrix b from scalar a. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
NormTraits< T >::SquaredNormType dot | ( | const MultiArrayView< 2, T, C1 > & | x, |
const MultiArrayView< 2, T, C2 > & | y | ||
) |
calculate the inner product of two matrices representing vectors. Typically, matrix x has a single row, and matrix y has a single column, and the other dimensions match. In addition, this function handles the cases when either or both of the two inputs are transposed (e.g. it can compute the dot product of two column vectors). A PreconditionViolation
exception is thrown when the shape conditions are violated.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
NormTraits< T >::SquaredNormType dot | ( | const MultiArrayView< 1, T, C1 > & | x, |
const MultiArrayView< 1, T, C2 > & | y | ||
) |
calculate the inner product of two vectors. The vector lengths must match.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
void cross | ( | const MultiArrayView< 1, T, C1 > & | x, |
const MultiArrayView< 1, T, C2 > & | y, | ||
MultiArrayView< 1, T, C3 > & | r | ||
) |
calculate the cross product of two vectors of length 3. The result is written into r.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
void cross | ( | const MultiArrayView< 2, T, C1 > & | x, |
const MultiArrayView< 2, T, C2 > & | y, | ||
MultiArrayView< 2, T, C3 > & | r | ||
) |
calculate the cross product of two matrices representing vectors. That is, x, y, and r must have a single column of length 3. The result is written into r.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
TemporaryMatrix< T > cross | ( | const MultiArrayView< 2, T, C1 > & | x, |
const MultiArrayView< 2, T, C2 > & | y | ||
) |
calculate the cross product of two matrices representing vectors. That is, x, and y must have a single column of length 3. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
void outer | ( | const MultiArrayView< 2, T, C1 > & | x, |
const MultiArrayView< 2, T, C2 > & | y, | ||
MultiArrayView< 2, T, C3 > & | r | ||
) |
calculate the outer product of two matrices representing vectors. That is, matrix x must have a single column, and matrix y must have a single row, and the other dimensions must match. The result is written into r.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
TemporaryMatrix< T > outer | ( | const MultiArrayView< 2, T, C1 > & | x, |
const MultiArrayView< 2, T, C2 > & | y | ||
) |
calculate the outer product of two matrices representing vectors. That is, matrix x must have a single column, and matrix y must have a single row, and the other dimensions must match. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
calculate the outer product of a matrix (representing a vector) with itself. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
TemporaryMatrix< T > outer | ( | const TinyVector< T, N > & | x | ) |
calculate the outer product of a TinyVector with itself. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
void smul | ( | const MultiArrayView< 2, T, C1 > & | a, |
T | b, | ||
MultiArrayView< 2, T, C2 > & | r | ||
) |
multiply matrix a with scalar b. The result is written into r. a and r must have the same shape.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
void smul | ( | T | a, |
const MultiArrayView< 2, T, C2 > & | b, | ||
MultiArrayView< 2, T, C3 > & | r | ||
) |
multiply scalar a with matrix b. The result is written into r. b and r must have the same shape.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
void mmul | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b, | ||
MultiArrayView< 2, T, C3 > & | r | ||
) |
perform matrix multiplication of matrices a and b. The result is written into r. The three matrices must have matching shapes.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > mmul | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b | ||
) |
perform matrix multiplication of matrices a and b. a and b must have matching shapes. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
void pmul | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b, | ||
MultiArrayView< 2, T, C3 > & | r | ||
) |
multiply two matrices a and b pointwise. The result is written into r. All three matrices must have the same shape.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > pmul | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b | ||
) |
multiply matrices a and b pointwise. a and b must have matching shapes. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator* | ( | const MultiArrayView< 2, T, C > & | a, |
PointWise< U > | b | ||
) |
multiply matrices a and b pointwise. a and b must have matching shapes. The result is returned as a temporary matrix.
Usage:
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator* | ( | const MultiArrayView< 2, T, C > & | a, |
T | b | ||
) |
multiply matrix a with scalar b. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator* | ( | T | a, |
const MultiArrayView< 2, T, C > & | b | ||
) |
multiply scalar a with matrix b. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TinyVector< T, N > operator* | ( | const Matrix< T, A > & | a, |
const TinyVectorBase< T, N, DATA, DERIVED > & | b | ||
) |
multiply matrix a with TinyVector b. a must be of size N x N
. Vector b and the result vector are interpreted as column vectors.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TinyVector< T, N > operator* | ( | const TinyVectorBase< T, N, DATA, DERIVED > & | a, |
const Matrix< T, A > & | b | ||
) |
multiply TinyVector a with matrix b. b must be of size N x N
. Vector a and the result vector are interpreted as row vectors.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator* | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b | ||
) |
perform matrix multiplication of matrices a and b. a and b must have matching shapes. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
void sdiv | ( | const MultiArrayView< 2, T, C1 > & | a, |
T | b, | ||
MultiArrayView< 2, T, C2 > & | r | ||
) |
divide matrix a by scalar b. The result is written into r. a and r must have the same shape.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
void pdiv | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b, | ||
MultiArrayView< 2, T, C3 > & | r | ||
) |
divide two matrices a and b pointwise. The result is written into r. All three matrices must have the same shape.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > pdiv | ( | const MultiArrayView< 2, T, C1 > & | a, |
const MultiArrayView< 2, T, C2 > & | b | ||
) |
divide matrices a and b pointwise. a and b must have matching shapes. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator/ | ( | const MultiArrayView< 2, T, C > & | a, |
PointWise< U > | b | ||
) |
divide matrices a and b pointwise. a and b must have matching shapes. The result is returned as a temporary matrix.
Usage:
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator/ | ( | const MultiArrayView< 2, T, C > & | a, |
T | b | ||
) |
divide matrix a by scalar b. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
TemporaryMatrix< T > operator/ | ( | T | a, |
const MultiArrayView< 2, T, C > & | b | ||
) |
Create a matrix whose elements are the quotients between scalar a and matrix b. The result is returned as a temporary matrix.
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespace: vigra::linalg
Find the index of the minimum element in a matrix.
The function returns the index in column-major scan-order sense, i.e. according to the order used by MultiArrayView::operator[]
. If the matrix is actually a vector, this is just the row or columns index. In case of a truly 2-dimensional matrix, the index can be converted to an index pair by calling MultiArrayView::scanOrderIndexToCoordinate()
Required Interface:
#include <vigra/matrix.hxx>
Namespace: vigra
Find the index of the maximum element in a matrix.
The function returns the index in column-major scan-order sense, i.e. according to the order used by MultiArrayView::operator[]
. If the matrix is actually a vector, this is just the row or columns index. In case of a truly 2-dimensional matrix, the index can be converted to an index pair by calling MultiArrayView::scanOrderIndexToCoordinate()
Required Interface:
#include <vigra/matrix.hxx>
Namespace: vigra
int argMinIf | ( | MultiArrayView< 2, T, C > const & | a, |
UnaryFunctor | condition | ||
) |
Find the index of the minimum element in a matrix subject to a condition.
The function returns -1
if no element conforms to condition. Otherwise, the index of the maximum element is returned in column-major scan-order sense, i.e. according to the order used by MultiArrayView::operator[]
. If the matrix is actually a vector, this is just the row or columns index. In case of a truly 2-dimensional matrix, the index can be converted to an index pair by calling MultiArrayView::scanOrderIndexToCoordinate()
Required Interface:
#include <vigra/matrix.hxx>
Namespace: vigra
int argMaxIf | ( | MultiArrayView< 2, T, C > const & | a, |
UnaryFunctor | condition | ||
) |
Find the index of the maximum element in a matrix subject to a condition.
The function returns -1
if no element conforms to condition. Otherwise, the index of the maximum element is returned in column-major scan-order sense, i.e. according to the order used by MultiArrayView::operator[]
. If the matrix is actually a vector, this is just the row or columns index. In case of a truly 2-dimensional matrix, the index can be converted to an index pair by calling MultiArrayView::scanOrderIndexToCoordinate()
Required Interface:
#include <vigra/matrix.hxx>
Namespace: vigra
linalg::TemporaryMatrix< T > pow | ( | MultiArrayView< 2, T, C > const & | v, |
T | exponent | ||
) |
Matrix point-wise power.
linalg::TemporaryMatrix< T > sqrt | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise sqrt.
linalg::TemporaryMatrix< T > exp | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise exp.
linalg::TemporaryMatrix< T > log | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise log.
linalg::TemporaryMatrix< T > log10 | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise log10.
linalg::TemporaryMatrix< T > sin | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise sin.
linalg::TemporaryMatrix< T > asin | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise asin.
linalg::TemporaryMatrix< T > cos | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise cos.
linalg::TemporaryMatrix< T > acos | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise acos.
linalg::TemporaryMatrix< T > tan | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise tan.
linalg::TemporaryMatrix< T > atan | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise atan.
linalg::TemporaryMatrix< T > round | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise round.
linalg::TemporaryMatrix< T > floor | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise floor.
linalg::TemporaryMatrix< T > ceil | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise ceil.
linalg::TemporaryMatrix< T > abs | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise abs.
Matrix point-wise square.
linalg::TemporaryMatrix< T > sign | ( | MultiArrayView< 2, T, C > const & | v | ) |
Matrix point-wise sign.
void columnStatistics | ( | ... | ) |
Compute statistics of every column of matrix A.
The result matrices must be row vectors with as many columns as A.
Declarations:
compute only the mean:
compute mean and standard deviation:
compute mean, standard deviation, and norm:
Usage:
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
void rowStatistics | ( | ... | ) |
Compute statistics of every row of matrix A.
The result matrices must be column vectors with as many rows as A.
Declarations:
compute only the mean:
compute mean and standard deviation:
compute mean, standard deviation, and norm:
Usage:
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
void covarianceMatrixOfColumns | ( | MultiArrayView< 2, T1, C1 > const & | features, |
MultiArrayView< 2, T2, C2 > & | covariance | ||
) |
Compute the covariance matrix between the columns of a matrix features.
The result matrix covariance must by a square matrix with as many rows and columns as the number of columns in matrix features.
#include <vigra/matrix.hxx>
Namespace: vigra
TemporaryMatrix< T > covarianceMatrixOfColumns | ( | MultiArrayView< 2, T, C > const & | features | ) |
Compute the covariance matrix between the columns of a matrix features.
The result is returned as a square temporary matrix with as many rows and columns as the number of columns in matrix features.
#include <vigra/matrix.hxx>
Namespace: vigra
void covarianceMatrixOfRows | ( | MultiArrayView< 2, T1, C1 > const & | features, |
MultiArrayView< 2, T2, C2 > & | covariance | ||
) |
Compute the covariance matrix between the rows of a matrix features.
The result matrix covariance must by a square matrix with as many rows and columns as the number of rows in matrix features.
#include <vigra/matrix.hxx>
Namespace: vigra
TemporaryMatrix< T > covarianceMatrixOfRows | ( | MultiArrayView< 2, T, C > const & | features | ) |
Compute the covariance matrix between the rows of a matrix features.
The result is returned as a square temporary matrix with as many rows and columns as the number of rows in matrix features.
#include <vigra/matrix.hxx>
Namespace: vigra
void prepareColumns | ( | ... | ) |
Standardize the columns of a matrix according to given DataPreparationGoals
.
For every column of the matrix A, this function computes mean, standard deviation, and norm. It then applies a linear transformation to the values of the column according to these statistics and the given DataPreparationGoals
. The result is returned in matrix res which must have the same size as A. Optionally, the transformation applied can also be returned in the matrices offset and scaling (see below for an example how these matrices can be used to standardize more data according to the same transformation).
The following DataPreparationGoals
are supported:
ZeroMean
UnitSum
UnitVariance
UnitNorm
ZeroMean | UnitVariance
ZeroMean | UnitNorm
Declarations:
Standardize the matrix and return the parameters of the linear transformation. The matrices offset and scaling must be row vectors with as many columns as A.
Only standardize the matrix.
Usage:
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
void prepareRows | ( | ... | ) |
Standardize the rows of a matrix according to given DataPreparationGoals
.
This algorithm works in the same way as prepareColumns() (see there for detailed documentation), but is applied to the rows of the matrix A instead. Accordingly, the matrices holding the parameters of the linear transformation must be column vectors with as many rows as A.
Declarations:
Standardize the matrix and return the parameters of the linear transformation. The matrices offset and scaling must be column vectors with as many rows as A.
Only standardize the matrix.
Usage:
#include <vigra/matrix.hxx> or
#include <vigra/linear_algebra.hxx>
Namespaces: vigra and vigra::linalg
bool leastSquares | ( | MultiArrayView< 2, T, C1 > const & | A, |
MultiArrayView< 2, T, C2 > const & | b, | ||
MultiArrayView< 2, T, C3 > & | x, | ||
std::string | method = "QR" |
||
) |
Ordinary Least Squares Regression.
Given a matrix A with m
rows and n
columns (with m >= n
), and a column vector b of length m
rows, this function computes the column vector x of length n
rows that solves the optimization problem
When b is a matrix with k
columns, x must also have k
columns, which will contain the solutions for the corresponding columns of b. Note that all matrices must already have the correct shape.
This function is just another name for linearSolve(), perhaps leading to more readable code when A is a rectangular matrix. It returns false
when the rank of A is less than n
. See linearSolve() for more documentation.
#include <vigra/regression.hxx>
Namespaces: vigra and vigra::linalg
bool weightedLeastSquares | ( | MultiArrayView< 2, T, C1 > const & | A, |
MultiArrayView< 2, T, C2 > const & | b, | ||
MultiArrayView< 2, T, C3 > const & | weights, | ||
MultiArrayView< 2, T, C4 > & | x, | ||
std::string | method = "QR" |
||
) |
Weighted Least Squares Regression.
Given a matrix A with m
rows and n
columns (with m >= n
), a vector b of length m
, and a weight vector weights of length m
with non-negative entries, this function computes the vector x of length n
that solves the optimization problem
where diag(weights)
creates a diagonal matrix from weights. The algorithm calls leastSquares() on the equivalent problem
where the square root of weights is just taken element-wise.
When b is a matrix with k
columns, x must also have k
columns, which will contain the solutions for the corresponding columns of b. Note that all matrices must already have the correct shape.
The function returns false
when the rank of the weighted matrix A is less than n
.
#include <vigra/regression.hxx>
Namespaces: vigra and vigra::linalg
bool ridgeRegression | ( | MultiArrayView< 2, T, C1 > const & | A, |
MultiArrayView< 2, T, C2 > const & | b, | ||
MultiArrayView< 2, T, C3 > & | x, | ||
double | lambda | ||
) |
Ridge Regression.
Given a matrix A with m
rows and n
columns (with m >= n
), a vector b of length m
, and a regularization parameter lambda >= 0.0
, this function computes the vector x of length n
that solves the optimization problem
This is implemented by means of singularValueDecomposition().
When b is a matrix with k
columns, x must also have k
columns, which will contain the solutions for the corresponding columns of b. Note that all matrices must already have the correct shape.
The function returns false
if the rank of A is less than n
and lambda == 0.0
.
#include <vigra/regression.hxx>
Namespaces: vigra and vigra::linalg
bool weightedRidgeRegression | ( | MultiArrayView< 2, T, C1 > const & | A, |
MultiArrayView< 2, T, C2 > const & | b, | ||
MultiArrayView< 2, T, C3 > const & | weights, | ||
MultiArrayView< 2, T, C4 > & | x, | ||
double | lambda | ||
) |
Weighted ridge Regression.
Given a matrix A with m
rows and n
columns (with m >= n
), a vector b of length m
, a weight vector weights of length m
with non-negative entries, and a regularization parameter lambda >= 0.0
this function computes the vector x of length n
that solves the optimization problem
where diag(weights)
creates a diagonal matrix from weights. The algorithm calls ridgeRegression() on the equivalent problem
where the square root of weights is just taken element-wise. This solution is computed by means of singularValueDecomposition().
When b is a matrix with k
columns, x must also have k
columns, which will contain the solutions for the corresponding columns of b. Note that all matrices must already have the correct shape.
The function returns false
if the rank of A is less than n
and lambda == 0.0
.
#include <vigra/regression.hxx>
Namespaces: vigra and vigra::linalg
bool ridgeRegressionSeries | ( | MultiArrayView< 2, T, C1 > const & | A, |
MultiArrayView< 2, T, C2 > const & | b, | ||
MultiArrayView< 2, T, C3 > & | x, | ||
Array const & | lambda | ||
) |
Ridge Regression with many lambdas.
This executes ridgeRegression() for a sequence of regularization parameters. This is implemented so that the singularValueDecomposition() has to be executed only once. lambda must be an array conforming to the std::vector
interface, i.e. must support lambda.size()
and lambda[k]
. The columns of the matrix x will contain the solutions for the corresponding lambda, so the number of columns of the matrix x must be equal to lambda.size()
, and b must be a columns vector, i.e. cannot contain several right hand sides at once.
The function returns false
when the matrix A is rank deficient. If this happens, and one of the lambdas is zero, the corresponding column of x will be skipped.
#include <vigra/regression.hxx>
Namespaces: vigra and vigra::linalg
Least Angle Regression.
#include <vigra/regression.hxx>
Namespaces: vigra and vigra::linalg
Declarations:
This function implements Least Angle Regression (LARS) as described in
B.Efron, T.Hastie, I.Johnstone, and R.Tibshirani: "Least Angle Regression", Annals of Statistics 32(2):407-499, 2004.
It is an efficient algorithm to solve the L1-regularized least squares (LASSO) problem
and the L1-regularized non-negative least squares (NN-LASSO) problem
where A is a matrix with m
rows and n
columns (often with m < n
), b a vector of length m
, and a regularization parameter s >= 0.0. L1-regularization has the desirable effect that it causes the solution x to be sparse, i.e. only the most important elements in x (called the active set) have non-zero values. The key insight of the LARS algorithm is the following: When the solution vector x is considered as a function of the regularization parameter s, then x(s) is a piecewise linear function, i.e. a polyline in n-dimensional space. The knots of the polyline x(s) are located precisely at those values of s where one variable enters or leaves the active set and can be efficiently computed.
Therefore, leastAngleRegression() returns the entire solution path as a sequence of knot points, starting at
The sequences of active sets and corresponding variable weights are returned in activeSets and solutions respectively. That is, activeSets[i]
is an ArrayVector<int> containing the indices of the variables that are active at the i-th knot, and solutions
is a Matrix<T> containing the weights of those variables, in the same order (see example below). Variables not contained in activeSets[i]
are zero at this solution.
The behavior of the algorithm can be adapted by LeastAngleRegressionOptions:
n
solutions. Usage:
Required Interface:
T
must be numeric type (compatible to double) Array1 a1;
a1.push_back(ArrayVector<int>());
Array2 a2;
a2.push_back(Matrix<T>());
Non-negative Least Squares Regression.
Given a matrix A with m
rows and n
columns (with m >= n
), and a column vector b of length m
rows, this function computes a column vector x of length n
with non-negative entries that solves the optimization problem
Both b and x must be column vectors (i.e. matrices with 1
column). Note that all matrices must already have the correct shape. The solution is computed by means of leastAngleRegression() with non-negativity constraint.
#include <vigra/regression.hxx>
Namespaces: vigra and vigra::linalg
Declarations:
© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de) |