We extend the scale-dependent Gaussian Streaming Model (GSM) to produce analytical predictions for the anisotropic redshift-space correlation function for biased tracers in modified gravity models. Employing the Convolution Lagrangian Perturbation Theory (CLPT) re-summation scheme, with a local Lagrangian bias schema provided by the peak-background split formalism, we predict the necessary ingredients that enter the GSM, the real-space halo pairwise velocity and the pairwise velocity dispersion. We further consider effective field theory contributions to the pairwise velocity dispersion in order to model correctly its large scale behavior. We apply our method on two widely-considered modified gravity models, the chameleon-screened f(R) Hu-Sawicki model and the nDGP Vainshtein model and compare our predictions against state-of-the-art N-body simulations for these models. We demonstrate that the GSM approach to predict the monopole and the quadrupole of the redshift-space correlation function for halos, gives very good agreement with the simulation data, for a wide range of screening mechanisms, levels of screening and halo masses at z=0.5 and z=1. Our work shows the applicability of the GSM, for cosmologies beyond GR, demonstrating that it can serve as a powerful predictive tool for the next stage of cosmological surveys like DESI, Euclid, LSST and WFIRST.