  We consider the nonparametric multivariate isotonic regression problem, where
the regression function is assumed to be nondecreasing with respect to each
predictor. Our goal is to construct a Bayesian credible interval for the
function value at a given interior point with assured limiting frequentist
coverage. We put a prior on unrestricted step-functions, but make inference
using the induced posterior measure by an "immersion map" from the space of
unrestricted functions to that of multivariate monotone functions. This allows
maintaining the natural conjugacy for posterior sampling. A natural immersion
map to use is a projection via a distance, but in the present context, a block
isotonization map is found to be more useful. The approach of using the induced
"immersion posterior" measure instead of the original posterior to make
inference provides a useful extension of the Bayesian paradigm, particularly
helpful when the model space is restricted by some complex relations. We
establish a key weak convergence result for the posterior distribution of the
function at a point in terms of some functional of a multi-indexed Gaussian
process that leads to an expression for the limiting coverage of the Bayesian
credible interval. Analogous to a recent result for univariate monotone
functions, we find that the limiting coverage is slightly higher than the
credibility, the opposite of a phenomenon observed in smoothing problems.
Interestingly, the relation between credibility and limiting coverage does not
involve any unknown parameter. Hence by a recalibration procedure, we can get a
predetermined asymptotic coverage by choosing a suitable credibility level
smaller than the targeted coverage, and thus also shorten the credible
intervals.
