kliff.uq.mcmc_utils

kliff.uq.mcmc_utils.mser(chain, dmin=1, dstep=10, dmax=-1, full_output=False)[source]

Estimate the equilibration time using marginal standard error rule (MSER).

This is done by calculating the standard error (square) of chain_d, where chain_d contains the last n-d element of the chain (n is the total number of iterations for each chain), for progresively larger d values, starting from dmin upto dmax, incremented by dstep. The SE values are stored in a list. Then we search the minimum element in the list and return the index of that element.

Parameters:
  • chain (ndarray) – (nsteps,) Array containing the time series.

  • dmin (Optional[int]) – Index where to start the search in the time series.

  • dstep (Optional[int]) – How much to increment the search is done.

  • dmax (Optional[int]) – Index where to stop the search in the time series.

  • full_output (Optional[bool]) – A flag to return the list of squared standard error.

Return type:

Union[int, dict]

Returns:

Estimate of the equilibration time using MSER. If full_output=True, then a dictionary containing the estimated equilibration time and the list of squared standard errors will be returned.

kliff.uq.mcmc_utils.autocorr(chain, *args, **kwargs)[source]

Use emcee package to estimate the autocorrelation length.

Parameters:
  • chain (ndarray) – (nwalkers, nsteps, ndim,) Chains from the MCMC simulation. Note that the burn-in time needs to be discarded prior to this calculation

  • args – Additional positional and keyword arguments of emcee.autocorr.integrated_time.

  • kwargs – Additional positional and keyword arguments of emcee.autocorr.integrated_time.

Return type:

ndarray

Returns:

Estimate of the autocorrelation length for each parameter.

kliff.uq.mcmc_utils.rhat(chain, time_axis=1, return_WB=False)[source]

Compute the value of \hat{r} proposed by Brooks and Gelman [BrooksGelman1998].

If the samples come from PTMCMC simulation, then the chain needs to be from one of the temperature only.

Parameters:
  • chain (ndarray) – The MCMC chain as a ndarray, preferrably with the shape (nwalkers, nsteps, ndim,). However, the shape can also be (nsteps, nwalkers, ndim,), but the argument time_axis needs to be set to 0.

  • time_axis (Optional[int]) – Axis in which the time series is stored (0 or 1). For emcee results, the time series is stored in axis 0, but for ptemcee for a given temperature, the time axis is 1.

  • return_WB (Optional[bool]) – A flag to return covariance matrices within and between chains.

Return type:

Union[float, Tuple[float, ndarray, ndarray]]

Returns:

The value of rhat. if return_WB=True, also returns matrices of covariance within and between the chains.

References

[BrooksGelman1998]

Brooks, S.P., Gelman, A., 1998. General Methods for Monitoring Convergence of Iterative Simulations. Journal of Computational and Graphical Statistics 7, 434455. https://doi.org/10.1080/10618600.1998.10474787