In the Wikipedia article on Monte Carlo integration, they give an example of stratified sampling for integrating the function
$$ f(x,y) = \begin{cases} 1 \;,\; x^2+y^2 < 1 \\ 0 \;,\; x^2+y^2 \geq 1\end{cases} $$
over the unit square $x,y \in [-1, 1]$. This integral basically represents the area of the unit circle.
So I was wondering whether stratified sampling is even useful for such volume estimations, also in higher dimensions. The standard error for the Monte Carlo integral of binary functions (i.e. either $1$ or $0$), which represents the volume of the region where $f \equiv 1$, can be expressed in the following way:
$$ \begin{aligned} \sigma_N^2 &= \frac{1}{N-1}\sum_{i=1}^{N}\left[f_i - \langle f\rangle \right]^2 \\ &= \frac{1}{N-1}\left[\sum_i f_i^2 - 2N\langle f\rangle^2 + N\langle f\rangle^2 \right] \\ &= \frac{1}{N-1}\left[\sum_i f_i^2 - N\langle f\rangle^2\right] \\ &= \frac{N}{N-1}\left[\langle f\rangle - \langle f\rangle^2\right] \\ \Rightarrow \mathrm{standard\,error\,of\,}\langle f\rangle &= \frac{\sigma_N}{\sqrt{N}} = \sqrt{\frac{\langle f\rangle - \langle f\rangle^2}{N-1}} \end{aligned} $$
where $f_i$ is the $i$-th Monte Carlo sample and $\sum_i f_i = N\langle f\rangle$ and $f_i^2 = f_i$ has been used.
The function $\langle f\rangle - \langle f\rangle^2$ behaves in the following way, so the maximum error is obtained for $\langle f\rangle = \frac{1}{2}$:
The error relative to the actual estimate, i.e.
$$ \frac{\mathrm{standard\,error}}{\langle f\rangle} = \sqrt{\frac{\frac{1}{\langle f\rangle} - 1}{N-1}} $$
is largest for $\langle f\rangle \ll 1$.
How does stratified sampling help in this case to reduce the standard error of the volume estimate?