Estimating $f(\mathbb{E}[X])$ with a guaranteed error performance Given "black-box" sample access to a random variable $X$**, are there results that give an algorithm that approximates $f(\mathbb{E}[X])$ with a user-specified error bound, ideally using as few i.i.d. samples from $X$ as possible? (See the section "Update (May 11)".)
In this question:

*

*$X$ is a random variable that does not take on a single value with probability 1.

*$f(x)$ is a known function belonging to a given class of functions.

The algorithm should—

*

*Ensure the expected (absolute) error and/or mean squared error of the estimate is within a user-specified error tolerance ($\epsilon$), or if that is not possible,

*return an estimate that is within a user-specified tolerance ($\epsilon$) on the absolute error or relative error with probability greater than 1 minus $\delta$, where $\delta$ is user-specified (relative error means $|\hat\mu/f(\mathbb{E}[X])-1|$ where $\hat\mu$ is the estimate).

Before the Update
In this section, $X$ lies in the interval [0, 1].
The classes of functions $f$ that I care about are:

*

*C1: Continuous functions that map [0, 1] to [0, 1] and are polynomially bounded (meaning that $f(x)$ and $1-f(x)$ are both bounded from below by $\min(x^n, (1-x)^n)$ for some integer $n$) (Keane and O'Brien 1994).

*C2: Continuous functions that map [0, 1] to [0, 1] and are not necessarily polynomially bounded.

*C3: Piecewise continuous functions that map [0, 1] to [0, 1] (notably, each piece's domain must be a non-trivial interval and all the pieces must cover all of [0, 1]).

Note the following:

*

*The following algorithms are different from the algorithms asked for here:

*

*The Gamma Bernoulli Approximation scheme (Huber 2017) as well as the algorithm of Kunsch et al. 2019 both estimate the mean of $X$, namely $\mathbb{E}[X]$, rather than a function of that mean, namely $f(\mathbb{E}[X])$.

*Some algorithms produce unbiased estimates of $f(\mathbb{E}[X])$, but not with a user-specified error (Jacob and Thiery 2015).



*Algorithms like those being asked for here are especially useful because they can help build so-called "approximate Bernoulli factories", or algorithms that approximately sample the probability $f(\lambda)$ given a coin with probability of heads of $\lambda$.

*I suspect that the algorithm's performance will depend on the "smoothness" of $f(X)$, including its modulus of continuity (see also (Holtz et al. 2011)).

Update (May 11)
The responses so far helped me greatly. I now have the following more focused questions on this topic:
Let $f(x)$ be a continuous function.

*

*Given the setting at the top of the question, is there an algorithm, besides the algorithm of Kunsch et al. (2019), that can approximate $\mathbb{E}[X]$ (or $f(\mathbb{E}[X])$) with either a high probability of a "small" absolute error or one of a "small" relative error, when the distribution of $X$ is unbounded, and additional assumptions on the distribution of $X$ apply, such as—

*

*being unimodal (having one peak) and symmetric (mirrored on each side of the peak), and/or

*following a geometric distribution, and/or

*having decreasing or nonincreasing probabilities?

Notice that merely having finite moments is not enough (Theorem 3.4, Kunsch et al.).  Here, the accuracy tolerances for small error and high probability are user-specified.  My article on estimation algorithms already gives a relative-error algorithm for the geometric distribution in a note.


*How can the method in the answer by "guy" be adapted to discontinuous functions $g$, so that the algorithm finds $g(\mathbb{E}[X])$ with either a high probability of a "small" absolute error or one of a "small" relative error at all points in [0, 1] except at a "negligible" area around $g$'s discontinuities? Is it enough to replace $g$ with a continuous function $f$ that equals $g$ everywhere except at that "negligible" area? Here, the accuracy tolerances for small error, high probability, and "negligible" area are user-specified. Perhaps the tolerance could be defined as the integral of absolute differences between $f$ and $g$ instead of "negligible area"; in that case, how should the continuous $f$ be built? (In this question, $X$ lies in the interval [0, 1].)
Motivation
The answers to this question will help improve my article on randomized estimation algorithms, which describes practical estimation algorithms in a way that programmers can easily implement them.
References

*

*Keane, M. S., and O'Brien, G. L., "A Bernoulli factory", ACM Transactions on Modeling and Computer Simulation 4(2), 1994.

*Huber, M., 2017. A Bernoulli mean estimate with known relative error distribution. Random Structures & Algorithms, 50(2), pp.173-182. (preprint in arXiv:1309.5413v2 [math.ST], 2015).

*Kunsch, Robert J., Erich Novak, and Daniel Rudolf. "Solvable integration problems and optimal sample size selection." Journal of Complexity 53 (2019): 40-67. Also in https://arxiv.org/pdf/1805.08637.pdf.

*Pierre E. Jacob. Alexandre H. Thiery. "On nonnegative unbiased estimators." Ann. Statist. 43 (2) 769 - 784, April 2015.

*Holtz, O., Nazarov, F., Peres, Y., "New Coins from Old, Smoothly", Constructive Approximation 33 (2011).


** And possibly a second source of randomness (such as unbiased random bits). However, the random bits could also come from the random variable $X$ itself, assuming $X$ doesn't take on the same value with probability 1, via randomness extraction techniques such as the von Neumann extractor (1951). On the other hand, if $X$ can be degenerate, this second source will be necessary.
 A: C1: Yes, by the argument for C2.
C2: Yes, this is possible. The key result here is that any continuous function $f$ on $[0,1]$ is also uniformly continuous: for any $\epsilon$, there is some $\gamma$ such that if $|x_1-x_2|<\gamma$ then $|f(x_1)-f(x_2)|<\epsilon$. So it is enough to have a probability of $1-\delta$ of approximating $E[X]$ to within $\gamma$, and then applying $f$ will give a probability of $1-\delta$ of approximating $f(E[X])$ to within $\epsilon$.
C3:
No algorithm will work for all functions of the class C3. If we had an algorithm, we could use it on a discontinuous function $f$ to tell whether $E[X]>\frac12$, but there are always borderline cases that won't work with the desired probability.
In more detail, suppose the user specifies:

*

*$f(x) = 1[x>\frac12]$, i.e. $1$ if $x>\frac12$, or $0$ otherwise.

*$\epsilon=\frac15$ as the accuracy -- so getting $f$ to an accuracy of $\epsilon$ is the same as getting $f$ exactly

*$\delta=\frac14$, i.e. a probability of $\frac34$ of achieving this accuracy

Suppose the algorithm replies that either:

*

*$n$ observations are required, or

*with probability of at least $1-\frac\delta2$, at most $n$ observations are required.

Now consider a Bernoulli variable $X$ which is equally likely to be 0 or 1.
Let $S$ be the sample mean of $X$ after $n$ observations. Then $$P\left[S=\frac12\right]=\binom{n}{n/2}2^{-n}\simeq\sqrt{\frac2{\pi n}}$$
So there is close to a 50% chance that $S>\frac12$, and close to a 50% chance that $S<\frac12$. There is no way to have a $\frac34$ chance of getting the right answer, or even a $\frac34-\frac18$ chance of getting the right answer when we allow a $\frac18$ possibility of not enough data.
With a slightly longer argument, we could show that the same pattern holds when $X$ is Bernoulli and $E[X]$ is anything in the range $\frac12\pm\sqrt{n}/16$.
So even if we allow a $\frac\delta2$ probability that $n$ observations are insufficient, there is still no way to get $f(E[X])$ which will be accurate to within $\epsilon$ the needed $1-\delta$ of the time.
A: Elaborating on my comment to the OP, I'll give an explicit answer when $f$ is continuous.
Let $\psi(\epsilon) = \sup\{\psi : |x - y| < \psi \Longrightarrow |f(x) - f(y)| < \epsilon\}$. Because $f$ is continuous on $[0,1]$ it is also uniformly continuous, implying that $\psi(\epsilon) \downarrow 0$ as $\epsilon \downarrow 0$.
By Hoeffding's inequality we have
$$
  \Pr(|f(\bar X) - f(E(X))| \ge \epsilon)
  \le \Pr(|\bar X - E(X)| \ge \psi(\epsilon))
  \le 2\exp\{-2n\psi(\epsilon)^2\},
$$
where $\bar X = \frac{1}{n} \sum_i X_i$. We can make this less than $\delta$ by taking
$$
n \ge \frac{\log(2/\delta)}{2 \psi(\epsilon)^2}.
$$
In the special case where $f$ is Lipschitz continuous, i.e., $|f(x) - f(y)| \le L |x - y|$, we have $\psi(\epsilon) \ge \epsilon / L$, and the bound simplifies to
$$
  n \ge \frac{L^2 \log(2 / \delta)}{2 \epsilon^2}.
$$
This can be weakened to allow for non-Lipschitz functions; for example, you get a similar bound by assuming $|f(x) - f(y)| \le L |x - y|^\alpha$ for some $\alpha \le 1$ (i.e., that $f$ is Holder continuous), which should cover just about every continuous function you might care about, at least on $[0,1]$.
