  We study the effect of two Modified Gravity (MG) theories, $f(R)$ and nDGP,
on three probes of large-scale structure, the real space power spectrum
estimator $Q_0$, bispectrum and voids, and validate fast approximate COLA
simulations against full $N$-body simulations for the prediction of these
probes. We find that using the first three even multipoles of the redshift
space power spectrum to estimate $Q_0$ is enough to reproduce the MG boost
factors of the real space power spectrum for both halo and galaxy catalogues.
By analysing the bispectrum and reduced bispectrum of Dark Matter (DM), we show
that the strong MG signal present in the DM bispectrum is mainly due to the
enhanced power spectrum. We warn about adopting screening approximations in
simulations as this neglects non-linear contributions that can source a
significant component of the MG bispectrum signal at the DM level, but we argue
that this is not a problem for the bispectrum of galaxies in redshift space
where the signal is dominated by the non-linear galaxy bias. Finally, we
perform void-finding on our galaxy mock catalogues by the ZOBOV watershed
algorithm. To apply a linear model for Redshift-Space Distortion (RSD) in the
void-galaxy cross-correlation function, we first examine the effects of MG on
the void profiles entering into the RSD model. We find relevant MG signals in
the integrated-density, velocity dispersion and radial velocity profiles in the
nDGP theory. Fitting the RSD model for the linear growth rate, we recover the
linear theory prediction in an nDGP model, which is larger than the
$\Lambda$CDM prediction at the $3 \sigma$ level. In $f(R)$ theory we cannot
naively compare the results of the fit with the linear theory prediction as
this is scale-dependent, but we obtain results that are consistent with the
$\Lambda$CDM prediction.
