The ΛCDM model used in modern cosmology has major theoretical problems. To attempt to solve these many modified gravity theories have been proposed. We calculate velocity statistics (streaming pair-wise velocities and the velocity correlation function) for different models (Hu- Sawicki f(R) and symmetron) using data from n-body simulations. We find that the streaming velocities show the clearest deviations from GR. The F5 model deviates quite strongly, while for symmetron models the time of fifth force turn on is more important than the force range. We have not accounted for all the sources of variance and covariance, but suggest running many simulations with different realisations of the initial matter distribution to solve this. For future work we suggest calculating the statistics while binning the halos by their environments, and comparing our simulations to observations.