Quantile Monte Carlo standard error

Inputs: $M\ge 2$ equal-length scalar chains with $N\ge 2$ observations each, a quantile level $\alpha \in (0,1)$, and the fixed levels $u_{-}=.16$ and $u_{+}=.84$ used for a one-standard-error approximation. Output: The approximate error interval and MCSE below, or not estimable.

1. If $N<4$, return not estimable. This procedure calls Split multi-chain autocorrelation effective sample size, whose input chains must have length at least $4$.
2. Pool the draws and compute ${\widehat{q}}_{\alpha }=Q_S(\alpha )$ by Left empirical quantile of finite draws. Form the indicator chains $I_{m,n}=\mathbf{1}\{x_{m,n}\le {\widehat{q}}_{\alpha }\}$, Blom-normalize their pooled values by Blom rank-normal scores for pooled draws, and compute $E_{\alpha }$ with Split multi-chain autocorrelation effective sample size. If the indicator chains are constant or $E_{\alpha }$ is not estimable, return not estimable.
3. Set $a=E_{\alpha }\alpha +1$ and $b=E_{\alpha }(1-\alpha )+1$. For $u\in (0,1)$, define the lower Beta quantile $Q_{a,b}(u)=\inf \{t\in [0,1]:F_{a,b}(t)\ge u\}$, where $F_{a,b}(t)=\mathrm{Beta}(a,b)((0,t])$ and Beta distribution supplies the probability measure. Let $v_{-}=Q_{a,b}(u_{-})$ and $v_{+}=Q_{a,b}(u_{+})$. For integers $n\ge 1$, the sets $E_n=(0,1\mathrm{/}n]$ in $(0,1)$ decrease to the empty set. By The real numbers, $\mathbb{R}$ is the Cauchy-sequence real field, and The Cauchy-sequence reals have the least-upper-bound property says it is a complete ordered field. Thus For every $\epsilon >0$ in a complete ordered field there is a natural $n\ge 1$ with $1\mathrm{/}n<\epsilon$ gives, for every $p\in (0,1)$, an $n$ with $1\mathrm{/}n<p$. The sets have finite Beta measure because the Beta law is a probability measure. Thus Continuity from above when one set has finite measure gives $F_{a,b}(1\mathrm{/}n)=\mathrm{Beta}(a,b)(E_n)\to 0$. Since $F_{a,b}$ is nondecreasing by Measures are monotone and $F_{a,b}(1)=1$, for each $u\in \{u_{-},u_{+}\}\subset (0,1)$ there is an $n$ with $F_{a,b}(1\mathrm{/}n)<u$; no $t\le 1\mathrm{/}n$ then belongs to the defining level set, while $t=1$ does. Hence $0<1\mathrm{/}n\le Q_{a,b}(u)\le 1$, so $v_{-},v_{+}\in (0,1]$ and both calls to $Q_S$ lie in its domain. Also $u_{-}<u_{+}$ makes the upper level set a subset of the lower one, so $v_{-}\le v_{+}$. The order-statistic definition makes $Q_S$ nondecreasing, hence $Q_S(v_{-})\le Q_S(v_{+})$.
4. Return the approximate Monte Carlo error interval $[Q_S(v_{-}),Q_S(v_{+})]$ and the estimate $\stackrel{}{\mathrm{MCSE}}({\widehat{q}}_{\alpha })=\frac{Q_S(v_{+})-Q_S(v_{-})}{2}\mathrm{.}$

All Beta shape parameters are positive because $E_{\alpha }>0$ and $0<\alpha <1$. This is a finite-run approximation, not a theorem that the interval has a prescribed coverage probability. Its statistical interpretation requires conditions for the relevant quantile and indicator-chain limit behavior; none are asserted here.