We discuss the calculation of the double occupancy using dynamical mean-field theory in $\it finite$ dimensions. The double occupancy can be determined from the susceptibility of the auxiliary impurity model or from the lattice susceptibility. The former method typically overestimates, whereas the latter underestimates the double occupancy. We illustrate this for the square-lattice Hubbard model. We propose an approach for which both methods lead to identical results by construction and which resolves this ambiguity. This self-consistent dual boson scheme results in a double occupancy that is numerically close to benchmarks available in the literature.