  Roughly every 2-10 minutes, a pair of stellar mass black holes merge
somewhere in the Universe. A small fraction of these mergers are detected as
individually resolvable gravitational-wave events by advanced detectors such as
LIGO and Virgo. The rest contribute to a stochastic background. We derive the
statistically optimal search strategy for a background of unresolved binaries.
Our method applies Bayesian parameter estimation to all available data. Using
Monte Carlo simulations, we demonstrate that the search is both "safe" and
effective: it is not fooled by instrumental artefacts such as glitches, and it
recovers simulated stochastic signals without bias. Given realistic
assumptions, we estimate that the search can detect the binary black hole
background with about one day of design sensitivity data versus $\approx 40$
months using the traditional cross-correlation search. This framework
independently constrains the merger rate and black hole mass distribution,
breaking a degeneracy present in the cross-correlation approach. The search
provides a unified framework for population studies of compact binaries, which
is cast in terms of hyper-parameter estimation. We discuss a number of
extensions and generalizations including: application to other sources (such as
binary neutron stars and continuous-wave sources), simultaneous estimation of a
continuous Gaussian background, and applications to pulsar timing.
