In this paper we discuss a Schwinger-Dyson (SD) approach for determining the time evolution of the unequal time correlation functions of a nonequilibrium classical field theory, where the classical system is described by an initial density matrix at time
t=0. We focus on $\ensuremath{\lambda}{\ensuremath{\varphi}}^{4}$ field theory in
1+1 space-time dimensions where we can perform exact numerical simulations by sampling an ensemble of initial conditions specified by the initial density matrix. We discuss two approaches. The first, the bare vertex approximation (BVA), is based on ignoring vertex corrections to the SD equations in the auxiliary field formalism relevant for
1/N expansions. The second approximation is a related approximation made to the SD equations of the original formulation in terms of $\ensuremath{\varphi}$ alone. We compare these SD approximations as well as a Hartree approximation with exact numerical simulations. We find that both approximations based on the SD equations yield good agreement with exact numerical simulations and cure the late time oscillation problem of the Hartree approximation. We also discuss the relationship between the quantum and classical SD equations.