Patterns in static

Apophenia

Macros | Functions
Calculate moments (mean, var, kurtosis) for the data in a gsl_vector.

Macros

#define apop_mean(in)
 
#define apop_var(in)   apop_vector_var(in)
 

Functions

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)
 

Detailed Description

These functions simply take in a GSL vector and return its mean, variance, or kurtosis; the covariance functions take two GSL vectors as inputs.

See also
db_moments

For apop_vector_var_m(vector, mean), mean is the mean of the vector. This saves the trouble of re-calculating the mean if you've already done so. E.g.,

gsl_vector *v;
double mean, var;
//Allocate v and fill it with data here.
mean = apop_vector_mean(v);
var = apop_vector_var_m(v, mean);
printf("Your vector has mean %g and variance %g\n", mean, var);

Macro Definition Documentation

#define apop_mean (   in)

Returns the mean of the elements of the vector v.

An alias for apop_vector_mean. Returns the mean of the data in the given vector.

#define apop_var (   in)    apop_vector_var(in)

An alias for apop_vector_var. Returns the variance of the data in the given vector.

Function Documentation

double apop_vector_correlation ( const gsl_vector *  ina,
const gsl_vector *  inb 
)

Returns the correlation coefficient of two vectors. It's just $ {\hbox{cov}(a,b)\over \sqrt(\hbox{var}(a)) * \sqrt(\hbox{var}(b))}.$

double apop_vector_kurtosis ( const gsl_vector *  in)

Returns the sample kurtosis (divide by $n-1$) of the data in the given vector. Corrections are made to produce an unbiased result as per Appendix M (PDF) of Modeling with data.

  • This does not normalize the output: the kurtosis of a ${\cal N}(0,1)$ is $3 \sigma^4$, not three, one, or zero.
double apop_vector_kurtosis_pop ( gsl_vector const *  v,
gsl_vector const *  weights 
)

Returns the population kurtosis ( $\sum_i (x_i - \mu)^4/n)$) of the data in the given vector, with an optional weighting.

Parameters
vThe data vector
weightsThe weight vector. If NULL, assume equal weights.
Returns
The weighted kurtosis. No sample adjustment given weights.
  • Some people like to normalize the kurtosis by dividing by variance squared, or by subtracting three; those things are not done here, so you'll have to do them separately if need be.
  • This function uses the Designated initializers syntax for inputs.
double apop_vector_skew ( const gsl_vector *  in)

Returns an unbiased estimate of the sample skew (population skew times y $n^2/(n^2-1)$) of the data in the given vector.

double apop_vector_skew_pop ( gsl_vector const *  v,
gsl_vector const *  weights 
)

Returns the population skew $(\sum_i (x_i - \mu)^3/n))$ of the data in the given vector. Observations may be weighted.

Parameters
vThe data vector
weightsThe weight vector. Default: equal weights for all observations.
Returns
The weighted skew. No sample adjustment given weights.
  • Some people like to normalize the skew by dividing by variance $^{3/2}$; that's not done here, so you'll have to do so separately if need be.
  • Apophenia tries to be smart about reading the weights. If weights sum to one, then the system uses w->size as the number of elements, and returns the usual sum over $n-1$. If weights > 1, then the system uses the total weights as $n$. Thus, you can use the weights as standard weightings or to represent elements that appear repeatedly.
double apop_vector_var_m ( const gsl_vector *  in,
const double  mean 
)

Returns the variance of the data in the given vector, given that you've already calculated the mean.

Parameters
inthe vector in question
meanthe mean, which you've already calculated using apop_vector_mean.

Autogenerated by doxygen on Sun Oct 26 2014 (Debian 0.999b+ds3-2).