The most plausible cause is strong dependence between the latent variables and covariance hyperparameters in their joint posterior distribution. In the usual conditionally independent Gaussian process classification model, put , , and . The probit regression likelihood is , while has a multivariate normal distribution with covariance matrix .
For distinct inputs and positive length scales, is positive definite. Its Gaussian probability density function contributes
The unit-rate exponential prior distribution on therefore gives
Its scale is tightly tied to the current magnitude and shape of . Changing changes the eigenvectors and eigenvalues of its covariance matrix; a latent vector typical under the old covariance may be very atypical under a substantially different one. Holding fixed can thus restrict the range of a hyperparameter update, and holding the hyperparameters fixed restricts the next latent-function update. Even exact sampling of each full conditional distribution can move slowly along a narrow ridge of the joint posterior distribution.
Binary observations provide limited information about the absolute size of large correctly signed latent values. This can leave substantial scale uncertainty and strengthen the dependence. Long correlation lengths can also make the covariance matrix nearly singular and the latent coordinates highly dependent. These are plausible explanations for a large mixing time, not proofs that every data set causes slow convergence. A blocked Gibbs sampler removes dependence within its chosen blocks, but not dependence between them. Repeated identical inputs should be represented by one shared latent value; otherwise the nominal Gaussian matrix is singular and the displayed inverse formula must be replaced by a representation on its support.