![]() |
Object-oriented Scientific Computing Library: Version 0.910
|
File containing statistics template functions. More...
#include <o2scl/err_hnd.h>
#include <gsl/gsl_sort.h>
Go to the source code of this file.
This file contains several function templates for computing statistics of vectors of double-precision data. It includes mean, median, variance, standard deviation, covariance, correlation, and other functions.
No additional range checking is done on the vectors.
Definition in file vec_stats.h.
Functions | |
Vector functions | |
template<class vec_t > | |
double | vector_mean (size_t n, const vec_t &data) |
Compute the mean of the first n elements of a vector. | |
template<class vec_t > | |
double | vector_variance_fmean (size_t n, const vec_t &data, double mean) |
Compute variance with specified mean known in advance. | |
template<class vec_t > | |
double | vector_variance (size_t n, const vec_t &data, double mean) |
Compute the variance with specified mean. | |
template<class vec_t > | |
double | vector_variance (size_t n, const vec_t &data) |
Compute the variance. | |
template<class vec_t > | |
double | vector_stddev_fmean (size_t n, const vec_t &data, double mean) |
Standard deviation with specified mean known in advance. | |
template<class vec_t > | |
double | vector_stddev (size_t n, const vec_t &data) |
Standard deviation with specified mean. | |
template<class vec_t > | |
double | vector_stddev (size_t n, const vec_t &data, double mean) |
Standard deviation with specified mean. | |
template<class vec_t > | |
double | vector_absdev (size_t n, const vec_t &data, double mean) |
Absolute deviation from the specified mean. | |
template<class vec_t > | |
double | vector_absdev (size_t n, const vec_t &data) |
Absolute deviation from the computed mean. | |
template<class vec_t > | |
double | vector_skew (size_t n, const vec_t &data, double mean, double stddev) |
Skewness with specified mean and standard deviation. | |
template<class vec_t > | |
double | vector_skew (size_t n, const vec_t &data) |
Skewness with computed mean and standard deviation. | |
template<class vec_t > | |
double | vector_kurtosis (size_t n, const vec_t &data, double mean, double stddev) |
Kurtosis with specified mean and standard deviation. | |
template<class vec_t > | |
double | vector_kurtosis (size_t n, const vec_t &data) |
Kurtosis with computed mean and standard deviation. | |
template<class vec_t > | |
double | vector_lag1_autocorr (size_t n, const vec_t &data, double mean) |
Lag-1 autocorrelation. | |
template<class vec_t > | |
double | vector_lag1_autocorr (size_t n, const vec_t &data) |
Lag-1 autocorrelation. | |
template<class vec_t > | |
double | vector_covariance (size_t n, const vec_t &data1, const vec_t &data2, double mean1, double mean2) |
Compute the covariance of two vectors. | |
template<class vec_t > | |
double | vector_covariance (size_t n, const vec_t &data1, const vec_t &data2) |
Compute the covariance of two vectors. | |
template<class vec_t > | |
double | vector_correlation (size_t n, const vec_t &data1, const vec_t &data2) |
Pearson's correlation. | |
template<class vec_t , class vec2_t > | |
double | vector_pvariance (size_t n1, const vec_t &data1, size_t n2, const vec2_t &data2) |
Pooled variance. | |
template<class vec_t > | |
double | vector_quantile_sorted (size_t n, const vec_t &data, const double f) |
Quantile from sorted data (ascending only) | |
template<class vec_t > | |
double | vector_median_sorted (size_t n, const vec_t &data) |
Return the median of sorted (ascending or descending) data. | |
template<class vec_t , class vec2_t , class vec3_t > | |
double | vector_chi_squared (size_t n, const vec_t &obs, const vec2_t &exp, const vec3_t &err) |
Compute the chi-squared statistic. | |
Weighted vector functions | |
template<class vec_t > | |
double | wvector_mean (size_t n, const vec_t &data, const vec_t &weights) |
Compute the mean of weighted data. | |
template<class vec_t > | |
double | wvector_factor (size_t n, const vec_t &weights) |
Compute a normalization factor for weighted data. | |
template<class vec_t , class vec2_t > | |
double | wvector_variance_fmean (size_t n, const vec_t &data, const vec2_t &weights, double wmean) |
Compute the variance of a weighted vector with a mean known in advance. | |
template<class vec_t , class vec2_t > | |
double | wvector_variance (size_t n, const vec_t &data, const vec2_t &weights, double wmean) |
Compute the variance of a weighted vector with specified mean. | |
template<class vec_t , class vec2_t > | |
double | wvector_variance (size_t n, const vec_t &data, const vec2_t &weights) |
Compute the variance of a weighted vector where mean is computed automatically. | |
template<class vec_t , class vec2_t > | |
double | wvector_stddev_fmean (size_t n, const vec_t &data, const vec2_t &weights, double wmean) |
Compute the standard deviation of a weighted vector with a mean known in advance. | |
template<class vec_t , class vec2_t > | |
double | wvector_stddev (size_t n, const vec_t &data, const vec2_t &weights) |
Compute the standard deviation of a weighted vector where mean is computed automatically. | |
template<class vec_t , class vec2_t > | |
double | wvector_stddev (size_t n, const vec_t &data, const vec2_t &weights, double wmean) |
Compute the standard deviation of a weighted vector with specified mean. | |
template<class vec_t , class vec2_t > | |
double | wvector_sumsq (size_t n, const vec_t &data, const vec2_t &weights, double wmean) |
Compute the weighted sum of squares of data about the specified weighted mean. | |
template<class vec_t , class vec2_t > | |
double | wvector_sumsq (size_t n, const vec_t &data, const vec2_t &weights) |
Compute the weighted sum of squares of data about the weighted mean. | |
template<class vec_t , class vec2_t > | |
double | wvector_absdev (size_t n, const vec_t &data, const vec2_t &weights, double wmean) |
Compute the absolute deviation of data about a specified mean. | |
template<class vec_t , class vec2_t > | |
double | wvector_absdev (size_t n, const vec_t &data, const vec2_t &weights) |
Compute the absolute deviation of data about a specified mean. | |
template<class vec_t , class vec2_t > | |
double | wvector_skew (size_t n, const vec_t &data, const vec2_t &weights, double wmean, double wsd) |
Compute the skewness of data with specified mean and standard deviation. | |
template<class vec_t , class vec2_t > | |
double | wvector_skew (size_t n, const vec_t &data, const vec2_t &weights) |
Compute the skewness of data with specified mean and standard deviation. | |
template<class vec_t , class vec2_t > | |
double | wvector_kurtosis (size_t n, const vec_t &data, const vec2_t &weights, double wmean, double wsd) |
Compute the kurtosis of data with specified mean and standard deviation. | |
template<class vec_t , class vec2_t > | |
double | wvector_kurtosis (size_t n, const vec_t &data, const vec2_t &weights) |
Compute the kurtosis of data with specified mean and standard deviation. |
double vector_mean | ( | size_t | n, |
const vec_t & | data | ||
) |
This function produces the same results as gsl_stats_mean()
.
If n
is zero, this will return zero.
Definition at line 55 of file vec_stats.h.
double vector_variance_fmean | ( | size_t | n, |
const vec_t & | data, | ||
double | mean | ||
) |
This function computes
where the value of is given in
mean
.
This function produces the same results as gsl_stats_variance_with_fixed_mean()
.
Definition at line 76 of file vec_stats.h.
double vector_variance | ( | size_t | n, |
const vec_t & | data, | ||
double | mean | ||
) |
This function computes
where the value of is given in
mean
.
This function produces the same results as gsl_stats_variance_m
.
If n
is 0 or 1, this function will call the error handler.
Definition at line 100 of file vec_stats.h.
double vector_variance | ( | size_t | n, |
const vec_t & | data | ||
) |
This function computes
where is the mean computed with vector_mean().
This function produces the same results as gsl_stats_variance
.
If n
is 0 or 1, this function will call the error handler.
Definition at line 124 of file vec_stats.h.
double vector_stddev_fmean | ( | size_t | n, |
const vec_t & | data, | ||
double | mean | ||
) |
This function computes
where the value of is given in
mean
.
This function produces the same results as gsl_stats_sd_with_fixed_mean()
.
If n
is zero, this function will return zero without calling the error handler.
Definition at line 151 of file vec_stats.h.
double vector_stddev | ( | size_t | n, |
const vec_t & | data | ||
) |
This function computes
where is the mean computed with vector_mean().
This function produces the same results as gsl_stats_sd()
.
If n
is 0 or 1, this function will call the error handler.
Definition at line 170 of file vec_stats.h.
double vector_stddev | ( | size_t | n, |
const vec_t & | data, | ||
double | mean | ||
) |
This function computes
where the value of is given in
mean
.
This function produces the same results as gsl_stats_sd_m()
.
If n
is 0 or 1, this function will call the error handler.
Definition at line 196 of file vec_stats.h.
double vector_absdev | ( | size_t | n, |
const vec_t & | data, | ||
double | mean | ||
) |
This function computes
where the value of is given in
mean
.
This function produces the same results as gsl_stats_absdev_m()
.
If n
is zero, this function will return zero without calling the error handler.
Definition at line 222 of file vec_stats.h.
double vector_absdev | ( | size_t | n, |
const vec_t & | data | ||
) |
This function computes
where the value of is mean as computed from vector_mean().
This function produces the same results as gsl_stats_absdev()
.
If n
is zero, this function will return zero without calling the error handler.
Definition at line 250 of file vec_stats.h.
double vector_skew | ( | size_t | n, |
const vec_t & | data, | ||
double | mean, | ||
double | stddev | ||
) |
This function computes
where the values of and
are given in
mean
and stddev
.
This function produces the same results as gsl_stats_skew_m_sd()
.
If n
is zero, this function will return zero without calling the error handler.
Definition at line 271 of file vec_stats.h.
double vector_skew | ( | size_t | n, |
const vec_t & | data | ||
) |
This function computes
where the values of and
are computed using vector_mean() and vector_stddev().
This function produces the same results as gsl_stats_skew()
.
If n
is zero, this function will return zero without calling the error handler.
Definition at line 297 of file vec_stats.h.
double vector_kurtosis | ( | size_t | n, |
const vec_t & | data, | ||
double | mean, | ||
double | stddev | ||
) |
This function computes
where the values of and
are given in
mean
and stddev
.
This function produces the same results as gsl_stats_kurtosis_m_sd()
.
If n
is zero, this function will return zero without calling the error handler.
Definition at line 320 of file vec_stats.h.
double vector_kurtosis | ( | size_t | n, |
const vec_t & | data | ||
) |
This function computes
where the values of and
are computed using vector_mean() and vector_stddev().
This function produces the same results as gsl_stats_kurtosis()
.
If n
is zero, this function will return zero without calling the error handler.
Definition at line 346 of file vec_stats.h.
double vector_lag1_autocorr | ( | size_t | n, |
const vec_t & | data, | ||
double | mean | ||
) |
This function computes
This function produces the same results as gsl_stats_lag1_autocorrelation_m()
.
If n
is zero, this function will call the error handler.
Definition at line 369 of file vec_stats.h.
double vector_lag1_autocorr | ( | size_t | n, |
const vec_t & | data | ||
) |
This function computes
This function produces the same results as gsl_stats_lag1_autocorrelation()
.
If n
is zero, this function will call the error handler.
Definition at line 404 of file vec_stats.h.
double vector_covariance | ( | size_t | n, |
const vec_t & | data1, | ||
const vec_t & | data2, | ||
double | mean1, | ||
double | mean2 | ||
) |
This function computes
where and
are specified in
mean1
and mean2
, respectively.
This function produces the same results as gsl_stats_covariance_m()
.
If n
is zero, this function will return zero without calling the error handler.
Definition at line 426 of file vec_stats.h.
double vector_covariance | ( | size_t | n, |
const vec_t & | data1, | ||
const vec_t & | data2 | ||
) |
This function computes
where and
are the averages of
data1
and data2
and are computed automatically using vector_mean().
This function produces the same results as gsl_stats_covariance()
.
If n
is zero, this function will return zero without calling the error handler.
Definition at line 455 of file vec_stats.h.
double vector_correlation | ( | size_t | n, |
const vec_t & | data1, | ||
const vec_t & | data2 | ||
) |
This function computes the Pearson correlation coefficient between data1
and data2
.
This function produces the same results as gsl_stats_correlation()
.
If n
is zero, this function will call the error handler.
Definition at line 488 of file vec_stats.h.
double vector_pvariance | ( | size_t | n1, |
const vec_t & | data1, | ||
size_t | n2, | ||
const vec2_t & | data2 | ||
) |
This function produces the same results as gsl_stats_pvariance()
.
If n
is zero, this function will return zero without calling the error handler.
Definition at line 543 of file vec_stats.h.
double vector_quantile_sorted | ( | size_t | n, |
const vec_t & | data, | ||
const double | f | ||
) |
This function returns the quantile f
of data which has already been sorted in ascending order. The quantile, , is found by interpolation using
where and
.
This function produces the same results as gsl_stats_quantile_from_sorted_data()
.
No checks are made to ensure the data is sorted, or to ensure that . If
n
is zero, this function will return zero without calling the error handler.
Definition at line 570 of file vec_stats.h.
double vector_median_sorted | ( | size_t | n, |
const vec_t & | data | ||
) |
This function returns the median of sorted data (either ascending or descending), assuming the data has already been sorted. When the data set has an odd number of elements, the median is the value of the element at index , otherwise, the median is taken to be the average of the elements at indices
and
.
This function produces the same results as gsl_stats_median_from_sorted_data()
.
No checks are made to ensure the data is sorted. If n
is zero, this function will return zero without calling the error handler.
Definition at line 598 of file vec_stats.h.
double vector_chi_squared | ( | size_t | n, |
const vec_t & | obs, | ||
const vec2_t & | exp, | ||
const vec3_t & | err | ||
) |
This function computes
where are the observed values,
are the expected values, and
are the errors.
Definition at line 622 of file vec_stats.h.
double wvector_mean | ( | size_t | n, |
const vec_t & | data, | ||
const vec_t & | weights | ||
) |
This function computes
This function produces the same results as gsl_stats_wmean()
.
Definition at line 650 of file vec_stats.h.
double wvector_factor | ( | size_t | n, |
const vec_t & | weights | ||
) |
This function is used internally in wvector_variance(size_t n, vec_t &data, const vec2_t &weights, double wmean) and wvector_stddev(size_t n, vec_t &data, const vec2_t &weights, double wmean) .
Definition at line 672 of file vec_stats.h.
double wvector_variance_fmean | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights, | ||
double | wmean | ||
) |
This function computes
This function produces the same results as gsl_stats_wvariance_with_fixed_mean()
.
Definition at line 702 of file vec_stats.h.
double wvector_variance | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights, | ||
double | wmean | ||
) |
This function produces the same results as gsl_stats_wvariance_m()
.
Definition at line 725 of file vec_stats.h.
double wvector_variance | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights | ||
) |
This function produces the same results as gsl_stats_wvariance()
.
Definition at line 742 of file vec_stats.h.
double wvector_stddev_fmean | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights, | ||
double | wmean | ||
) |
This function produces the same results as gsl_stats_wsd_with_fixed_mean()
.
Definition at line 756 of file vec_stats.h.
double wvector_stddev | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights | ||
) |
This function produces the same results as gsl_stats_wsd()
.
Definition at line 768 of file vec_stats.h.
double wvector_stddev | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights, | ||
double | wmean | ||
) |
This function produces the same results as gsl_stats_wsd_m()
.
Definition at line 781 of file vec_stats.h.
double wvector_sumsq | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights, | ||
double | wmean | ||
) |
This function produces the same results as gsl_stats_wtss_m()
.
Definition at line 797 of file vec_stats.h.
double wvector_sumsq | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights | ||
) |
This function produces the same results as gsl_stats_wtss()
.
Definition at line 818 of file vec_stats.h.
double wvector_absdev | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights, | ||
double | wmean | ||
) |
This function produces the same results as gsl_stats_wabsdev_m()
.
Definition at line 831 of file vec_stats.h.
double wvector_absdev | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights | ||
) |
This function produces the same results as gsl_stats_wabsdev()
.
Definition at line 852 of file vec_stats.h.
double wvector_skew | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights, | ||
double | wmean, | ||
double | wsd | ||
) |
This function produces the same results as gsl_stats_wskew_m_sd()
.
Definition at line 866 of file vec_stats.h.
double wvector_skew | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights | ||
) |
This function produces the same results as gsl_stats_wskew()
.
Definition at line 888 of file vec_stats.h.
double wvector_kurtosis | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights, | ||
double | wmean, | ||
double | wsd | ||
) |
This function produces the same results as gsl_stats_wkurtosis_m_sd()
.
Definition at line 901 of file vec_stats.h.
double wvector_kurtosis | ( | size_t | n, |
const vec_t & | data, | ||
const vec2_t & | weights | ||
) |
This function produces the same results as gsl_stats_wkurtosis()
.
Definition at line 923 of file vec_stats.h.
Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).