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, wherechain_dcontains the lastelement of the chain (n is the total number of iterations for each chain), for progresively larger d values, starting from
dminuptodmax, incremented bydstep. 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
emceepackage 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 calculationargs – 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
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