We analyze the probability distribution and spatial correlations around a stationary state of a general reaction-diffusion system. The stochastic description is based on a multivariate master equation. We use a WKB expansion of the probability density and determine the leading term, the stochastic potential, to second order in deviations from a homogeneous stationary state. For a system below, but near, a Turing instability, the spatial correlations become long range and display a macroscopic structure that will emerge after the bifurcation. We derive explicit expressions for both the probability density and the correlation function. For systems close to the instability, the correlation function is approximately an exponentially damped cosine function. We derive explicit expression for the correlation length and the amplitude of the correlation function; they are inversely proportional to the square root of the largest eigenvalue of the deterministic system. Our approach differs from earlier work in that systems with many chemical species are treated, asymptotic approximations are derived, correlations are given a geometrical picture in terms of eigenvectors of the Jacobian of an associated Hamiltonian system, and higher order terms in the stochastic potential are possible to obtain analytically (which we do not pursue). Results of the theoretical analysis are applied to the Sel＇kov model. Exact and approximate solutions at lowest order agree well.