In geophysical imaging, uncertainty quantification is crucial for decision making. 4D seismic imaging aims to accurately recover changes that take place within a reservoir. These changes are typically characterized by their magnitude and their extent. We perform a Bayesian inversion using a Metropolis Hastings algorithm to sample our posterior distribution of 4D velocity models given observed data. To model the 4D change we use a discrete cosine transformation, and attempt to recover the lowest frequency coefficients, so that we can model realistic changes with only a few degrees of freedom. Unlike most of uncertainty quantification methodologies that use expensive forward solvers, we speed up our computations by using a numerically exact local acoustic solver.