![]() |
|
Macros | |
#define | Check_vw |
Functions | |
long double | apop_vector_sum (const gsl_vector *in) |
double | apop_vector_skew (const gsl_vector *in) |
double | apop_vector_kurtosis (const gsl_vector *in) |
double | apop_vector_skew_pop (gsl_vector const *v, gsl_vector const *weights) |
double | apop_vector_kurtosis_pop (gsl_vector const *v, gsl_vector const *weights) |
double | apop_vector_var_m (const gsl_vector *in, const double mean) |
double | apop_vector_correlation (const gsl_vector *ina, const gsl_vector *inb) |
double | apop_vector_distance (const gsl_vector *ina, const gsl_vector *inb, const char metric, const double norm) |
void | apop_vector_normalize (gsl_vector *in, gsl_vector **out, const char normalization_type) |
void | apop_matrix_normalize (gsl_matrix *data, const char row_or_col, const char normalization) |
long double | apop_matrix_sum (const gsl_matrix *m) |
double | apop_matrix_mean (const gsl_matrix *data) |
void | apop_matrix_mean_and_var (const gsl_matrix *data, double *mean, double *var) |
apop_data * | apop_data_summarize (apop_data *indata) |
double | apop_vector_mean (gsl_vector const *v, gsl_vector const *weights) |
double | apop_vector_var (gsl_vector const *v, gsl_vector const *weights) |
double | apop_vector_cov (const gsl_vector *v1, const gsl_vector *v2, const gsl_vector *weights) |
apop_data * | apop_data_covariance (const apop_data *in) |
apop_data * | apop_data_correlation (const apop_data *in) |
double | apop_kl_divergence (apop_model *from, apop_model *to, int draw_ct, gsl_rng *rng) |
long double | apop_multivariate_gamma (double a, int p) |
long double | apop_multivariate_lngamma (double a, int p) |
int | apop_matrix_is_positive_semidefinite (gsl_matrix *m, char semi) |
void | vfabs (double *x) |
double | apop_matrix_to_positive_semidefinite (gsl_matrix *m) |
Basic moments and some distributions.
#define Check_vw |
Returns the matrix of correlation coefficients relating each column with each other.
in | A data matrix: rows are observations, columns are variables. If you give me a weights vector, I'll use it. |
out->error='a' | Allocation error. |
Returns the sample variance/covariance matrix relating each column of the matrix to each other column.
in | An apop_data set. If the weights vector is set, I'll take it into account. |
out->error='a' | Allocation error. |
double apop_kl_divergence | ( | apop_model * | from, |
apop_model * | to, | ||
int | draw_ct, | ||
gsl_rng * | rng | ||
) |
Kullback-Leibler divergence.
This measure of the divergence of one distribution from another has the form . Notice that it is not a distance, because there is an asymmetry between
and
, so one can expect that
.
from | the ![]() NULL ) |
to | the ![]() NULL ) |
draw_ct | If I do the calculation via random draws, how many? (Default = 1e5) |
rng | A gsl_rng . If NULL or number of threads is greater than 1, I'll take care of the RNG; see apop_rng_get_thread. (Default = NULL ) |
This function can take empirical histogram-type models (apop_pmf) or continuous models like apop_loess or apop_normal.
If there is a PMF (I'll try from
first, under the presumption that you are measuring the divergence of data from an observed data distribution), then I'll step through it for the points in the summation.
GSL_NEGINF
. If apop_opts.verbose >=1
I print a message as well.If neither distribution is a PMF, then I'll take draw_ct
random draws from to
and evaluate at those points.
apop_opts.verbose = 3
for observation-by-observation info.int apop_matrix_is_positive_semidefinite | ( | gsl_matrix * | m, |
char | semi | ||
) |
Test whether the input matrix is positive semidefinite.
A covariance matrix will always be PSD, so this function can tell you whether your matrix is a valid covariance matrix.
Consider the 1x1 matrix in the upper left of the input, then the 2x2 matrix in the upper left, on up to the full matrix. If the matrix is PSD, then each of these has a positive determinant. This function thus calculates determinants for an
x
matrix.
m | The matrix to test. If NULL , I will return zero—not PSD. |
semi | If anything but 's', check for positive definite, not semidefinite. (default 's') |
See also apop_matrix_to_positive_semidefinite, which will change the input to something PSD.
void apop_matrix_normalize | ( | gsl_matrix * | data, |
const char | row_or_col, | ||
const char | normalization | ||
) |
Normalize each row or column in the given matrix, one by one.
Basically just a convenience fn to iterate through the columns or rows and run apop_vector_normalize for you.
data | The data set to normalize. |
row_or_col | Either 'r' or 'c'. |
normalization | see apop_vector_normalize. |
double apop_matrix_to_positive_semidefinite | ( | gsl_matrix * | m | ) |
First, this function passes tests, but is under development.
It takes in a matrix and converts it to the `closest' positive semidefinite matrix.
m | On input, any matrix; on output, a positive semidefinite matrix. |
Adapted from the R Matrix package's nearPD, which is Copyright (2007) Jens Oehlschlägel [and is GPL].
long double apop_multivariate_gamma | ( | double | a, |
int | p | ||
) |
The multivariate generalization of the Gamma distribution.
Because is undefined for
, this function returns
NAN
when takes on one of those values.
See also apop_multivariate_lngamma, which is more numerically stable in most cases.
long double apop_multivariate_lngamma | ( | double | a, |
int | p | ||
) |
The log of the multivariate generalization of the Gamma; see also apop_multivariate_gamma.
double apop_vector_cov | ( | const gsl_vector * | v1, |
const gsl_vector * | v2, | ||
const gsl_vector * | weights | ||
) |
Find the sample covariance of a pair of vectors, with an optional weighting. This only makes sense if the weightings are identical, so the function takes only one weighting vector for both.
v1,v2 | The data vectors |
weights | The weight vector. Default: equal weights for all elements. |
double apop_vector_mean | ( | gsl_vector const * | v, |
gsl_vector const * | weights | ||
) |
Find the mean, weighted or unweighted.
v | The data vector |
weights | The weight vector. Default: assume equal weights. |
void apop_vector_normalize | ( | gsl_vector * | in, |
gsl_vector ** | out, | ||
const char | normalization_type | ||
) |
This function will normalize a vector, either such that it has mean zero and variance one, or ranges between zero and one, or sums to one.
in | A gsl_vector which you have already allocated and filled. NULL input gives NULL output. (No default) |
out | If normalizing in place, NULL . If not, the address of a gsl_vector . Do not allocate. (default = NULL .) |
normalization_type | 'p': normalized vector will sum to one. E.g., start with a set of observations in bins, end with the percentage of observations in each bin. (the default) 'r': normalized vector will range between zero and one. Replace each X with (X-min) / (max - min). 's': normalized vector will have mean zero and variance one. Replace each X with ![]() ![]() 'm': normalize to mean zero: Replace each X with ![]() |
Example
double apop_vector_var | ( | gsl_vector const * | v, |
gsl_vector const * | weights | ||
) |
Find the sample variance of a vector, weighted or unweighted.
v | The data vector |
weights | The weight vector. If NULL, assume equal weights. |
w->size
as the number of elements, and returns the usual sum over