Tomographer  v5.4
Tomographer C++ Framework Documentation
Tomographer::BinningAnalysis< Params, LoggerType_ > Class Template Reference

Simple binning analysis for determining error bars. More...

#include <tomographer/mhrw_bin_err.h>

+ Inheritance diagram for Tomographer::BinningAnalysis< Params, LoggerType_ >:
+ Collaboration diagram for Tomographer::BinningAnalysis< Params, LoggerType_ >:

Public Types

enum  
 Constants for error bar convergence analysis. More...
 
typedef Params::ValueType ValueType
 Type of the value(s) for which we are calculating error bars. See BinningAnalysisParams::ValueType.
 
typedef Params::CountIntType CountIntType
 Type used to count the number of samples. See BinningAnalysisParams::CountIntType.
 
typedef Params::SamplesArray SamplesArray
 Type of the internal buffer for the raw samples. See BinningAnalysisParams::SamplesArray.
 
typedef Params::BinSumArray BinSumArray
 Type used to store the sum of values. See BinningAnalysisParams::BinSumArray.
 
typedef Params::BinSumSqArray BinSumSqArray
 Type used to store the sum of squares of values at each binning level. See BinningAnalysisParams::BinSumSqArray.
 
typedef LoggerType_ LoggerType
 Type of the logger we will be logging debugging messages to.
 
- Public Types inherited from Tomographer::Tools::EigenAlignedOperatorNewProvider
typedef EigenAlignedOperatorNewProvider OperatorNewProviderType
 

Public Member Functions

 BinningAnalysis (int num_track_values_, int num_levels_, LoggerType &logger_)
 Constructor. More...
 
void reset ()
 Reset this object and start again as if freshly constructed. More...
 
template<typename Derived >
void processNewValues (const Eigen::DenseBase< Derived > &vals)
 Process new raw samples. More...
 
template<typename CalcValType , ENABLED_IF( NumTrackValuesCTime==1) >
void processNewValue (const CalcValType val)
 Process a new value (if we're tracking a single function only) More...
 
Eigen::Index effectiveSampleSize () const
 Number of samples required to flush all levels of the binning analysis. More...
 
CountIntType getNumFlushes () const
 Return the number of times the collected samples were flushed. More...
 
template<ENABLED_IF( StoreBinSums) >
auto getBinMeans () const
 Get the average of each tracked value observed. More...
 
auto getBinSqmeans () const
 Get the raw average of the squared values observed, for each binning level. More...
 
template<ENABLED_IF( StoreBinSums) >
const BinSumArraygetBinSum () const
 Get the sum of each tracked value observed. More...
 
const BinSumSqArraygetBinSumsq () const
 Get the raw sums of the squared values observed, at each binning level. More...
 
template<typename Derived >
BinSumSqArray calcErrorLevels (const Eigen::ArrayBase< Derived > &means) const
 Calculate the error bars of samples at different binning levels. More...
 
template<typename Derived >
BinSumArray calcErrorLastLevel (const Eigen::ArrayBase< Derived > &means) const
 Calculate the error bar of samples (from the last binning level). More...
 
template<ENABLED_IF( StoreBinSums) >
BinSumSqArray calcErrorLevels () const
 Calculate the error bars of samples at different binning levels. More...
 
template<ENABLED_IF( StoreBinSums) >
BinSumArray calcErrorLastLevel () const
 Calculate the error bar of samples (from the last binning level). More...
 
Eigen::ArrayXi determineErrorConvergence (const Eigen::Ref< const BinSumSqArray > &error_levels) const
 Attempt to determine if the error bars have converged. More...
 

Public Attributes

const Tools::StaticOrDynamic< int,(NumTrackValuesCTime==Eigen::Dynamic), NumTrackValuesCTimenumTrackValues
 The number of functions being tracked/analyzed. More...
 
const Tools::StaticOrDynamic< int,(NumLevelsCTime==Eigen::Dynamic), NumLevelsCTimenumLevels
 The number of levels in the binning analysis. More...
 
const Tools::StaticOrDynamic< CountIntType,(SamplesSizeCTime==Eigen::Dynamic),(SamplesSizeCTime< 0 ? 0 :SamplesSizeCTime)> samplesSize
 The size of our samples buffer. (avoid using, might change in the future) More...
 

Static Public Attributes

static constexpr int NumTrackValuesCTime = Params::NumTrackValuesCTime
 Number of values we are tracking/analyzing, if known at compile-time or Eigen::Dynamic. See BinningAnalysisParams::NumTrackValuesCTime.
 
static constexpr int NumLevelsCTime = Params::NumLevelsCTime
 Number of binning levels in our binning analysis, if known at compile-time or Eigen::Dynamic. See BinningAnalysisParams::NumLevelsCTime.
 
static constexpr int NumLevelsPlusOneCTime = Params::NumLevelsPlusOneCTime
 Number of binning levels in our binning analysis plus one, if known at compile-time or Eigen::Dynamic. See BinningAnalysisParams::NumLevelsPlusOneCTime.
 
static constexpr int SamplesSizeCTime = Params::SamplesSizeCTime
 Size of the buffer which holds the raw sequence, if of fixed size and known at compile-time. See BinningAnalysisParams::SamplesSizeCTime.
 
static constexpr bool StoreBinSums = Params::StoreBinSums
 Whether we should store the bin sums, or whether they can be provided independently by the user. See BinningAnalysisParams::StoreBinSums.
 

Detailed Description

template<typename Params, typename LoggerType_>
class Tomographer::BinningAnalysis< Params, LoggerType_ >

Simple binning analysis for determining error bars.

This class implements the theoretical procedure described in the theory page Binning Analysis.

This class in fact can perform binning analysis in parallel on several different functions (or values you're integrating). This is useful, for example, to determine error bars and run the binning analysis independently on each bin of a histogram (see ValueHistogramWithBinningMHRWStatsCollector).

Several values here may be either specified at compile-time or at run-time, inspired by Eigen's mechanism using Eigen::Dynamic. For example, the number of independent values to track may be fixed at compile-time in the NumTrackValues template parameter of BinningAnalysisParams which you specify as the Params parameter here; however if the number of values will be determined dynamically at run-time, the NumTrackValues template parameter should be set to Eigen::Dynamic.

Raw samples are added to the analysis by calling processNewValues() with a list of values, one new value per independent function integration that we are tracking. You should call this function repeatedly, once per sample. This class then takes care of adding those raw samples to an internal samples buffer and flushing this buffer whenever is needed.

The final results of the binning analysis can be obtained by calling calcErrorLevels() to get the error bars at each binning level. Note that the analysis will only include all the flushed samples, i.e. all samples up to a multiple of samplesSize(). If additional samples were given which didn't fill the samples buffer, they are ignored for the error analysis but they are included in getBinMeans() and getBinSum().

To calculate the error levels, this class also needs to keep track of the sum of all samples to calculate their mean. However, this is often done independently, e.g. with a ValueHistogramMHRWStatsCollector. In order to avoid calculating this twice, you may tell this BinningAnalysis class NOT to worry about calculating the means of the samples and that you will provide the sample means yourself. This is specified using the StoreBinSums template parameter to BinningAnalysisParams. In this case, you'll need to use the variants of calcErrorLevels() and calcErrorLastLevel() with one argument, where you specify the bin means yourself.

You may also access some other information such as the sum of the squares of the samples at different binning levels (see getBinSqmeans()).

Template Parameters
Paramsis a BinningAnalysisParams with template parameters. This type specifies such types and values as the type of the values we're integrating, the number of independent functions we're integrating, and so on.
LoggerTypea logger type as usual, for us to provide information about what we're doing and how the values look like.
Todo:
Further optimizations. For now, the raw samples are stored until they can be flushed up to the last level. Change this to do what ALPS' implementation does: flushing whenever it can, using a much shorter memory buffer, see https://alps.comp-phys.org/svn/alps1/trunk/alps/src/alps/alea/simplebinning.h .

Definition at line 361 of file mhrw_bin_err.h.

Member Enumeration Documentation

§ anonymous enum

template<typename Params, typename LoggerType_>
anonymous enum

Constants for error bar convergence analysis.

Enumerator
UNKNOWN_CONVERGENCE 

Unable to determine whether the error bars have converged.

CONVERGED 

The error bars appear to have converged.

NOT_CONVERGED 

The error bars don't seem to have converged.

Definition at line 442 of file mhrw_bin_err.h.

Constructor & Destructor Documentation

§ BinningAnalysis()

template<typename Params, typename LoggerType_>
Tomographer::BinningAnalysis< Params, LoggerType_ >::BinningAnalysis ( int  num_track_values_,
int  num_levels_,
LoggerType logger_ 
)
inline

Constructor.

You must specify also the number of values that we will be tracking independently (numTrackValues) and the number of binning levels (numLevels). If compile-time values have been provided as template parameters not equal to Eigen::Dynamic, these values MUST be equal to their compile-time given counterpart values.

Specify also a reference to a logger for logging messages. See Logging and Loggers.

Definition at line 520 of file mhrw_bin_err.h.

Member Function Documentation

§ calcErrorLastLevel() [1/2]

template<typename Params, typename LoggerType_>
template<typename Derived >
BinSumArray Tomographer::BinningAnalysis< Params, LoggerType_ >::calcErrorLastLevel ( const Eigen::ArrayBase< Derived > &  means) const
inline

Calculate the error bar of samples (from the last binning level).

Return a vector of numTrackValues elements, where the i -th item corresponds to the error bar of the i -th value determined from binning level numLevels.

If the error bars converged (see determineErrorConvergence()), this should be a good estimate of the error bars on the corresponding values.

Use this variant of the function if this class doesn't store the bin means. If so, you need to provide the value of the means explicitly to the parameter means.

Definition at line 798 of file mhrw_bin_err.h.

§ calcErrorLastLevel() [2/2]

template<typename Params, typename LoggerType_>
template<ENABLED_IF( StoreBinSums) >
BinSumArray Tomographer::BinningAnalysis< Params, LoggerType_ >::calcErrorLastLevel ( ) const
inline

Calculate the error bar of samples (from the last binning level).

Return a vector of numTrackValues elements, where the i -th item corresponds to the error bar of the i -th value determined from binning level numLevels.

If the error bars converged (see determineErrorConvergence()), this should be a good estimate of the error bars on the corresponding values.

Use this variant of the function if this class stores the bin means. If this is not the case, you will need to call the variant calcErrorLastLevel(const Eigen::ArrayBase<Derived> & means) with the values of the means.

Definition at line 837 of file mhrw_bin_err.h.

§ calcErrorLevels() [1/2]

template<typename Params, typename LoggerType_>
template<typename Derived >
BinSumSqArray Tomographer::BinningAnalysis< Params, LoggerType_ >::calcErrorLevels ( const Eigen::ArrayBase< Derived > &  means) const
inline

Calculate the error bars of samples at different binning levels.

Return an array of shape (numTrackValues, numLevels) where element (i,k) corresponds to the error of the i -th value at binning level k. Binning level k=0 corresponds to the naive error bar from of the samples (no binning).

Use this variant of the function if this class doesn't store the bin means. If so, you need to provide the value of the means explicitly to the parameter means.

Todo:
this should be optimizable, using directly bin_sumsq and not effectively repeating the powersOfTwo constants...

Definition at line 753 of file mhrw_bin_err.h.

§ calcErrorLevels() [2/2]

template<typename Params, typename LoggerType_>
template<ENABLED_IF( StoreBinSums) >
BinSumSqArray Tomographer::BinningAnalysis< Params, LoggerType_ >::calcErrorLevels ( ) const
inline

Calculate the error bars of samples at different binning levels.

Return an array of shape (numTrackValues, numLevels) where element (i,k) corresponds to the error of the i -th value at binning level k. Binning level k=0 corresponds to the naive error bar from of the samples (no binning).

Use this variant of the function if this class stores the bin means. If this is not the case, you will need to call the variant calcErrorLevels(const Eigen::ArrayBase<Derived> & means) with the values of the means.

Definition at line 819 of file mhrw_bin_err.h.

§ determineErrorConvergence()

template<typename Params, typename LoggerType_>
Eigen::ArrayXi Tomographer::BinningAnalysis< Params, LoggerType_ >::determineErrorConvergence ( const Eigen::Ref< const BinSumSqArray > &  error_levels) const
inline

Attempt to determine if the error bars have converged.

Call this method after calculating the error bars for each level with calcErrorLevels(). Use the return value of that function to feed in the input here.

Returns
an array of integers, of length numTrackValues, each set to one of CONVERGED, NOT_CONVERGED or UNKNOWN_CONVERGENCE.

Contains code inspired by ALPS project, see https://alps.comp-phys.org/svn/alps1/trunk/alps/src/alps/alea/simplebinning.h.

Definition at line 854 of file mhrw_bin_err.h.

§ effectiveSampleSize()

template<typename Params, typename LoggerType_>
Eigen::Index Tomographer::BinningAnalysis< Params, LoggerType_ >::effectiveSampleSize ( ) const
inline

Number of samples required to flush all levels of the binning analysis.

Make sure the random walk gathers a number of samples which is a multiple of effectiveSampleSize(), to make sure that the binning analysis has taken into account each and every sample. Otherwise, the binning analysis could ignore the last few samples, and, for instance, you could have an error bar exactly equal to zero (if no samples up to then were observed), but a bin count with one sample. This really bad for a fit, for instance, where the error bar is used to calculate the relevant \( \chi^2 \) quantity.

Since
Added in Tomographer 5.3.

Definition at line 647 of file mhrw_bin_err.h.

§ getBinMeans()

template<typename Params, typename LoggerType_>
template<ENABLED_IF( StoreBinSums) >
auto Tomographer::BinningAnalysis< Params, LoggerType_ >::getBinMeans ( ) const
inline

Get the average of each tracked value observed.

Definition at line 689 of file mhrw_bin_err.h.

§ getBinSqmeans()

template<typename Params, typename LoggerType_>
auto Tomographer::BinningAnalysis< Params, LoggerType_ >::getBinSqmeans ( ) const
inline

Get the raw average of the squared values observed, for each binning level.

The vector bin_sqmeans.col(0) contains the raw average of the squares of the raw values observed, bin_sqmeans.col(1) the raw average of the squares of the values averaged 2 by 2 (i.e. at the first binning level), and so on.

Definition at line 703 of file mhrw_bin_err.h.

§ getBinSum()

template<typename Params, typename LoggerType_>
template<ENABLED_IF( StoreBinSums) >
const BinSumArray& Tomographer::BinningAnalysis< Params, LoggerType_ >::getBinSum ( ) const
inline

Get the sum of each tracked value observed.

This is only available if the StoreBinSums template parameter was set to true.

Definition at line 729 of file mhrw_bin_err.h.

§ getBinSumsq()

template<typename Params, typename LoggerType_>
const BinSumSqArray& Tomographer::BinningAnalysis< Params, LoggerType_ >::getBinSumsq ( ) const
inline

Get the raw sums of the squared values observed, at each binning level.

The vector bin_sumsq.col(0) contains the raw sum of the squares of the raw values observed, bin_sumsq.col(1) the raw sum of the squares of the values averaged 2 by 2 (i.e. at the first binning level), and so on.

Definition at line 737 of file mhrw_bin_err.h.

§ getNumFlushes()

template<typename Params, typename LoggerType_>
CountIntType Tomographer::BinningAnalysis< Params, LoggerType_ >::getNumFlushes ( ) const
inline

Return the number of times the collected samples were flushed.

This corresponds to the number of values of which the most coarse-grained binned averaging consists of.

Definition at line 683 of file mhrw_bin_err.h.

§ processNewValue()

template<typename Params, typename LoggerType_>
template<typename CalcValType , ENABLED_IF( NumTrackValuesCTime==1) >
void Tomographer::BinningAnalysis< Params, LoggerType_ >::processNewValue ( const CalcValType  val)
inline

Process a new value (if we're tracking a single function only)

Use this variant of the function if we're tracking a single function only, so that you don't have to specify a 1-element "array" to processNewValues().

Definition at line 626 of file mhrw_bin_err.h.

§ processNewValues()

template<typename Params, typename LoggerType_>
template<typename Derived >
void Tomographer::BinningAnalysis< Params, LoggerType_ >::processNewValues ( const Eigen::DenseBase< Derived > &  vals)
inline

Process new raw samples.

Call this function whenever you have a new sample with corresponding values. The argument should evaluate to a column vector; each element of the vector corresponds to the value of a function we're tracking. The length of the vector must equal numTrackValues().

This will add the samples to the internal raw sample buffer, and flush the buffer as required.

Definition at line 562 of file mhrw_bin_err.h.

§ reset()

template<typename Params, typename LoggerType_>
void Tomographer::BinningAnalysis< Params, LoggerType_ >::reset ( )
inline

Reset this object and start again as if freshly constructed.

After calling reset(), you may use this object as if you had freshly constructed it. All values are reset to zero and the samples buffer is emptied.

Definition at line 542 of file mhrw_bin_err.h.

Member Data Documentation

§ numLevels

template<typename Params, typename LoggerType_>
const Tools::StaticOrDynamic<int, (NumLevelsCTime==Eigen::Dynamic), NumLevelsCTime> Tomographer::BinningAnalysis< Params, LoggerType_ >::numLevels

The number of levels in the binning analysis.

The number may be obtained by calling binning_analysis.numLevels(). This will either return the compile-time fixed value NumTrackValuesCTime, or the value which was set dynamically at runtime.

Eigen::Dynamic is never returned. See Tools::StaticOrDynamic.

Definition at line 428 of file mhrw_bin_err.h.

§ numTrackValues

template<typename Params, typename LoggerType_>
const Tools::StaticOrDynamic<int, (NumTrackValuesCTime==Eigen::Dynamic), NumTrackValuesCTime> Tomographer::BinningAnalysis< Params, LoggerType_ >::numTrackValues

The number of functions being tracked/analyzed.

The number may be obtained by calling binning_analysis.numTrackValues(). This will either return the compile-time fixed value NumTrackValuesCTime, or the value which was set dynamically at runtime.

Eigen::Dynamic is never returned. See Tools::StaticOrDynamic.

Definition at line 417 of file mhrw_bin_err.h.

§ samplesSize

template<typename Params, typename LoggerType_>
const Tools::StaticOrDynamic<CountIntType, (SamplesSizeCTime==Eigen::Dynamic), (SamplesSizeCTime < 0 ? 0 : SamplesSizeCTime)> Tomographer::BinningAnalysis< Params, LoggerType_ >::samplesSize

The size of our samples buffer. (avoid using, might change in the future)

See BinningAnalysisParams::SamplesSizeCTime.

Warning
Avoid using this property, this might change in the future. [I may decide to change the algorithm, not flushing all after the last-level bin is full but half-flushing as we go, to save memory space and hopefully gain in terms of memory cache]

Definition at line 439 of file mhrw_bin_err.h.


The documentation for this class was generated from the following file: