Solution (source code)

= Solution

The <Einstein-de Sitter universe> has $a\propto t^{2/3}$ and $H=\dot a/a\propto a^{-3/2}$. Normalizing $a_0=1$ therefore gives $H(z)=H_0(1+z)^{3/2}$. Differentiating $1+z=1/a$ gives
$$
\frac{dt}{dz}=-\frac1{(1+z)H(z)}=-\frac1{H_0(1+z)^{5/2}}.
$$
In the local comoving rest frame, a <photon> crosses proper path length $d\ell=c|dt|$. Galaxies with proper <number density> $n(z)$ and geometric <cross section> $\sigma_G=\pi r_G^2$ contribute the mean number of intersections $d\tau=n(z)\sigma_Gd\ell$. The infinitesimal probability of an intersection is the same to first order, giving
$$
\boxed{dP=d\tau=\frac{c\pi r_G^2n(z)}{H_0(1+z)^{5/2}}\,dz.}
$$
The radius is a physical radius, and $n(z)$ is a proper <number density> evaluated at epoch $z$, rather than a density already multiplied by $dz$. For a finite path with overlaps, the integrated mean count $\tau$ need not itself be the probability of at least one intersection.