We study the dynamics of small inhomogeneities in an expanding universe collapsing to form bound structures using full solutions of the Einstein-Vlasov (N -body) equations. We compare these to standard Newtonian N -body solutions using quantities defined with respect to fiducial observers in order to bound relativistic effects. We focus on simplified initial conditions containing a limited range of length scales, but vary the inhomogeneities from small magnitude, where the Newtonian and GR calculations agree quite well, to large magnitude, where the background metric receives an order one correction. For large inhomogeneities, we find that the collapse of overdensities tends to happen faster in Newtonian calculations relative to fully general-relativistic ones. Even in this extreme regime, the differences in the spacetime evolution outside the regions of large gravitational potential and velocity are small. For standard cosmological values, we corroborate the robustness of Newtonian N -body simulations to model large scale perturbations and the related cosmic variance in the local expansion rate.