We show that any 3+1-dimensional Milne model is future nonlinearly, asymptotically stable in the set of solutions to the Einstein-Vlasov system. For the analysis of the Einstein equations we use the constant-mean-curvature-spatial-harmonic gauge. For the distribution function the proof makes use of geometric L^2-estimates based on the Sasaki-metric. The resulting estimates on the energy momentum tensor are then upgraded by employing the natural continuity equation for the energy density.