Recession velocity, redshift and different cosmological models? Equation 1 (page 5) in “Expanding Confusion: common misconceptions of cosmological horizons and the superluminal expansion of the universe” by Davis and Lineweaver gives the general relativistic relation between velocity and cosmological redshift:$$v_{rec}\left(t,z\right)=\frac{c}{R_{0}}\dot{R}\left(t\right)\int_{0}^{z}\frac{dz^{\prime}}{H\left(z^{\prime}\right)}.$$In Figure 1 (page 7) they use this equation to show the relationship between velocity and redshift for a range of Fridemann-Robertson-Walker models (using different values of $\Omega_{m}$ and $\Omega_{\Lambda}$). They assume $\dot{R}\left(t\right)=\dot{R}_{0}$.

How do they do that? In other words, what is the relationship between redshift and the density parameters $\Omega_{m}$, $\Omega_{\Lambda}$ and $\Omega_{r}$ (which they appear to ignore)?
 A: The thing you are missing is perhaps
$$H(z) = H_0\left((1-\Omega_\Lambda -\Omega_m)(1+z)^2 + (\Omega_\Lambda + \Omega_m)\frac{\rho(z)}{\rho_0}\right)^{1/2}$$
Where the densities $\rho$ dependend on the matter content, so they must be split into different epochs (of matter, radiation, etc. domination) to do the full integral but they just. This equation can be derived from the Friedmann equations. Once you plug in the values for the quantities measured today, $0$, you can compute the integral and obtain the velocity. You can find more details on the following link (https://ned.ipac.caltech.edu/level5/Peacock/Peacock3_2.html) or in any standard book on cosmology such as Mukhanov's "Physical Foundations of Cosmology" (2005), Chapter 2.
EDIT: The formula above holds for universes that are not spatially flat, $k\neq 0$. 
Starting from the Friedmann equation:
$$H^2(t) = \frac{8\pi G}{3}\rho(t)- \frac{k}{a^2} + \frac{\Lambda}{3}$$
It turns out that for the "Expanding Confusion" document, they take spatial curvature to be cero (flat universe) so that the density parameters add up to 1 (so you can always eliminate one of them, in your case to eliminate $\Omega_r$) and the universe to be composed only of matter, radiation and a cosmological constant so that with the usual definitions
$$\rho_{crit} = \frac{3H_0^2}{8\pi G}\qquad \textrm{and}\qquad \Omega_X = \frac{\rho_X}{\rho_{crit}},$$
with $0$ denoting the values today and knowing how the different components of the universe behave with respect to the scale factor allows us to rewrite the Friedmann equation in terms of the scale factor. Recall for matter $\rho_m\propto a^{-3}$, for radiation $\rho_r\propto a^{-4}$ and for dark energy we assume consant in this case, then
$$H(a) = H_0\sqrt{\Omega_\Lambda + \Omega_m a^{-3} + \Omega_r a^{-4}}$$
with $a_0 =1$. Now using the relation between scale factor and redshift
$$\frac{a_0}{a(t)}=1+z$$
in the previous formula, this is were the problems arise, using $\Omega_r = 1- \Omega_m - \Omega_\Lambda$ and a scaling of radiation as $a^{-2}$ you get their result, however using the correct scaling for radiation you get:
$$H(z) = H_0 (1+z)\left( 1 + \Omega_m z + \Omega_\Lambda\left(\frac{1}{(1+z)^2}-1\right) + \color{red}{2z\Omega_r + z^2\Omega_r} \right)^{1/2}$$
I hope this helps, but now I am also curios about how they get that formula...
A: I eventually realised that in order to answer my question I needed to understand the derivation of D&L's Equation 25: $$H\left(z\right)=H_{0}\left(1+z\right)\left[1+\varOmega_{m}z+\varOmega_{\varLambda}\left(\frac{1}{\left(1+z\right)^{2}}-1\right)\right]^{1/2}.$$
This is the Hubble parameter parameterised by redshift not time.  The critical density at present is$$\rho_{c}=\frac{3H_{0}^{2}}{8\pi G}.$$
Therefore $$H_{0}^{2}=\frac{8\pi G\rho_{c}}{3}$$
and $$\varOmega=\frac{\rho}{\rho_{c}}=\frac{8\pi G\rho}{3H_{0}^{2}}.$$
The Friedmann equation is
$$H^{2}\left(t\right)=\frac{8\pi G}{3}\rho-\frac{kc^{2}}{R^{2}}.$$
Divide by $H_{0}^{2}$ $$\frac{H^{2}}{H_{0}^{2}}=\frac{8\pi G}{3H_{0}^{2}}\rho-\frac{kc^{2}}{R^{2}H_{0}^{2}}$$
$$\frac{H^{2}}{H_{0}^{2}}=\varOmega-\frac{kc^{2}}{R^{2}H_{0}^{2}}.$$Where $\varOmega=\varOmega_{m}+\varOmega_{\varLambda}+\varOmega_{r}$. To find $-kc^{2}$ set$$-kc^{2}=H_{0}^{2}R_{0}^{2}-\frac{R_{0}^{2}8\pi G\rho}{3}$$ and get$$\frac{H^{2}}{H_{0}^{2}}=\varOmega+\frac{H_{0}^{2}R_{0}^{2}}{H_{0}^{2}R^{2}}-\frac{R_{0}^{2}8\pi G\rho}{3R^{2}H_{0}^{2}}$$
$$\frac{H^{2}}{H_{0}^{2}}=\varOmega+\frac{R_{0}^{2}}{R^{2}}\left(1-\frac{8\pi G\rho}{3H_{0}^{2}}\right)$$ $$\frac{H^{2}}{H_{0}^{2}}=\varOmega+\frac{R_{0}^{2}}{R^{2}}\left(1-\varOmega\right).$$
Now write as a function of redshift $z$. 
First, $$\frac{R_{0}^{2}}{R^{2}}=\left(1+z\right)^{2}$$giving$$H^{2}\left(z\right)=H_{0}^{2}\left[\varOmega+\left(1+z\right)^{2}\left(1-\varOmega\right)\right].$$
And$$\rho\left(z\right)=\rho_{m}\left(1+z\right)^{3}+\rho_{\varLambda}+\rho_{r}\left(1+z\right)^{4}$$
giving$$H^{2}\left(z\right)=H_{0}^{2}\left[\varOmega_{m}\left(1+z\right)^{3}+\varOmega_{\varLambda}+\varOmega_{r}\left(1+z\right)^{4}+\left(1+z\right)^{2}\left(1-\varOmega\right)\right].$$To simplify, let $\varOmega_{r}=0$ (as do Davis and Lineweaver)$$H^{2}\left(z\right)=H_{0}^{2}\left[\varOmega_{m}\left(1+z\right)^{3}+\varOmega_{\varLambda}+\left(1+z\right)^{2}\left(1-\varOmega_{m}-\varOmega_{\varLambda}\right)\right]$$ $$H^{2}\left(z\right)=H_{0}^{2}\left(1+z\right)^{2}\left[\varOmega_{m}\left(1+z\right)+\frac{\varOmega_{\varLambda}}{\left(1+z\right)^{2}}+1-\varOmega_{m}-\varOmega_{\varLambda}\right]$$ $$H^{2}\left(z\right)=H_{0}^{2}\left(1+z\right)^{2}\left[1+\varOmega_{m}z+\varOmega_{\varLambda}\left(\frac{1}{\left(1+z\right)^{2}}-1\right)\right]$$
$$H\left(z\right)=H_{0}\left(1+z\right)\left[1+\varOmega_{m}z+\varOmega_{\varLambda}\left(\frac{1}{\left(1+z\right)^{2}}-1\right)\right]^{1/2}.$$
