  We present a new class of statistical error reduction techniques for
Monte-Carlo simulations. Using covariant symmetries, we show that correlation
functions can be constructed from inexpensive approximations without
introducing any systematic bias in the final result. We introduce a new class
of covariant approximation averaging techniques, known as all-mode averaging
(AMA), in which the approximation takes account of contributions of all
eigenmodes through the inverse of the Dirac operator computed from the
conjugate gradient method with a relaxed stopping condition. In this paper we
compare the performance and computational cost of our new method with
traditional methods using correlation functions and masses of the pion,
nucleon, and vector meson in $N_f=2+1$ lattice QCD using domain-wall fermions.
This comparison indicates that AMA significantly reduces statistical errors in
Monte-Carlo calculations over conventional methods for the same cost.
