exhibit significant negative bias [@damerdji:1995; @vats:flegal:2018]. To overcome this, [@vats:flegal:2018] propose the lugsail BM estimator that offsets the negative bias in the opposite direction. For $r \geq 1$ and $0 \leq c < 1$, the lugsail BM estimator for the $k$th chain is: $$\hat{\Sigma}_{k,L} = \dfrac{1}{1 - c}\hat{\Sigma}_{k,b} - \dfrac{c}{1 - c}\hat{\Sigma}_{k, \lfloor b/r \rfloor} \; ,$$ where $\hat{\Sigma}_{k, \lfloor b/r \rfloor}$ is the BM estimator with batch size $\lfloor b/r \rfloor$. Setting $r = 1$ yields the estimator in .

Replicated batch means estimator {#eq:rbm}
================================

[@arg:and:2006] constructed a univariate RBM estimator for $\phi$-mixing processes by essentially concatenating the $m$ Markov chains. We construct a multivariate RBM estimator and prove critical results for their application to MCMC, under weaker mixing conditions. For a batch size of $b$, the multivariate RBM estimator is $$\label{eq:RBM}
    \hat{\Sigma}_{R, b} = \dfrac{b}{am-1}\sum_{k=1}^{m}\sum_{l=1}^{a}\left(\bar{Y}_{k,l} - \hat{\mu}\right)\left(\bar{Y}_{k,l} - \hat{\mu} \right)^T \; .$$ The idea behind RBM is to pool all the $am$ batch mean vectors and measure the squared deviation from the overall mean, $\hat{\mu}$. Thus, $\hat{\Sigma}_{R, b}$ acknowledges the high variability in the batch mean vectors when each of the Markov chains is in different parts of the state space. Regular RBM can also suffer from 