Lensed event sampling (Analytical formulation)

Initial setup

Note: I will interchangeably use terms like observable events and detectable events.

Define all the parameters involved.

Source parameters: $\theta \in \{$ $m_1$ (mass of the heavier one), $m_2$ (mass of the lighter one), $\iota$ (inclination-angle), $\phi$ (phase-of-coalescence), $\psi$ (polarization-angle), $ra$ (right-ascension), $dec$ (declination) $\}$ and $z_s$ : red-shift of the source.
Lens parameters: $\theta_L \in \{$ $\sigma$ (velocity-dispersion), $q$ (axis-ratio), $\psi$ (axis-rotation), $\gamma$ (spectral-index), $[\gamma_1,\gamma_2]$ (external-shear), $[e_1,e_2]$ (ellipticity), $\beta$ (source position) $\}$
$z_L$ : red-shift of the galaxy lens
image param: $\{$ $\mu_i$ (magnification), $dt_i$ (time-delay), $n_i$ (morse-phase) $\}$. There is subscript $i$ because there can be multiple images.

Given $d N^L_{obs}(z_s)$ is the number of lensed GW detectable events from sources at red-shift $z_s$ in a spherical shell of thickness $d z_s$, then, let rate of lensing (number of lensed events happening per unit time) is given by,

\begin{equation} \begin{split} \mathcal{R}_L &= \int_{z_{min}}^{z_{max}} \frac{d N^L_{obs}(z_s)}{d \tau} \\ &= \int_{z_{min}}^{z_{max}} \frac{d N^L_{obs}(z_s)}{d \tau\;dV_c} \frac{dV_c}{dz_s} dz_s \end{split} \tag{1} \end{equation}

$\frac{d N^L_{obs}(z_s)}{d \tau \;dV_c}$ is the observed merger rate density at source-frame, and $\frac{dV_c}{dz_s}$ is the differential co-moving volume at red-shift $z_s$. After taking care of time-dilation, the expression looks,

\begin{equation} \begin{split} \mathcal{R}_L &= \int_{z_{min}}^{z_{max}} \frac{d N^L_{obs}(z_s)}{d t\;dV_c} \frac{1}{1+z_s} \frac{dV_c}{dz_s} dz_s \\ &= \int_{z_{min}}^{z_{max}} R^L_{obs}(z_s) dz_s \end{split} \tag{2} \end{equation}

Observed merger rate density of lensed GW events at source red-shift $z_s$ (detector-frame) is given by $R^L_{obs}(z_s) = \frac{d N^L_{obs}(z_s)}{d t\;dV_c}$. And, let $R(z_s)$ be the rate of the merger of unlensed events (source frame) at $z_s$, regardless of whether it is detectable or not.

\begin{equation} \begin{split} \mathcal{R}_L &= \int_{z_{min}}^{z_{max}} R^L_{obs}(z_s) dz_s \\ &= \int_{z_{min}}^{z_{max}} R(z_s) P(obs,SL|z_s) dz_s \end{split} \tag{3} \end{equation}

$P(obs,SL|z_s)$ is the probability of observing a strong lensing event at red-shift $z_s$.

Optical depth

Strong lensing probability (where does it come from?)

\begin{equation} \begin{split} \mathcal{R}_L &= \int_{z_{min}}^{z_{max}} R(z_s) P(obs|z_s, SL) P(SL|z_s) dz_s \end{split} \tag{4} \end{equation}

Probability of observing an event given that it is located at redshift $z_s$ and it’s strongly lensed: $P(obs|z_s, SL)$. Strong lensing probability with source at redshift $z_s$ (optical depth): $P(SL|z_s)$. Now, using Bayes’ theorem,

\begin{equation} \begin{split} P(SL|z_s) &= \frac{P(z_s|SL) P(SL)}{P(z_s)} \\ P(z_s) P(SL|z_s) &= P(z_s|SL) P(SL) \\ \frac{R(z_s)}{N_1} P(SL|z_s) &= P(z_s|SL) P(SL) \\ P(SL|z_s) &= P(z_s|SL) P(SL) \frac{N_1}{R(z_s)} \end{split} \tag{5} \end{equation}

Normalizing factor: $N_1 = \int_{z_{min}}^{z_{max}} R(z_s) dz_s$.

Taking care of normalization

Similarly, when lensing condition applied, let $N_2 = \int_{z_{min}}^{z_{max}} R(z_s) P(SL|z_s) dz_s$.

\begin{equation} \begin{split} P(SL) &= \int_{z_{min}}^{z_{max}} P(SL|z_s) P(z_s) dz_s \\ &= \int_{z_{min}}^{z_{max}} P(SL|z_s) \frac{R(z_s)}{N_1} dz_s \\ &= \frac{N_2}{N_1} \end{split} \tag{6} \end{equation}

Now putting together, equation 5 becomes,

\begin{equation} \begin{split} \frac{R(z_s)}{N_1} P(SL|z_s) &= P(z_s|SL) \frac{N_2}{N_1} \\ R(z_s) P(SL|z_s) &= N_2 P(z_s|SL) \end{split} \tag{7} \end{equation}

Optical depth in rate equation

Replace the above result in the integrand of Equation 4. This also takes care of the normalizing factor. Note that $P(z_s|SL)$ is a normalized pdf of source red-shifts, $z_s$, conditioned on strong lensing.

\begin{equation} \begin{split} \mathcal{R}_L &= N_2 \int_{z_{min}}^{z_{max}} P(z_s|SL) P(obs|z_s, SL) dz_s \\ & \text{consider } \int \rightarrow \int_{zl}\int_{\beta}\int_{\theta}\int_{\theta_L} \\ \mathcal{R}_L &= N_2 \int P(z_s|SL) P(obs|\theta, \theta_L, \beta, z_s, SL) \\ & \;\; P(\beta|\theta_L, z_s, SL) P(\theta_L|z_s, SL) P(\theta) d\beta d\theta d\theta_L dz_s \end{split} \tag{8} \end{equation}

For $P(z_s|SL)= \frac{R(z_s)}{N_2} P(SL|z_s)$, from equation 7, I have considered ‘optical depth’ ($\tau(z_s) = P(SL|z_s)$) is a function of $z_s$ only. Otherwise, we need to consider the cross-section ($P(SL|z_s, \theta_L)$) which will not be discussed here. Below shows how to get $P(SL|z_s)$, i.e., the probability of strong lensing of source at $z_s$. $dN(z_l)$ is the number of galaxy lenses at red-shift $z_l$ (in $d z_l$).

General formulation of optical depth expression

\begin{equation} \begin{split} P(SL|z_s) &= \int^{z_s}_0 \frac{P(SL|\theta_L, z_s)}{4\pi} dN(z_l) \\ P(SL|z_s) &= \int \frac{P(SL|z_s, z_l, \sigma, q)}{4\pi} \frac{dN(z_l)}{dz_l d\sigma dq} dz_l d\sigma dq \\ P(SL|z_s) &= \int \frac{P(SL|z_s, z_l, \sigma, q)}{4\pi} \frac{dN(z_l)}{dV_c d\sigma dq} \frac{dV_c}{dz_l} dz_l d\sigma dq \\ & \text{writing the cross-section $P(SL|z_s, z_l, \sigma, q)$ as $\phi$} \\ P(SL|z_s) &= \int \frac{\phi}{4\pi} \frac{dN(z_l)}{dV_c d\sigma dq} \frac{dV_c}{dz_l} dz_l d\sigma dq \\ \end{split} \tag{9} \end{equation}

Optical depth for SIS lens

Consider the SIS case of Haris et al. 2018. Take $\phi$ as $\phi_{SIS}$.

\begin{equation} \begin{split} P(SL|z_s) &= \int \frac{\phi_{SIS}}{4\pi} \frac{dN(z_l)}{dV_c d\sigma} \frac{dV_c}{dz_l} dz_l d\sigma \\ P(SL|z_s) &= \int \frac{\pi \theta_E^2}{4\pi} \big< n \big>_{\sigma\in P(\sigma)} P(\sigma) \frac{dV_c}{dz_l} dz_l d\sigma \end{split} \tag{10} \end{equation}

The cross-section of the SIS lens is $\pi \theta_E^2$, where $\theta_E$ is the Einstein radius. Haris has considered the number density of the lens, $<n>_{\sigma\in P(\sigma)}$ and the PDF of velocity dispersion, $P(\sigma)$ is independent of $z_l$. Take $<n>_{\sigma\in P(\sigma)}=n_o=8\times 10^{-3} h^3 Mpc^{-3}$.

\begin{equation} \begin{split} P(SL|z_s) &= \int \frac{\theta_E^2}{4} n_o P(\sigma) \frac{dV_c}{dz_l} dz_l d\sigma \\ P(SL|z_s) &= \int_0^{z_s} \Phi_{SIS}(z_l) dz_l \\ & \text{where, $\Phi_{SIS}(z_l)= \int \frac{\theta_E^2}{4} n_o P(\sigma) \frac{dV_c}{dz_l}d\sigma$ }. \\ & \;\; \Phi_{SIS}(z_l)= \left< \frac{\theta_E^2}{4} n_o \frac{dV_c}{dz_l}\right>_{\sigma\in P(\sigma)} \\ & \text{Note: $\theta_E$ and $\frac{dV_c}{dz_l}$ are functions of $z_l$.} \end{split} \tag{11} \end{equation}

Optical depth for SIE lens

Consider the SIE case with $\sigma$ distribution dependent on $z_l$. The expression for optical depth reads,

\begin{equation} \begin{split} P(SL|z_s) &= \int \frac{\phi_{SIE}}{4\pi} \frac{dN(z_l)}{dV_c d\sigma dq} \frac{dV_c}{dz_l} dz_l d\sigma dq \\ P(SL|z_s) &= \int \frac{\phi_{SIS}}{4\pi} \frac{\phi^{SIE}_{CUT}(q)}{\pi} \big< n \big>_{\sigma, q\in P(\sigma, q)} P(q|\sigma, z_l) P(\sigma|z_l) \frac{dV_c}{dz_l} dz_l d\sigma dq \\ P(SL|z_s) &= \int \frac{\phi_{SIS}}{4\pi} \frac{\phi^{SIE}_{CUT}(q)}{\pi} n_o P(q|\sigma, z_l) P(\sigma|z_l) \frac{dV_c}{dz_l} dz_l d\sigma dq \\ & \;\;\;\;\; \text{$\frac{\phi^{SIE}_{CUT}(q)}{\pi}$ will be found through interpolation.} \end{split} \tag{12a} \end{equation}

\begin{equation} \begin{split} P(SL|z_s) &= \int_0^{z_s} \Phi_{SIE}(z_l) dz_l \\ & \text{where, $\Phi_{SIE}(z_l)= \int \frac{\phi_{SIS}}{4\pi} \frac{\phi^{SIE}_{CUT}(q)}{\pi} n_o P(q|\sigma, z_l) P(\sigma|z_l) \frac{dV_c}{dz_l} d\sigma dq$ }. \\ & \;\; \Phi_{SIE}(z_l)= \left< \frac{\phi_{SIS}}{4\pi} \frac{\phi^{SIE}_{CUT}(q)}{\pi} n_o \frac{dV_c}{dz_l} \right>_{q\in P(q|\sigma,z_l), \sigma\in P(\sigma|z_l)} \\ & \text{If $\sigma$ is independent of $z_l$, then} \\ & \;\; \Phi_{SIE}(z_l)= \left< \frac{\phi_{SIS}}{4\pi} \frac{\phi^{SIE}_{CUT}(q)}{\pi} n_o \frac{dV_c}{dz_l} \right>_{q\in P(q|\sigma), \sigma\in P(\sigma)} \end{split} \tag{12b} \end{equation}

$\phi^{SIE}_{CUT}(q)$ is derive from the expression given in Fei Xu et al. 2022.

Rate equation (statistical formulation)

The final equation of the observed rate of lensed events is shown below. Note that, $z_s$ is sampled from its prior distribution and then rejection sampled with respect to optical depth.

\begin{equation} \begin{split} \mathcal{R}_L &= N_2 \bigg< P(obs|\theta, \theta_L, \beta, z_s, SL) \bigg>_{z_s\in P(z_s|SL), \theta\in P(\theta), \theta_L\in P(\theta_L|z_s, SL), \beta\in P(\beta|\theta_L, z_s, SL)} \end{split} \tag{13} \end{equation}

$P(obs|\theta, \theta_L, \beta, z_s, SL)$ checks whether the event is detectable or not.

Where the sampling priors can be further simplified as follows,

\begin{equation} \begin{split} P(z_s|SL) &= P(SL|z_s)\,P(z_s) \\ P(\theta_L|z_s, SL) &= P(SL|z_s, \theta_L)\, P(\theta_L|z_s) \\ P(\beta|z_s, \theta_L, SL) &= P(SL|z_s, \theta_L, \beta)\, P(\beta|z_s, \theta_L) \end{split} \tag{15} \end{equation}

This allows $z_s$ to be sampled from the astrophysical prior, $P(z_s)$, and then later rejection sample with respect to optical depth, $P(SL|z_s)$. The same is the case for $\theta_L$ ($z_l,\sigma,q$). The strong lensing condition is applied through rejection sampling with respect to $\phi^{SIE}_{CUT}(q)$ ($\propto \theta_E^2\,\phi^{SIE}_{CUT}$). For the source position, $\beta$, it is sampled within the caustic and then checked whether it has 2 or more images or not.

The order of sampling and rate calculation steps in LeR are listed below.

Sample $z_s$ from $\frac{R(z_s)}{N_1}$. And apply rejection sampling with optical depth, $P(SL|z_s)$. Other source parameters are sampled separately, $P(\theta)$.
$z_l$ from $P(z_l|z_s)$.
$\sigma$ together from $P(\sigma|z_l, SL)$.
$q$ from $P(q|\sigma)$.
Calculation of Einstein radius and apply lensing condition to the sampled lens parameters, $P(SL|z_s, z_l, \sigma, q) \propto \theta_E^2\,\phi^{SIE}_{CUT}$.
Other lens parameters ($e_1$, $e_2$, $\gamma_1$, $\gamma_2$, $\gamma$) are sampled independent of the SL condition, $P(e_1,e_2,\gamma_1,\gamma_2,\gamma)$. But, this will be rejection sampled later along with the image position.
Draw image position, $\beta$, from within the caustic boundary and solve the lens equation. Accept it if it results in 2 or more images, otherwise resample $\beta$.
Sometimes (once in 100-200 thousand), the strong lensing condition cannot be satisfied. For these particular events, resample lens parameters and draw image positions, i.e. repeat steps 2-7.
Calculate the magnification, $\mu_i$, time-delay, $dt_i$ and morse phase, $n_i$ for each of the lensed event.
Modify the luminosity distance, $D_l$ to $D_l^{eff}=D_l/\sqrt{|\mu_i|}$, and geocent_time to $t_{eff}=t_{geocent}+dt_i$.
Calculate SNR with gwsnr
Apply the SNR threshold and check whether the event is detectable or not.
Calculate rate of lensed events, $\mathcal{R}_L$ using equation 13.