The most direct approach to this problem starts with the Dyson-Schwinger equations for field correlators and truncates this hierarchy at a given order eg. through a 1/N expansion. Beyond mean fields such equations become non-local in time. This is a huge constraint on the possibility of their solution.
To circumvent this problem it was proposed that ab initio equal time field correlators could be used. The (Dyson-Schwinger) equation of motion for the the generating functional for equal time Green's functions can be written down for full, connected or 1PI field correlators. If these generating functionals have a series expansion an approximation can be achieved by truncating the series at a given order. Analogously another truncation can be achieved by using collecting terms of order 1/N to a certain order.
In order to study the properties of these aproximations and their validity against an exact quantum system I discussed the time evolution of the coupled hierarchy of equal time Green's functions for a system of N anharmonic oscillators.
This can be compared to the quantum roll for these oscillators. Under the assumption of radial symmetry, which is equivalent to a N dimensional single oscillator problem, the Schroedinger equation for the system is solved.
The truncated Green's functions beyond the second order truncation, have the property that they do not correspond to a positive definite probability. I discussed some of their general properties.
I also considered what possible variational truncations might be successful by projecting the exact answer on a particular set of basis wave functions.
A preprint is in preparation and the corresponding link will be added as soon as possible.