BAT  0.9.4
The Bayesian analysis toolkit
 All Classes Namespaces Functions Variables Enumerations
BCMath Namespace Reference

Some useful mathematic functions. More...

Functions

double LogGaus (double x, double mean, double sigma, bool norm)
 
double LogPoisson (double x, double par)
 
double ApproxBinomial (int n, int k, double p)
 
double LogApproxBinomial (int n, int k, double p)
 
double LogBinomFactor (int n, int k)
 
double ApproxLogFact (double x)
 
double LogFact (int n)
 
void CacheFactorial (unsigned int n)
 
double LogNoverK (int n, int k)
 
int Nint (double x)
 
double LogBreitWignerNonRel (double x, double mean, double Gamma, bool norm)
 
double LogBreitWignerRel (double x, double mean, double Gamma)
 
double LogChi2 (double x, int n)
 
double LogVoigtian (double x, double sigma, double gamma)
 
double LogGammaPDF (double x, double alpha, double beta)
 
double LogLogNormal (double x, double mean, double sigma)
 
double SplitGaussian (double *x, double *par)
 
double chi2 (double *x, double *par)
 
void RandomChi2 (std::vector< double > &randoms, int K)
 
TH1D * ECDF (const std::vector< double > &data)
 
std::vector< int > longestRuns (const std::vector< bool > &bitStream)
 
std::vector< double > longestRunsChi2 (const std::vector< double > &yMeasured, const std::vector< double > &yExpected, const std::vector< double > &sigma)
 
double longestRunFrequency (unsigned longestObserved, unsigned int nTrials)
 
double Rvalue (const std::vector< double > &chain_means, const std::vector< double > &chain_variances, const unsigned &chain_length, const bool &strict) throw (std::invalid_argument, std::domain_error)
 
double longestRunFrequency (unsigned int longestObserved, unsigned int nTrials)
 
p value methods
double CorrectPValue (const double &pvalue, const unsigned &npar, const unsigned &nobservations) throw (std::domain_error)
 
double FastPValue (const std::vector< unsigned > &observed, const std::vector< double > &expected, unsigned nIterations, unsigned seed) throw (std::invalid_argument)
 

Detailed Description

Some useful mathematic functions.

Author
Daniel Kollar
Kevin Kröninger
Jing Liu
Frederik Beaujean
Version
1.0
Date
08.2008 A namespace which encapsulates some mathematical functions necessary for BAT.

Function Documentation

double BCMath::ApproxBinomial ( int  n,
int  k,
double  p 
)

Calculates Binomial probability using approximations for factorial calculations if calculation for number greater than 20 required using the BCMath::ApproxLogFact function.

Definition at line 73 of file BCMath.cxx.

double BCMath::ApproxLogFact ( double  x)

Calculates natural logarithm of the n-factorial (n!) using Srinivasa Ramanujan approximation log(n!) = n*log(n) - n + log(n*(1.+4.*n*(1.+2.*n)))/6. + log(PI)/2. if n > 20. If n <= 20 it uses BCMath::LogFact to calculate it exactly.

Definition at line 125 of file BCMath.cxx.

void BCMath::CacheFactorial ( unsigned int  n)

Cache factorials for first

  • n integers. The cache is filled upon first call of LogFact().

Definition at line 169 of file BCMath.cxx.

double BCMath::CorrectPValue ( const double &  pvalue,
const unsigned &  npar,
const unsigned &  nobservations 
)
throw (std::domain_error
)

Correct a p value by transforming to a chi^2 with dof=nobservations, then back to a pvalue with dof reduced by number of fit parameters.

Parameters
pvalueThe p value to correct.
nparThe number of fit parameters.
nobservationsThe number of data points.
Returns
corrected p value

Definition at line 665 of file BCMath.cxx.

TH1D * BCMath::ECDF ( const std::vector< double > &  data)

Calculate the empirical cumulative distribution function for one dimensional data vector. For consistency, the ECDF of value smaller than the minimum observed (underflow bin) is zero, and for larger than maximum (overflow bin) it is one.

Parameters
datathe observations
Returns
histogram with normalized ECDF

Definition at line 333 of file BCMath.cxx.

double BCMath::FastPValue ( const std::vector< unsigned > &  observed,
const std::vector< double > &  expected,
unsigned  nIterations = 1e5,
unsigned  seed = 0 
)
throw (std::invalid_argument
)

Calculate the p value using fast MCMC for a histogram and the likelihood as test statistic. The method is explained in the appendix of http://arxiv.org/abs/1011.1674

Parameters
observedThe counts of observed events
expectedThe expected number of events (Poisson means)
nIterationsControls number of pseudo data sets
seedSet to nonzero value for reproducible results.
Returns
The p value

Definition at line 690 of file BCMath.cxx.

double BCMath::LogApproxBinomial ( int  n,
int  k,
double  p 
)

Calculates natural logarithm of the Binomial probability using approximations for factorial calculations if calculation for number greater than 20 required using the BCMath::ApproxLogFact function.

Definition at line 80 of file BCMath.cxx.

double BCMath::LogBinomFactor ( int  n,
int  k 
)

Calculates natural logarithm of the Binomial factor "n over k" using approximations for factorial calculations if calculation for number greater than 20 required using the BCMath::ApproxLogFact function. Even for large numbers the calculation is performed precisely, if n-k < 5

Definition at line 101 of file BCMath.cxx.

double BCMath::LogBreitWignerNonRel ( double  x,
double  mean,
double  Gamma,
bool  norm = false 
)

Calculates the logarithm of the non-relativistic Breit-Wigner distribution.

Definition at line 225 of file BCMath.cxx.

double BCMath::LogChi2 ( double  x,
int  n 
)

Calculates the logarithm of chi square function: chi2(double x; size_t n)

Definition at line 244 of file BCMath.cxx.

double BCMath::LogFact ( int  n)

Calculates natural logarithm of the n-factorial (n!)

Definition at line 135 of file BCMath.cxx.

double BCMath::LogGammaPDF ( double  x,
double  alpha,
double  beta 
)

Returns the log of the Gamma PDF.

Definition at line 273 of file BCMath.cxx.

double BCMath::LogGaus ( double  x,
double  mean = 0,
double  sigma = 1,
bool  norm = false 
)

Calculate the natural logarithm of a gaussian function with mean and sigma. If norm=true (default is false) the result is multiplied by the normalization constant, i.e. divided by sqrt(2*Pi)*sigma.

Definition at line 33 of file BCMath.cxx.

double BCMath::LogLogNormal ( double  x,
double  mean = 0,
double  sigma = 1 
)

Return the log of the log normal distribution

Definition at line 279 of file BCMath.cxx.

double BCMath::LogNoverK ( int  n,
int  k 
)

Calculates natural logarithm of the Binomial factor "n over k".

Definition at line 176 of file BCMath.cxx.

double BCMath::LogPoisson ( double  x,
double  par 
)

Calculate the natural logarithm of a poisson distribution.

Definition at line 57 of file BCMath.cxx.

double BCMath::LogVoigtian ( double  x,
double  sigma,
double  gamma 
)

Calculates the logarithm of normalized voigtian function: voigtian(double x, double sigma, double gamma)

voigtian is a convolution of the following two functions: gaussian(x) = 1/(sqrt(2*pi)*sigma) * exp(x*x/(2*sigma*sigma) and lorentz(x) = (1/pi)*(gamma/2) / (x*x + (gamma/2)*(gamma/2))

it is singly peaked at x=0. The width of the peak is decided by sigma and gamma, so they should be positive.

Definition at line 262 of file BCMath.cxx.

double BCMath::longestRunFrequency ( unsigned int  longestObserved,
unsigned int  nTrials 
)

Find the sampling probability that, given n independent Bernoulli trials with success rate = failure rate = 1/2, the longest run of consecutive successes is greater than the longest observed run. Key idea from Burr, E.J. & Cane, G. Longest Run of Consecutive Observations Having a Specified Attribute. Biometrika 48, 461-465 (1961).

Parameters
longestObservedactual longest run
nTrialsnumber of independent trials
Returns
frequency
std::vector< int > BCMath::longestRuns ( const std::vector< bool > &  bitStream)

Find the longest runs of zeros and ones in the bit stream

Parameters
bitStreaminput sequence of boolean values
Returns
runs first entry the longest zeros run, second entry the longest ones run

Definition at line 388 of file BCMath.cxx.

std::vector< double > BCMath::longestRunsChi2 ( const std::vector< double > &  yMeasured,
const std::vector< double > &  yExpected,
const std::vector< double > &  sigma 
)

Find the longest success/failure run in set of norm. distributed variables. Success = observation >= expectation. Runs are weighted by the total chi^2 of all elements in the run

Parameters
yMeasuredthe observations
yExpectedthe expected values
sigmathe theoretical uncertainties on the expectations
Returns
runs first entry the max. weight failure run, second entry the max. success run

Definition at line 441 of file BCMath.cxx.

int BCMath::Nint ( double  x)

Returns the nearest integer of a double number.

Definition at line 204 of file BCMath.cxx.

void BCMath::RandomChi2 ( std::vector< double > &  randoms,
int  K 
)

Get N random numbers distributed according to chi square function with K degrees of freedom

Definition at line 318 of file BCMath.cxx.

double BCMath::Rvalue ( const std::vector< double > &  chain_means,
const std::vector< double > &  chain_variances,
const unsigned &  chain_length,
const bool &  strict = true 
)
throw (std::invalid_argument,
std::domain_error
)

Compute the R-value according to Gelman-Rubin, GR1992 : Gelman, A. and Rubin, D.B. Inference from Iterative Simulation Using Multiple Sequences Statistical Science, Vol. 7, No. 4 (Nov. 1992), pp. 457-472

Parameters
chain_means
chain_variances
chain_length
strictIf true, use the algorithm laid forth in the paper, else use a relaxed version which generally leads to smaller R-values.
Returns
R-value

Definition at line 551 of file BCMath.cxx.

double BCMath::SplitGaussian ( double *  x,
double *  par 
)

Gaussian with different uncertainties to either side of the mode. Interface to use with TF1.

Parameters
xRandom variable
parMean, sigmadown, sigmaup
Returns
density f(x)

Definition at line 296 of file BCMath.cxx.