def apply_inverse(self, y):
        """
        Self-consistently apply the inverse of the computed kernel matrix to
        some vector or matrix of samples. This method subtracts the mean,
        sorts the samples, then returns the samples in the correct (unsorted)
        order.

        :param y: ``(nsamples, )`` or ``(nsamples, K)``
            The vector (or matrix) of sample values.

        """
        self.recompute(quiet=False)
        r = np.array(y, dtype=np.float64, order="F")
        r = self._check_dimensions(r, check_dim=False)

        # Broadcast the mean function
        m = [slice(None)] + [np.newaxis for _ in range(len(r.shape) - 1)]
        r -= self._call_mean(self._x)[m]

        # Do the solve
        if len(r.shape) == 1:
            b = self.solver.apply_inverse(r, in_place=True).flatten()
        else:
            b = self.solver.apply_inverse(r, in_place=True)
        return b
Self-consistently apply the inverse of the computed kernel matrix to
        some vector or matrix of samples. This method subtracts the mean,
        sorts the samples, then returns the samples in the correct (unsorted)
        order.

        :param y: ``(nsamples, )`` or ``(nsamples, K)``
            The vector (or matrix) of sample values.