Gromacs  2018.8
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Typedefs | Enumerations | Functions
statistics.h File Reference
#include <cstdio>
#include "gromacs/utility/real.h"
+ Include dependency graph for statistics.h:

Description

Declares simple statistics toolbox.

Authors
David van der Spoel david.nosp@m..van.nosp@m.dersp.nosp@m.oel@.nosp@m.icm.u.nosp@m.u.se

Typedefs

typedef struct gmx_stats * gmx_stats_t
 Abstract container type.
 

Enumerations

enum  {
  estatsOK, estatsNO_POINTS, estatsNO_MEMORY, estatsERROR,
  estatsINVALID_INPUT, estatsNOT_IMPLEMENTED, estatsNR
}
 Error codes returned by the routines.
 
enum  {
  elsqWEIGHT_NONE, elsqWEIGHT_X, elsqWEIGHT_Y, elsqWEIGHT_XY,
  elsqWEIGHT_NR
}
 Enum for statistical weights.
 
enum  { ehistoX, ehistoY, ehistoNR }
 Enum determining which coordinate to histogram.
 

Functions

gmx_stats_t gmx_stats_init ()
 Initiate a data structure. More...
 
void gmx_stats_free (gmx_stats_t stats)
 Destroy a data structure. More...
 
int gmx_stats_remove_outliers (gmx_stats_t stats, double level)
 Remove outliers from a straight line, where level in units of sigma. Level needs to be larger than one obviously. More...
 
int gmx_stats_add_point (gmx_stats_t stats, double x, double y, double dx, double dy)
 Add a point to the data set. More...
 
int gmx_stats_add_points (gmx_stats_t stats, int n, real *x, real *y, real *dx, real *dy)
 Add a series of datapoints at once. The arrays dx and dy may be NULL in that case zero uncertainties will be assumed. More...
 
int gmx_stats_get_point (gmx_stats_t stats, real *x, real *y, real *dx, real *dy, real level)
 Delivers data points from the statistics. More...
 
int gmx_stats_get_ab (gmx_stats_t stats, int weight, real *a, real *b, real *da, real *db, real *chi2, real *Rfit)
 Fit the data to y = ax + b, possibly weighted, if uncertainties have been input. da and db may be NULL. More...
 
int gmx_stats_get_a (gmx_stats_t stats, int weight, real *a, real *da, real *chi2, real *Rfit)
 Fit the data to y = ax, possibly weighted, if uncertainties have have been input. da and db may be NULL. More...
 
int gmx_stats_get_corr_coeff (gmx_stats_t stats, real *R)
 Get the correlation coefficient. More...
 
int gmx_stats_get_rmsd (gmx_stats_t stats, real *rmsd)
 Get the root mean square deviation. More...
 
int gmx_stats_get_npoints (gmx_stats_t stats, int *N)
 Get the number of points. More...
 
int gmx_stats_get_average (gmx_stats_t stats, real *aver)
 Computes and returns the average value. More...
 
int gmx_stats_get_sigma (gmx_stats_t stats, real *sigma)
 Computes and returns the standard deviation. More...
 
int gmx_stats_get_error (gmx_stats_t stats, real *error)
 Computes and returns the standard error. More...
 
int gmx_stats_get_ase (gmx_stats_t stats, real *aver, real *sigma, real *error)
 Pointers may be null, in which case no assignment will be done. More...
 
int gmx_stats_dump_xy (gmx_stats_t stats, FILE *fp)
 Dump the x, y, dx, dy data to a text file. More...
 
int gmx_stats_make_histogram (gmx_stats_t stats, real binwidth, int *nbins, int ehisto, int normalized, real **x, real **y)
 Make a histogram of the data present. More...
 
const char * gmx_stats_message (int estats)
 Return message belonging to error code. More...
 
int lsq_y_ax (int n, real x[], real y[], real *a)
 Fit a straight line y=ax thru the n data points x, y, return the slope in *a. More...
 
int lsq_y_ax_b (int n, real x[], real y[], real *a, real *b, real *r, real *chi2)
 Fit a straight line y=ax+b thru the n data points x, y. More...
 
int lsq_y_ax_b_xdouble (int n, double x[], real y[], real *a, real *b, real *r, real *chi2)
 Fit a straight line y=ax+b thru the n data points x, y. More...
 
int lsq_y_ax_b_error (int n, real x[], real y[], real dy[], real *a, real *b, real *da, real *db, real *r, real *chi2)
 Fit a straight line y=ax+b thru the n data points x, y. More...
 

Function Documentation

int gmx_stats_add_point ( gmx_stats_t  stats,
double  x,
double  y,
double  dx,
double  dy 
)

Add a point to the data set.

Parameters
[in]statsThe data structure
[in]xThe x value
[in]yThe y value
[in]dxThe error in the x value
[in]dyThe error in the y value
Returns
error code
int gmx_stats_add_points ( gmx_stats_t  stats,
int  n,
real x,
real y,
real dx,
real dy 
)

Add a series of datapoints at once. The arrays dx and dy may be NULL in that case zero uncertainties will be assumed.

Parameters
[in]statsThe data structure
[in]nNumber of points
[in]xThe array of x values
[in]yThe array of y values
[in]dxThe error in the x value
[in]dyThe error in the y value
Returns
error code
int gmx_stats_dump_xy ( gmx_stats_t  stats,
FILE *  fp 
)

Dump the x, y, dx, dy data to a text file.

Parameters
[in]statsThe data structure
[in]fpFile pointer
Returns
error code
void gmx_stats_free ( gmx_stats_t  stats)

Destroy a data structure.

Parameters
statsThe data structure
int gmx_stats_get_a ( gmx_stats_t  stats,
int  weight,
real a,
real da,
real chi2,
real Rfit 
)

Fit the data to y = ax, possibly weighted, if uncertainties have have been input. da and db may be NULL.

Parameters
[in]statsThe data structure
[in]weighttype of weighting
[out]aslope
[out]dasigma in a
[out]chi2normalized quality of fit
[out]Rfitcorrelation coefficient
Returns
error code
int gmx_stats_get_ab ( gmx_stats_t  stats,
int  weight,
real a,
real b,
real da,
real db,
real chi2,
real Rfit 
)

Fit the data to y = ax + b, possibly weighted, if uncertainties have been input. da and db may be NULL.

Parameters
[in]statsThe data structure
[in]weighttype of weighting
[out]aslope
[out]bintercept
[out]dasigma in a
[out]dbsigma in b
[out]chi2normalized quality of fit
[out]Rfitcorrelation coefficient
Returns
error code
int gmx_stats_get_ase ( gmx_stats_t  stats,
real aver,
real sigma,
real error 
)

Pointers may be null, in which case no assignment will be done.

Parameters
[in]statsThe data structure
[out]averAverage value
[out]sigmaStandard deviation
[out]errorStandard error
Returns
error code
int gmx_stats_get_average ( gmx_stats_t  stats,
real aver 
)

Computes and returns the average value.

Parameters
[in]statsThe data structure
[out]averAverage value
Returns
error code
int gmx_stats_get_corr_coeff ( gmx_stats_t  stats,
real R 
)

Get the correlation coefficient.

Parameters
[in]statsThe data structure
[out]Rthe correlation coefficient between the data (x and y) as input to the structure.
Returns
error code
int gmx_stats_get_error ( gmx_stats_t  stats,
real error 
)

Computes and returns the standard error.

Parameters
[in]statsThe data structure
[out]errorStandard error
Returns
error code
int gmx_stats_get_npoints ( gmx_stats_t  stats,
int *  N 
)

Get the number of points.

Parameters
[in]statsThe data structure
[out]Nnumber of data points
Returns
error code
int gmx_stats_get_point ( gmx_stats_t  stats,
real x,
real y,
real dx,
real dy,
real  level 
)

Delivers data points from the statistics.

Should be used in a while loop. Variables for either pointer may be NULL, in which case the routine can be used as an expensive point counter. Return the data points one by one. Return estatsOK while there are more points, and returns estatsNOPOINTS when the last point has been returned. If level > 0 then the outliers outside level*sigma are reported only.

Parameters
[in]statsThe data structure
[out]xThe array of x values
[out]yThe array of y values
[out]dxThe error in the x value
[out]dyThe error in the y value
[in]levelsigma level (see above)
Returns
error code
int gmx_stats_get_rmsd ( gmx_stats_t  stats,
real rmsd 
)

Get the root mean square deviation.

Parameters
[in]statsThe data structure
[out]rmsdthe root mean square deviation between x and y values.
Returns
error code
int gmx_stats_get_sigma ( gmx_stats_t  stats,
real sigma 
)

Computes and returns the standard deviation.

Parameters
[in]statsThe data structure
[out]sigmaStandard deviation
Returns
error code
gmx_stats_t gmx_stats_init ( )

Initiate a data structure.

Returns
the data structure
int gmx_stats_make_histogram ( gmx_stats_t  stats,
real  binwidth,
int *  nbins,
int  ehisto,
int  normalized,
real **  x,
real **  y 
)

Make a histogram of the data present.

Uses either binwidth to determine the number of bins, or nbins to determine the binwidth, therefore one of these should be zero, but not the other. If *nbins = 0 the number of bins will be returned in this variable. ehisto should be one of ehistoX or ehistoY. If normalized not equal to zero, the integral of the histogram will be normalized to one. The output is in two arrays, *x and *y, to which you should pass a pointer. Memory for the arrays will be allocated as needed. Function returns one of the estats codes.

Parameters
[in]statsThe data structure
[in]binwidthFor the histogram
[in]nbinsNumber of bins
[in]ehistoType (see enum above)
[in]normalizedsee above
[out]xsee above
[out]ysee above
Returns
error code
const char* gmx_stats_message ( int  estats)

Return message belonging to error code.

Parameters
[in]estatserror code
int gmx_stats_remove_outliers ( gmx_stats_t  stats,
double  level 
)

Remove outliers from a straight line, where level in units of sigma. Level needs to be larger than one obviously.

Parameters
[in]statsThe data structure
[in]levelThe sigma level
Returns
error code
int lsq_y_ax ( int  n,
real  x[],
real  y[],
real a 
)

Fit a straight line y=ax thru the n data points x, y, return the slope in *a.

Parameters
[in]nnumber of points
[in]xdata points x
[in]ydata point y
[out]aslope
Returns
error code
int lsq_y_ax_b ( int  n,
real  x[],
real  y[],
real a,
real b,
real r,
real chi2 
)

Fit a straight line y=ax+b thru the n data points x, y.

Parameters
[in]nnumber of points
[in]xdata points x
[in]ydata point y
[out]aslope
[out]bintercept
[out]rcorrelation coefficient
[out]chi2quality of fit
Returns
error code
int lsq_y_ax_b_error ( int  n,
real  x[],
real  y[],
real  dy[],
real a,
real b,
real da,
real db,
real r,
real chi2 
)

Fit a straight line y=ax+b thru the n data points x, y.

Parameters
[in]nnumber of points
[in]xdata points x
[in]ydata point y
[in]dyuncertainty in data point y
[out]aslope
[out]bintercept
[out]daerror in slope
[out]dberror in intercept
[out]rcorrelation coefficient
[out]chi2quality of fit
Returns
error code
int lsq_y_ax_b_xdouble ( int  n,
double  x[],
real  y[],
real a,
real b,
real r,
real chi2 
)

Fit a straight line y=ax+b thru the n data points x, y.

Parameters
[in]nnumber of points
[in]xdata points x
[in]ydata point y
[out]aslope
[out]bintercept
[out]rcorrelation coefficient
[out]chi2quality of fit
Returns
error code