![]() |
|
Functions | |
double | apop_det_and_inv (const gsl_matrix *in, gsl_matrix **out, int calc_det, int calc_inv) |
gsl_matrix * | apop_matrix_inverse (const gsl_matrix *in) |
double | apop_matrix_determinant (const gsl_matrix *in) |
apop_data * | apop_matrix_pca (gsl_matrix *data, int const dimensions_we_want) |
apop_data * | apop_dot (const apop_data *d1, const apop_data *d2, char form1, char form2) |
gsl_vector * | apop_numerical_gradient (apop_data *data, apop_model *model, double delta) |
This page describes some standard bits of linear algebra that Apophenia facilitates.
See also the printing functions, Assorted printing functions, and the Convenience functions.
double apop_det_and_inv | ( | const gsl_matrix * | in, |
gsl_matrix ** | out, | ||
int | calc_det, | ||
int | calc_inv | ||
) |
Calculate the determinant of a matrix, its inverse, or both, via LU decomposition. The in
matrix is not destroyed in the process.
in | The matrix to be inverted/determined. |
out | If you want an inverse, this is where to place the matrix to be filled with the inverse. Will be allocated by the function. |
calc_det | 0: Do not calculate the determinant. 1: Do. |
calc_inv | 0: Do not calculate the inverse. 1: Do. |
calc_det == 1
, then return the determinant. Otherwise, just returns zero. If calc_inv!=0
, then *out
is pointed to the matrix inverse. In case of difficulty, I will set *out=NULL
and return NaN
. A convenience function for dot products, which requires less prep and typing than the gsl_cblas_dgexx
functions.
Second, it makes some use of the semi-overloading of the apop_data structure. d1
may be a vector or a matrix, and the same for d2
, so this function can do vector dot matrix, matrix dot matrix, and so on. If d1
includes both a vector and a matrix, then later parameters will indicate which to use.
d1 | the left part of ![]() |
d2 | the right part of ![]() |
form1 | 't' or 'p': transpose or prime d1->matrix , or, if d1->matrix is NULL , read d1->vector as a row vector.'n' or 0: no transpose. (the default) 'v': ignore the matrix and use the vector. |
form2 | As above, with d2 . |
NULL
and the matrix has the dot product; if either or both are vectors, the vector has the output and the matrix is NULL
.out->error='a' | Allocation error. |
out->error='d' | dimension-matching error. |
out->error='m' | GSL math error. |
NULL | If you ask me to take the dot product of NULL, I return NULL. [May some day change.] |
d1
is a vector and d2
a matrix, then apop_dot(d1,d2,'t')
won't work, because 't'
now refers to d1
. Instead use apop_dot(d1,d2,.form2='t')
or apop_dot(d1,d2,0, 't')
Sample code:
double apop_matrix_determinant | ( | const gsl_matrix * | in | ) |
Find the determinant of a matrix. The in
matrix is not destroyed in the process.
See also apop_matrix_inverse , or apop_det_and_inv to do both at once.
in | The matrix to be determined. |
gsl_matrix* apop_matrix_inverse | ( | const gsl_matrix * | in | ) |
Inverts a matrix. The in
matrix is not destroyed in the process. You may want to call apop_matrix_determinant first to check that your input is invertible, or use apop_det_and_inv to do both at once.
in | The matrix to be inverted. |
apop_data* apop_matrix_pca | ( | gsl_matrix * | data, |
int const | dimensions_we_want | ||
) |
Principal component analysis: hand in a matrix and (optionally) a number of desired dimensions, and I'll return a data set where each column of the matrix is an eigenvector. The columns are sorted, so column zero has the greatest weight. The vector element of the data set gives the weights.
You also specify the number of elements your principal component space should have. If this is equal to the rank of the space in which the input data lives, then the sum of weights will be one. If the dimensions desired is less than that (probably so you can prepare a plot), then the weights will be accordingly smaller, giving you an indication of how much variation these dimensions explain.
data | The input matrix. (No default. If NULL , return NULL and print a warning iff apop_opts.verbose >= 1 .) I modify int in place so that each column has mean zero. |
dimensions_we_want | (default: the size of the covariance matrix, i.e. data->size2 ) The singular value decomposition will return this many of the eigenvectors with the largest eigenvalues. |
out->error=='a' | Allocation error. |
gsl_vector* apop_numerical_gradient | ( | apop_data * | data, |
apop_model * | model, | ||
double | delta | ||
) |
The GSL provides one-dimensional numerical differentiation; here's the multidimensional extension.
data | The data set to use for all evaluations. It remains constant throughout. |
model | The model, expressing the function whose derivative is sought. The gradient is taken via small changes along the model parameters. |
delta | The size of the differential. If you explicitly give me a delta , I'll use it. If delta is not specified, but model has method_settings of type apop_ml_params , then the delta element is used for the differential. Else, I use 1e-3. |