More Rigorous Error Analysis

[emerali]

Vanilla Jackknife/Bootstrap assume IID samples, and thus underestimate standard errors (this is true even for simple estimators which are just means). There are Bootstrap methods specifically for time-series data which could be useful (e.g. Stationary Bootstrap).
There may be a way to turn one of these into a "streaming" bootstrap

Also: autocorrelation time we're computing seems to be larger than what's computed by other packages. Should try calculating autocorrelation time using a binning analysis instead as it may be more "stable" than the FFT based method we're using.

Should compute autocorrelation time for all observables in order to get accurate error bars.

Useful Papers:

• "Everything you wanted to know about Data Analysis and Fitting but were afraid to ask": https://arxiv.org/abs/1210.3781
• "Efficient estimation of autocorrelation spectra": https://arxiv.org/abs/1810.05079

[emerali]

Update:

Simple Estimators:

• BinningAnalysis.jl's LogBinner allows us to perform a binning analysis online and it works pretty well. Will use this for production runs.

Complicated/Non-linear estimators:
Stationary Bootstrap (in principle) can be performed in a streaming/online fashion, provided you've already set the average block-length parameter. However, computing the optimal block-length requires access to the entire time-series.

Fortunately, it seems that SB isn't extremely sensitive to an incorrect block-length (see Politis and Romano (1994)). In addition, the optimal block-length scales as b ~ cbrt(N) where N is the full length of the time-series. So, we can likely do something like:

1. After equilibriation, run a single batch of MC steps (N samples), keeping all the measurement needed for the bootstrapping procedure.
2. Compute the optimal block-length from these measurements (if there's more than 1 quantity the estimator depends on, compute the block-length for each quantity, then take the median of the block-lengths).
3. Rescale the block-length: b' = cbrt(B) b where B is the number of batches
4. Initialise an online SB object, setting the average block length to b', and adding this first batch of data to it
5. run the rest of the MC simulation, putting new measurements into the SB object instead of keeping them in memory.