  Based on two independent samples X_1,...,X_m and X_{m+1},...,X_n drawn from
multivariate distributions with unknown Lebesgue densities p and q
respectively, we propose an exact multiple test in order to identify
simultaneously regions of significant deviations between p and q. The
construction is built from randomized nearest-neighbor statistics. It does not
require any preliminary information about the multivariate densities such as
compact support, strict positivity or smoothness and shape properties. The
properly adjusted multiple testing procedure is shown to be sharp-optimal for
typical arrangements of the observation values which appear with probability
close to one. The proof relies on a new coupling Bernstein type exponential
inequality, reflecting the non-subgaussian tail behavior of a combinatorial
process. For power investigation of the proposed method a reparametrized
minimax set-up is introduced, reducing the composite hypothesis "p=q" to a
simple one with the multivariate mixed density (m/n)p+(1-m/n)q as infinite
dimensional nuisance parameter. Within this framework, the test is shown to be
spatially and sharply asymptotically adaptive with respect to uniform loss on
isotropic H\"older classes. The exact minimax risk asymptotics are obtained in
terms of solutions of the optimal recovery.
