© The Authors 2023

1. Introduction

Strong gravitational lensing is a powerful tool for studying galaxies and cosmology. Strong lenses have been used to probe the mass structure of massive galaxies (Auger et al. 2010; Oguri et al. 2014; Sonnenfeld et al. 2015; Shajib et al. 2021), to detect substructure (Vegetti et al. 2012; Hezaveh et al. 2016; Nierenberg et al. 2020), to carry out detailed studies of magnified star-forming galaxies (Jones et al. 2013), and to measure the expansion rate of the Universe with time delays (see Treu & Marshall 2016, for a review).

Strong lenses, however, are a biased subset of the general population of galaxies and background sources. In order for a lens-source pair to be included in a strong lensing survey, the following conditions must be met: (1) the lens must produce multiple images of the source, (2) the multiple images must be detected and resolved, and (3) the system must be recognised as a lens by the survey.

Condition 1 poses constraints on the mass distribution of the lens, as well as on the geometry of the system. First of all, not all lenses are capable of producing multiple images, only those whose surface mass density, Σ(θ), is higher than the critical surface mass density for lensing, Σ_cr, at least one position, θ (Schneider et al. 1992). This condition excludes objects with a very diffuse mass distribution from the population of lenses. Second, objects with a higher concentration of mass allow for a larger portion of the source plane to be multiply imaged and are therefore over-represented in lens samples. Third, with fixed lens properties, sources that are located farther away along the line of sight are more likely to be lensed.

Condition 2 sets limits on the observable properties of the sources and, again, the lens mass distribution, in a way that depends on the survey characteristics. On the one hand, brighter sources are more likely to be detected when lensed. On the other hand, lensing magnification allows the detection of fainter objects with respect to the field. However, if the angular separation between the multiple images is too small, the system may not be classified as a lens. Since the image separation increases with lens mass, this means that the lens sample is biased against low-mass objects.

Condition 3 can introduce additional biases, depending on the efficiency of the lens finding process on which the survey is based. Most commonly, lenses with a smaller image separation are more difficult to identify even if detected and resolved, because of contamination from the lens galaxy. This introduces an additional bias towards lenses with a higher or more concentrated mass.

The cumulative effect of these conditions can be summarised as follows. The probability distribution, P_SL, of a sample of strong lenses from a given survey with selection criterion S​ is given by (Sonnenfeld 2022)

(1)

In this equation, ψ_g is the set of parameters describing the properties of foreground galaxies that are relevant for lensing, such as their redshift and mass distribution; ψ_s is the set of parameters describing background sources; P_g and P_s describe the parent distribution of foreground galaxies and background sources in the absence of lensing; and P_sel describes the probability of selecting a lens-source system with parameters ψ_g and ψ_s given the criterion, S​, used to define a lens. This last factor takes both physical and survey selection effects into account, that is, whether a lens with parameters ψ_g produces a strongly lensed image of a source with parameters ψ_s, and whether such an image can be detected and recognised as a strong lens.

The information needed for the left-hand side of Eq. (1) is directly accessible from strong lensing observations. If the main goal of a lensing survey is to characterise the properties of the strong lens population, then it can be accomplished by directly analysing this term. For many applications of strong lensing, however, the aim is to constrain the properties of the general galaxy or source population, P_g and P_s, which are coupled in a non-trivial way via the lens selection probability, P_sel. In order to obtain an unbiased estimate of either P_g or P_s, it is necessary to invert Eq. (1). In principle, this can be done with a Bayesian hierarchical formalism (Sonnenfeld 2022), but knowledge of the lens selection probability, P_sel, is required. This factor can be written as the following product:

(2)

where P_det is the probability of detecting a strong lensing event, and P_find(ψ_g, ψ_s|S, det) is the probability of correctly classifying it as such¹. The detection probability, P_det, can be obtained via simulations. The main challenge is characterising P_find: in most of the existing strong lensing surveys, the process of determining whether a system is included in a strong lens sample is typically a combination of several cuts, usually involving a non-trivial visual selection step.

For the above reasons, the problem of inverting Eq. (1) is a difficult one to tackle exactly. A few studies have attempted to explicitly account for strong lensing selection effects, usually by making ad hoc simplifying assumptions (Sonnenfeld et al. 2015, 2019; Oldham et al. 2017). Whether it is necessary to invert Eq. (1), however, depends on the severity of the strong lensing bias that needs to be corrected and on the accuracy requirements for the key quantities of interest.

In this paper we aim to quantify the strength of the strong lensing bias on a series of foreground galaxy and background source parameters. In particular, we aim to determine how strong lenses differ from the parent population of foreground galaxies and background sources in terms of (a) the radial mass structure of the lenses (i.e. their stellar and dark matter mass density profiles), (b) the ellipticity of the lenses, and (c) the size-magnitude distribution of the lensed sources. This determination depends on (1) how the lens detection probability, P_det, varies as a function of galaxy and source properties, (2) the efficiency of a survey in correctly classifying detected strong lenses (i.e. P_find), and (3) the shape of the galaxy and source parameter distribution, P_g and P_s. To understand this third point we can imagine the limiting case in which both P_g and P_s are Dirac delta functions (i.e. all lenses and sources are identical): in this limit, P_SL simply reduces to the product P_gP_s up to a scaling constant; this corresponds to a case in which there is no lensing bias.

Mandelbaum et al. (2009) carried out a thorough study of point (1): they quantified how the properties of a lens galaxy determine its probability of creating a lensing event with a point source. In this work we revisit the Mandelbaum et al. (2009) study, expanding it to the extended source case. We simulated individual lenses and examined how the lens detection probability varies with lens and source properties. In addition, we address points (2) and (3): we simulated large populations of strong lenses using empirical models, we simulated the lens detection and finding phase, and we quantified the lensing bias under various scenarios. We explored how the results change as a function of the efficiency of a lens survey in discovering small image separation lenses, and as a function of the scatter in mass parameters at fixed light.

Finally, we considered galaxy-quasar lenses. The strong lensing bias of lensed quasars was studied by Collett & Cunnington (2016) in the context of power-law lens mass models with external convergence. Here we address the question of how different galaxy-quasar lenses are from galaxy-galaxy lens samples. This point is relevant for time-delay cosmography studies, in which measurements of the time delay between the multiple images of a strongly lensed quasar are used to constrain the expansion rate of the Universe. Galaxy-galaxy lenses can in principle be used to help break some of the model degeneracies affecting these measurements (Birrer et al. 2020; Birrer & Treu 2021), but any difference between the two lens classes can introduce biases if not corrected for. With this study we aim to quantify this bias.

The structure of this paper is as follows. In Sect. 2 we introduce the basics of gravitational lensing. In Sect. 3 we study individual lens systems. In Sect. 4 we describe our simulations of lens surveys. In Sect. 5 we show the results of our analysis of the simulated lens survey data. We discuss the results in Sect. 6 and draw conclusions in Sect. 7.

The Python code used for the simulation and analysis of the lens sample can be found in a dedicated section of a GitHub repository².

2. Strong lensing basics

The lensing properties of an object with respect to a source depend solely on its dimensionless surface mass density distribution, κ(θ) (also referred to as the convergence). This is the ratio between the surface mass density and the critical surface mass density for lensing:

(3)

The latter quantity is defined as

(4)

where c is the speed of light, G the gravitational constant, and D_d, D_s, and D_ds are the angular diameter distances between the observer and the lens, the observer and the source, and the lens and the source, respectively.

Given a source at angular position β, images of it form at the positions θ that are solutions of the lens equation

(5)

where α is the deflection angle of the lens. This can be expressed in terms of the dimensionless surface mass density by means of the following integral over the whole sky:

(6)

The images of the background source are in general magnified in total flux and in size, while preserving the original surface brightness.

2.1. The axisymmetric lens

In the special case of an axisymmetric lens we can simplify the notation by considering a single coordinate axis with origin at the centre of the lens. We label the components of the image position, source position and deflection angle along this axis as θ, β, and α, respectively. The lens equation for an axisymmetric lens then becomes

(7)

and the expression for the deflection angle reduces to

(8)

This can also be expressed in terms of the total projected mass enclosed within a circle with angular radius equal to θ:

(9)

A very important quantity for describing the strength of a strong lens is the Einstein radius, θ_Ein. For an axisymmetric lens, this is the radius corresponding to the solution of Eq. (7) for a source placed at the same angular position as the lens centre (β = 0). The circle with radius equal to θ_Ein is known as the tangential critical curve. Images that form there have infinite magnification in the tangential direction. It can be shown that θ_Ein satisfies the following condition,

(10)

where is the average surface mass density within a radius equal to θ:

(11)

Axisymmetric lenses of the kind considered in this work typically produce either one or three images of a point source. Figure 1 helps visualise this property. Plotted in Fig. 1 is the quantity θ − α(θ) as a function of the position in the image plane θ, for a few lens models. According to Eq. (7), images form at the locations where this quantity equals the position of the source. Therefore, given a source position β, the number of images and their location can be determined by drawing a horizontal line at the value β on the vertical axis, and finding the points where this line intersects the θ − α(θ) curve.

Fig. 1. Lens equation of axisymmetric lenses. Solid lines represent the right-hand side of Eq. (7) for five lenses with different density profiles. Given a source at position β, its lensed images form at the values of θ where the solid line intersects a horizontal line located at β on the vertical axis. Dotted lines represent positions of the radial caustics. Sources located within the radial caustic of a given lens produce three images. The lens models used in this simulation consist of a stellar component and a dark matter halo, as described in Sect. 3.1. Their density profile is plotted in Fig. 4.

For small values of β (for sources close to the lens centre), the number of images that are produced is three: the source is strongly lensed. For large values of β, instead, only one image is formed. The value of β where the transition occurs is known as the radial caustic, which is marked by the horizontal dotted line in Fig. 1. As can be seen from Fig. 1, a source at this location is mapped to the stationary point of the θ − α(θ) curve. That location on the image plane is known as the radial critical curve. Images that form there have a formally infinite magnification in the radial direction³.

2.2. The elliptical lens

In this paper we focus mostly on lens galaxies with elliptical isodensity contours. Given a surface mass density profile Σ(R), a lens with an elliptical mass distribution can be obtained by replacing the radial coordinate with the circularised radius:

(12)

where x and y are Cartesian axes centred on the lens centre, with x pointing towards the major axis, and where q is the minor-to-major axis ratio.

Figure 2 shows the source-plane caustics of elliptical lenses with different values of the axis ratio. We used the software GLAFIC (Oguri 2021) to obtain the caustic curves. The outermost curves are radial caustics, while the inner ones are tangential ones. The most striking difference with respect to the axisymmetric case (blue curves in Fig. 2) is the fact that the tangential caustic is transformed from a point into a diamond-like curve. Sources located within the diamond produce five images (one of which is usually highly de-magnified). Sources that lie in the region enclosed between the two caustics are imaged three times. Sources outside the radial caustic are imaged only once, as in the axisymmetric case. The fact that the number of images changes by two at a caustic crossing is a general feature of gravitational lenses with no singularities (Schneider et al. 1992).

Fig. 2. Caustics of lenses with fixed radial structure and different ellipticity. Source-plane angular coordinates are in units of the Einstein radius. The outer curve is the radial caustic, while the inner diamond is the tangential caustic. Point sources located outside the caustic are not strongly lensed. Sources that lie in the region enclosed by the two caustics produce three images, while sources inside the tangential caustic produce five images. The lens model adopted for this experiment consists of a stellar component and a dark matter halo, as described in Sect. 3.1. The two components have the same ellipticity.

2.3. Lensing event definition: point sources

In order to compute the probability of a lensing event we must provide an exact definition for it. A necessary condition for a lens-source system to qualify as a strong lens is the presence of multiple images of the source. As we showed above, this requires the source to lie within the region enclosed by the radial caustic. In order to recognise a strong lens in a real survey, however, it is not sufficient for multiple images to exist: they must be detected. For this reason, given the detection limit for a point source, m_lim, we define as strong lensing event any lens-source system with at least two images brighter than m_lim. This is the same definition used by van de Ven et al. (2009), upon which the Mandelbaum et al. (2009) work is based. Labelling with m₂ the magnitude of the second-brightest image, then, in the absence of photometric noise the lens detection probability P_det becomes

(13)

2.4. Lensing event definition: extended sources

Defining a strong lensing event in the case of an extended source is less straightforward. In principle, we should require parts of the source to be multiply imaged. In practice, it is not always easy to determine whether a lensed image contains multiple images or not. This is because, when the source is extended, some of the images can be blended together. In real strong lensing surveys, it is common to find lens candidates in which the lensed source consists of only a single arc. In those cases it is difficult to establish whether the arc is a set of blending images or not, and the decision of including such systems in a strong lens sample is often arbitrary. Here we adopt the following working definition: an extended source is strongly lensed if the number of surface brightness peaks that are detected is larger than its intrinsic number of peaks in the absence of lensing.

We explain how this definition applies in practice with a few examples. For simplicity, we focus on the case of a source with a single surface brightness peak, such as a Sérsic profile. All of the sources considered in this work belong to this family of surface brightness models. We adopt the following procedure to determine the number of detected peaks. Given the observed image of a lensed source, we define its footprint as the ensemble of pixels where the source is detected with S/N > 2. The footprint is in general composed of multiple disconnected regions, corresponding to the different images. In order to only include images that can be clearly identified, we add the condition that the integrated S/N of each separate region must be S/N > 10. This condition has the effect of removing from the source footprint any isolated region consisting only of a very small number of pixels. In a real-world application, it would be very hard to classify such marginal detections as images. If, after applying this cut, the source footprint is spread over multiple separated regions, then the system is classified as a lens. If the source footprint consists instead of a single region, we iteratively increase the surface brightness threshold used to define the footprint and count the number of isolated regions with S/N > 10. The maximum number of detected regions defines the number of detected peaks.

Figure 3 shows a few examples of how this criterion can be used to classify lenses. The first column shows the caustic structure and source position. The second column shows an image of the lensed source. The third column shows the source footprint defined with the procedure described above. Pink pixels correspond to the largest footprint that maximises the number of detected images, while purple pixels belong to the 2σ detection footprint. In the first, second and sixth example, these two footprints coincide. In the third, fourth and fifth example, instead, the 2σ detection footprint consists of a single region, but an increase of the surface brightness threshold leads to the detection of additional images.

Fig. 3. Criterion used to classify lensed images of extended sources. Six examples are shown. First column: caustics (red curves) and source position (blue circle). The radius of the circle corresponds to the radius at which the surface brightness is equal to 2σ the sky background rms fluctuation for a single pixel. In other words, the blue circle delimits the area of the source that can be detected. Second column: mock image of the lensed source, with added noise. Third column: footprint of the source. The purple footprint is obtained with a 2σ detection criterion. The coral region is the largest footprint with the highest number of detected images.

Strictly speaking, our lens definition criterion fails for a perfect Einstein ring. In order to cover such a scenario, we also classify as strongly lensed any sources that produce a footprint with a hole.

We point out that, although our lens definition relies on peak detections, it is not necessary for the peak of the source surface brightness to lie within the caustics in order for the source to be strongly lensed. This can be observed in the fourth example of Fig. 3: although the source centroid lies outside of the caustic (as shown in the first column), the outskirts of the source overlap with the lens, and therefore multiple images are produced.

We can identify three different regimes in strong lensing of extended sources, depending on the relative size of the source and of the lens caustics. In the limit of a very small source, the image configurations that are produced are qualitatively similar to those that can be obtained in the point-source case. When the source and caustic size are comparable, the multiple images tend to be blended into arcs or rings. Most, if not all, of the galaxy-scale strong lenses known belong to these two categories. However, there is a third regime in strong lensing, corresponding to the case in which the source size is much bigger than the caustics, such as in the fifth example of Fig. 3. In this regime, the overall size and total flux of the source are roughly preserved, and the lens produces only a relatively minor perturbation on a localised region of the image. Strong lenses of this kind can be difficult to detect, especially if the region of the image subject to strong lensing overlaps with the light from the foreground galaxy. The ultimate limiting factor to the detection of strong lenses in the large source size regime, however, is the ability to spatially resolve the multiple images. This limit is set by the size of the point spread function (PSF).

All of the examples that we considered were based on sources with a single surface brightness peak. In general this is not necessarily the case, and the intrinsic number of surface brightness peaks of a strongly lensed source is not known a priori. When dealing with a real lens candidate, applying our lens definition criterion requires showing that the observed number of surface brightness peaks can be reproduced with a lens model in which the source has a smaller number of peaks. In samples of lenses defined via visual inspection this process is typically done implicitly, by identifying multiply imaged blobs that belong to the same source element.

To our knowledge, we are the first to propose a surface brightness peak-based definition of a strong lensing event. A popular alternative definition of a strong lens is one based on magnification: only images that are magnified by more than a given threshold are considered as strongly lensed (see e.g. Hilbert et al. 2007). The problem with a definition of this kind is that magnification cannot be determined unambiguously from observations, unless the intrinsic properties of the lensed source are known from other means (e.g. if the source is a standard candle or a standard ruler). Because of the mass-sheet degeneracy (Falco et al. 1985), it is possible to vary the magnification of a lensed image while keeping its observed properties fixed: this could lead to the paradox of two identical-looking lenses that are classified differently on the basis of the underlying magnification. Although we could still use a magnification-based definition for the sake of carrying out our experiments, it would then be difficult to apply our results to real data. Our definition of a strong lensing event, instead, is robust with respect to the mass-sheet degeneracy.

3. Individual lenses

In this section we study how the probability of a strong lensing event varies as a function of lens and source properties. In order to do so, it is useful to introduce the concept of strong lensing cross-section. Given a foreground galaxy with parameters ψ_g, a background source with parameters ψ_s, and a criterion S​ to define a strong lensing event, the strong lensing cross-section is defined as (Sonnenfeld 2022)

(14)

where β is the position of the background source, is the ensemble of source parameter except for the position, and the integral is carried out over the whole sky. The definition above is valid for both a point source and an extended source: although there is no unique way of defining the position of an extended source, the integral over the sky ensures that the result is independent of how the source position is defined. In the limit of low density of background sources, which is satisfied in all practical cases, the probability of a lensing event is proportional to σ_SL. The lensing cross-section defined via Eq. (14) depends solely on the lens detection probability P_det and does not take into account whether the lens can be correctly classified by a lens finder. This separate selection step is captured by the term P_find in Eq. (2). In this section we consider exclusively the detectability of a lens, and therefore focus only on σ_SL. In Sect. 4, when considering specific lens survey simulations, we introduce P_find.

We computed σ_SL in a series of different scenarios of increasing complexity. The model family adopted to describe the radial density profile of the lenses was the same in all of our experiments. We describe this in Sect. 3.1. In Sect. 3.2 we show calculations of the strong lensing cross-section in the case of axisymmetric lenses and point sources. In Sect. 3.3 we generalise the lens geometry to elliptical, while in Sect. 3.4 we replace point sources with extended sources.

3.1. Lens density profile

In this work we focus on massive early-type galaxies as lenses, as these make up the vast majority of known strong lenses. We describe their mass distribution with a model consisting of two concentric components, one describing the baryons and one for the dark matter. We assume that the baryonic component consists entirely of stars, thereby neglecting gas, which is known to contribute very little to the mass of early-type galaxies in the inner regions that are probed by strong lensing. We then assume that the stars follow a Sérsic profile, with projected surface mass density given by

(15)

In the above equation,

(16)

M_* is the total mass, R_e is the radius enclosing half of the total mass, n is the Sérsic index, Γ is the incomplete Gamma function, and b is given by (Ciotti & Bertin 1999)

(17)

Throughout this paper we indicate with θ_e the angular size of the half-light radius.

We fix the Sérsic index of the lenses to n = 4, corresponding to a de Vaucouleurs model. Although early-type galaxies are often described with a free Sérsic index, a de Vaucouleurs profile is able to reproduce their surface brightness distribution to a few percent in the radial range 1 kpc < R < 30 kpc (see e.g. Sonnenfeld 2020), which is the region that is most relevant for strong lensing. Finally, we assume that the light distribution of the stellar component follows its mass distribution exactly. That is, we do not allow for the presence of gradients in the stellar mass-to-light ratio. In Sect. 6.5 we discuss qualitatively what the implications of relaxing this assumption would be.

We describe the dark matter halo with a generalised Navarro, Frenk, and White (gNFW) profile. We first define it by its three-dimensional distribution, which for a spherically symmetric profile is

(18)

The parameter γ_DM is the inner density slope, ρ₀ is a normalisation parameter, while r_s is the scale radius. The logarithmic slope of the density profile transitions from −γ_DM to −3 at a radius r ≈ r_s. The projected surface mass density of a gNFW profile can be expressed in terms of the following integral (Wyithe et al. 2001):

(19)

Figure 4 shows the dimensionless surface mass density profile of Sérsic + gNFW models with various values of the inner dark matter density slope γ_DM and of the fraction of projected dark matter mass within the half-light radius, f_DM. All of the profiles have a dark matter scale radius equal to ten times R_e, and are normalised in such a way that the Einstein radius is equal to the half-mass radius of the stellar component (these assumptions are dropped later). Two main features emerge from Fig. 4. First, the baryons generally dominate the total density in the inner regions (θ < θ_Ein), while the dark matter is the main component at large radii. Second, models with different dark matter fractions and inner dark matter slopes can conspire to produce very similar total density profiles. This is the case, for example, of the (f_DM = 0.5, γ_DM = 1.5) and the (f_DM = 0.2, γ_DM = 2.0) models (green and blue lines in Fig. 4).

Fig. 4. Dimensionless surface mass density profile of Sérsic + gNFW composite models. The Sérsic index of the baryonic component is fixed to n = 4, and the scale radius of the dark matter component is fixed to ten times the half-light radius. All of the profiles are normalised in such a way that the Einstein radius is equal to the half-light radius. Solid lines show the total density profile, dotted lines the dark matter density profile, and dashed lines the baryonic density profile. The dashed blue, orange, and green lines are identical, as they correspond to profiles with the same fraction of baryonic mass within the half-light radius.

3.2. Axisymmetric lenses, point sources

Axisymmetric lenses with a density profile of the kind introduced in Sect. 3.1 can produce either one or three images of a point source. This can be seen in Fig. 1, which shows the lens equation for various values of the dark matter fraction and the inner dark matter slope. All of the lenses shown in this figure have the same Einstein radius, which is equal in size to the half-mass radius of the stellar component. The radius of the radial caustic, marked by the dotted lines in Fig. 1, is a strong function of the lens properties: it is largest in lenses with a smaller dark matter fraction or steeper dark matter slopes. As a result, the source plane area that is subject to strong lensing is an even stronger function of these properties, since it scales with the square of the radial caustic radius. In order to compute the lensing cross-section, however, we must take the magnification into account because it determines whether multiple images can be detected or not.

Figure 5 shows the magnification of the secondary image as a function of source position. The secondary image is the one located in the region between the radial and tangential critical curves, opposite to the source with respect to the lens centre. In most practical cases this is the second brightest image. As Fig. 5 shows, the magnification is very large for sources close to the lens centre (small values of β), decreases with increasing source position and then increases in close proximity to the radial caustic. While for the model with f_DM = 0.8 (purple curve) the magnification is above unity everywhere, other lens models can produce highly de-magnified secondary images for large values of β. Depending on the intrinsic brightness of the lensed source, these images may or may not be detected.

Fig. 5. Magnification of the secondary image as a function of the source position. The lenses are axisymmetric composite models. Their density profile is plotted in Fig. 4. Vertical dotted lines mark the position of the radial caustic.

Using the definition of Eq. (14), we computed the lensing cross-section of a set of axisymmetric lenses, with respect to point sources with different brightnesses. In particular, we considered model lenses with fixed Einstein radius, with angular half-mass radius θ_e fixed to θ_Ein, and varying values of the dark matter fraction and dark matter inner slope. The results are shown in Fig. 6. Each line corresponds to a source with a given intrinsic (i.e. unlensed) magnitude, m_s. The difference between this magnitude and the limiting magnitude of the survey is indicated as follows:

(20)

Fig. 6. Strong lensing cross-section of an axisymmetric lens and a point source. The lens galaxy is a composite model, introduced in Sect. 3.1, with an angular half-light radius equal to the Einstein radius. The left panel shows the effect of varying the slope γDM for fixed fDM; conversely, the right panel shows variations as a function of fDM with γDM fixed. The system is defined as a strong lens if at least two images are detected. Different lines correspond to the difference, Δm, between the source magnitude and the survey detection limit for a point source.

We can see a clear trend between Δm and the lensing cross-section, in both panels of Fig. 6: σ_SL is larger for brighter sources. The trend saturates below a certain Δm, for sufficiently small values of γ_DM or for sufficiently large values of f_DM. In these cases the lensing cross-section coincides with the full area enclosed within the caustic, and further increasing the source brightness does not result in an increased value of σ_SL.

At fixed source brightness, trends with the dark matter inner slope or dark matter fraction are generally weak. The lines of Fig. 6, however, have been computed by keeping the Einstein radius fixed while varying γ_DM or f_DM. This is achieved by adjusting other properties of the lens, such as the total mass of the baryonic or the dark matter component (see Fig. 4). In practice, when varying one ingredient of the lens density profile, the Einstein radius varies in response. To get a complete view of how the lensing cross-section depends on different lens properties, we also computed how σ_SL responds in absolute terms by varying one lens parameter at a time. Figure 7 shows σ_SL as a function of stellar mass, half-light radius, inner dark matter slope, and projected dark matter mass enclosed within an aperture of 5 kpc, M_DM, 5.

Fig. 7. Absolute value of the strong lensing cross-section as a function of various lens properties. The reference lens is a galaxy at z = 0.3 with log M* = 11.5, Re = 7 kpc, γDM = 1.5, log MDM, 5 = 11.0, rs = 100 kpc, and a source at z = 1.5 In each panel, only one property of the lens is varied, as indicated on the label of the horizontal axis. Each curve corresponds to a different value of Δm, in accordance with Fig. 6. The dashed line in each panel shows the quantity .

The lensing cross-section increases with increasing stellar mass and dark matter mass, decreases with increasing R_e for bright sources, while is only a weak function of γ_DM. The lack of a clear trend between σ_SL and γ_DM appears to be in contradiction with the result of Mandelbaum et al. (2009), who found a strong positive correlation between σ_SL and γ_DM. The origin of this discrepancy lies in the different ways in which γ_DM is varied in the two experiments. While we varied γ_DM at fixed M_DM, 5, Mandelbaum et al. (2009) kept fixed the virial mass of the dark matter halo. At fixed virial mass, increasing the inner dark matter slope results in a correspondingly larger dark matter mass in the inner regions, which naturally results in a larger lensing cross-section.

Figures 6 and 7 show that the trends between lens properties and the strong lensing cross-section can be different for sources with different brightnesses. The net effect in a strong lensing survey is an average over the source population, weighted by the source luminosity function. This implies that surveys that target different families of sources, with different luminosity functions (for instance, galaxies or quasars), can have different strong lensing biases. We explore this possibility in Sect. 5.

3.3. Elliptical lenses, point sources

We measured σ_SL for lenses with a fixed radial density profile and different ellipticities, with respect to point sources of different brightnesses. In particular, we set f_DM = 0.5, γ_DM = 1.5, r_s = 10R_e, θ_Ein = θ_e, and set the ellipticities of both the baryonic and dark matter components to be the same, with the same orientation of the major axis. This is the lens model used to produce the caustics plot of Fig. 2. Figure 2 suggests that the size of the source plane area subject to strong lensing, the area enclosed within the outermost caustic, does not vary strongly with the ellipticity of the lens. Therefore, we expect the strong lensing cross-section to be a weak function of ellipticity.

We carried out the computation of σ_SL by means of simulation: we generated a large number of point sources randomly distributed over a given area that includes the caustic, then used GLAFIC to solve the lens equation, find the number of images and their magnification. We then measured the fraction of sources that are strongly lensed according to the criterion of Eq. (13) and multiplied this value by the area over which sources are located. The resulting σ_SL is shown in Fig. 8 as a function of the minor-to-major axis ratio, q.

Fig. 8. Point source strong lensing cross-section as a function of lens axis ratio. Solid lines show the cross-section based on the lens event definition of Eq. (13) and dashed lines the cross-section based on the detection of four images (quad cross-section). Lines of different colours correspond to sources of different intrinsic magnitudes. The dashed blue line, corresponding to the brightest source magnitude, overlaps completely with the dashed orange line. The dashed purple line is zero: very faint sources cannot produce any detectable quads. The parameters of the lens density profile are fDM = 0.5, γDM = 1.5, rs = 10Re, and θe = θEin. The baryonic and dark matter components have the same ellipticity and direction of the major axis.

For bright sources the strong lensing cross-section is approximately constant with axis ratio. This is because, as pointed out earlier, in the bright source regime the cross-section is determined by the area enclosed within the radial caustic, which does not vary much with lens ellipticity. For faint sources we observe a larger variation with q, with a factor of two difference between the largest and smallest value of σ_SL at fixed source brightness.

Most of the sources that result in detectable lenses produce two detectable images. These are sources that are located in the region enclosed between the radial and the tangential caustic. If the source is located within the tangential caustic, however, four detectable images are usually created. Lenses with four visible images, usually referred to as quad lenses, are sometimes given a high priority in certain lensing studies, because they offer more constraints compared to double lenses. For instance, quads make up the majority of the lenses used so far in time-delay studies (Millon et al. 2020). For this reason we also computed an alternative lensing cross-section, in which the definition of a lensing event requires the detection of four images, instead of two. This is plotted in Fig. 8 with dashed lines. The cross-section for quads is a strong function of lens ellipticity, for bright sources. This is a consequence of the fact that the area enclosed within the tangential caustic, which is where a source needs to be in order to produce four images, increases with increased lens ellipticity, as Fig. 2 shows. For sources that are intrinsically fainter than the detection limit, however, the quad cross-section is extremely small, regardless of ellipticity.

3.4. Elliptical lenses, extended sources

In the case of an extended source, the complexity of the problem is increased due to the addition of a series of features: the source surface brightness distribution, with its radial profile, shape and orientation, and the PSF. Moreover, as we discussed in Sect. 2.4, there are different regimes in strong lensing of extended sources, depending on the relative size of the source and the caustics of the lens. For this reason, we split our analysis into two parts. First, we explore the small source size regime, where the source size is comparable to or smaller than the lens caustics. Then, we consider cases in which the source size is bigger than the lens caustics, which we refer to as the large source regime.

3.4.1. Small source sizes

For the sake of reducing the dimensionality of the analysis, we focused on circularly symmetric sources. We also fixed the surface brightness profile to an exponential disk (i.e. a Sérsic model with n = 1). We then took a lens model with the same parameters used in Sect. 3.3 and a minor-to-major axis ratio of q = 0.7. We simulated a large number of images of extended sources with GLAFIC and measured the fraction of them that results in a strong lens, according to the definition of Sect. 2.4. For the small source size experiment, we used pixels with a size equal to 1/20θ_Ein and convolved the images with a Moffat PSF with a full width at half maximum (FWHM) of two pixels and a β parameter of 5.0. We also assumed that the background noise is an uncorrelated Gaussian field. We carried out experiments with sources with different values of the total unlensed flux, f, and half-light radius, θ_e, s. The results are shown in Fig. 9. The total flux of the source, indicated in the legend, is measured in units of the sky background rms fluctuation within an area equal to the square of the Einstein radius.

Fig. 9. Strong lensing cross-section of an extended source: small source size regime (θe, s < θEin). The cross-section is plotted as a function of the ratio between the half-light radius of the source and the Einstein radius of the lens. The lens model is fixed to be the same as in Fig. 8, with axis ratio q = 0.7. The source is a circular exponential profile. Each solid line corresponds to a different value of the total unlensed flux of the source. The flux, f, is expressed in terms of the background noise rms fluctuation measured over an area equal to . The vertical dashed lines correspond to the maximum size for which a galaxy with a given flux can be detected in the absence of lensing. The horizontal dotted line indicates the cross-section for a very bright point source (i.e. the area enclosed within the caustics).

As in the point source case, the lensing cross-section increases with increasing total flux, at fixed source size. At fixed flux, σ_SL stays approximately constant with increasing source size until a given value, then drops rapidly for larger sizes. From a qualitative point of view, this behaviour can be observed also in the absence of lensing: increasing the size of a galaxy while keeping its flux fixed lowers its average surface brightness. If the surface brightness drops below the sky rms fluctuation level, then it becomes very difficult to detect it. In order to determine whether there are lensing-specific features in the σ_SL − θ_e, s relation of Fig. 9, we computed, for each source flux, the maximum half-light radius for which it can be detected in the absence of lensing. We used the same criterion as that of Sect. 2.4 to define a detection: we defined the source footprint as the ensemble of pixels that are 2σ above the background and required the total signal-to-noise ratio within the footprint to be larger than ten. The resulting limiting sizes are shown as vertical lines in Fig. 9. For each given total flux, the non-lensing size limit is similar to the value of θ_e, s at which the lensing cross-section drops.

This result suggests that, to the first approximation, lensing does not introduce a strong selection in source size. While perhaps surprising, this follows from the fact that gravitational lensing preserves surface brightness: a source that can be detected in the absence of lensing will produce images with the same surface brightness when lensed, which can be detected as well. In order to classify a source as strongly lensed, however, we require that multiple images are observed. Only sources that lie within a well-defined region give rise to a strong lensing configuration. If part of the source extends outside of this region, then only a fraction of its flux contributes to creating a set of strongly lensed images. This lowers the signal-to-noise ratio of the multiple images compared to the point source case. The result is that σ_SL starts to decrease with increasing source size at values of θ_e, s that are smaller than the no-lensing detection limit, as observed in Fig. 9.

At the brightest flux explored in the experiment (green line in Fig. 9), the lensing cross-section at small source sizes is larger than the area enclosed within the radial caustic (black dotted line in Fig. 9). This is because, when the source is very bright, it can give rise to multiple images even while its centroid lies outside of the radial caustic, as long as its surface brightness distribution extends into it. This also explains why σ_SL increases with increasing source size, before dropping to zero: the more extended the source, the farther away from the lens centre it can be while still producing multiple images.

3.4.2. Large source sizes

For the large source size case we fixed the surface brightness profile of the source and varied the Einstein radius of the lens. In particular, we set the source half-light radius to 20 pixels and adjusted its total flux in such a way that the 2σ detection footprint in the absence of lensing extends out to the half-light radius. Then, starting from the lens model used in the previous section, we progressively increased the critical surface mass density to reduce the Einstein radius down to values comparable to the pixel size. Figure 10 shows the resulting strong lensing cross-section as a function of the ratio between lens Einstein radius and source half-light radius. For values of the Einstein radius close to the size of the PSF, the lensing cross-section (blue solid line) is very small: this is because the perturbations caused by lensing are not well resolved. For larger values, σ_SL stays approximately constant, around values that are comparable to the source size (horizontal dotted line) and much larger than the area enclosed within the caustics (red dashed line). We conclude that, in the large source size regime, the main factor that determines the lensing cross-section is the area of the background source, provided that the Einstein radius of the lens is larger than the PSF.

Fig. 10. Strong lensing cross-section of an extended source: large source size regime (θe, s ≳ θEin). The cross-section is plotted as a function of the ratio between the Einstein radius of the lens and the half-light radius of the source. The source model is fixed to a circular Sérsic profile with n = 1, detected out to the half-light radius. The source half-light radius measures 20 pixels, and the PSF has an FWHM of two pixels. The lens model is varied: starting from the model of Fig. 9, the critical surface mass density is increased to produce progressively smaller values of θEin. The solid line shows the lensing cross-section, in units of the footprint area of the source in the absence of lensing, as a function of the ratio between the lens Einstein radius and the source half-light radius. The dashed line delimits the area enclosed within the lens caustics, and the dotted line shows the value of the cross-section corresponding to the area of the source. The dash-dotted line is the FWHM of the PSF.

4. Lens population simulations

In the previous section we showed how the lensing cross-section, which is closely related to the lens detection probability P_det, varies as a function of lens and source properties. From here on we focus on the effect that those trends have on populations of lenses. We addressed this question by simulating populations of foreground galaxies and background sources, selecting strong lenses among them, and comparing the properties of the strong lens sample with the parent galaxy population. We simulated a lens-based search (as opposed to a source-based one), in which strongly lensed images are searched among a stellar mass-selected sample of galaxies. Our simulations are based on empirical models, in which existing observations of the baryonic component of galaxies are complemented with a set of assumptions on the mass distribution of the lenses. In Sect. 4.1 we explain how we built our foreground galaxy sample, while in Sect. 4.2 we describe the simulation of the background source population. In Sect. 4.3 we describe how our mock observations of lenses are generated. In Sect. 4.4 we apply a further selection step, based on the angular size of the Einstein radius.

4.1. Foreground galaxies

Our foreground galaxy population consists of a volume-limited sample of early-type galaxies, complete above a minimum observed stellar mass of 10¹¹ M_⊙. We chose to focus on early-type galaxies because many strong lensing surveys have preferentially targeted this class of objects in their lens-finding phase (e.g. Bolton et al. 2006; Gavazzi et al. 2012; Sonnenfeld et al. 2018). This, in turn, had a dual motivation: first, early-type galaxies are among the most massive objects in the Universe, and therefore they are more likely to be lenses; second, their smooth surface brightness distribution and red colour makes it easier to detect arcs from strongly lensed star forming galaxies around them.

We describe lenses with elliptical versions of the two-component model introduced in Sect. 3.1. In the following sections we explain how their parameters are generated.

4.1.1. Stellar mass and redshift distribution

We generated lens galaxies over a finite redshift range, 0.1 < z < 0.7. We chose these lower and upper limits because the value of the critical surface mass density Σ_cs becomes very large outside of this range, and the fraction of galaxies that act as strong lenses drops substantially as a result.

We drew stellar masses from the stellar mass function of quiescent galaxies measured by Muzzin et al. (2013). In particular, we used the following comoving number density distribution:

(21)

We set Φ^* = 1.009 × 10⁻³ Mpc⁻³, α = −0.92 and : these are the best-fit values measured by Muzzin et al. (2013).

We refer to the stellar masses drawn from this distribution as the observed stellar masses, , as opposed to the true stellar masses. The observed stellar mass is meant to represent an estimate of M_* based on stellar population synthesis modelling, which is the method used by Muzzin et al. (2013) to measure the galaxy stellar mass function. Stellar population synthesis measurements, however, are subject to systematic uncertainties, since they have not been calibrated on galaxies whose stellar mass is known by other means. We quantify the discrepancy between the observed and true stellar mass by means of the stellar population synthesis mismatch parameter α_sps, defined as follows:

(22)

The most important source of systematic uncertainty in the measurement of the stellar mass is the assumption of the stellar initial mass function (IMF). Muzzin et al. (2013) assumed a Kroupa IMF (Kroupa 2001) to derive their measurements. This means that, in the absence of other systematic effects, a galaxy with a Kroupa IMF has a value of α_sps = 1.

For each galaxy in the sample, we randomly drew a value of log α_sps from the following distribution:

(23)

where the notation 𝒩(μ, σ²) indicates a Gaussian with mean μ and variance σ². We set the mean of log α_sps to 0.1, as this is an intermediate value among estimates of α_sps from the literature (Conroy & van Dokkum 2012; Cappellari et al. 2013; Smith et al. 2015; Sonnenfeld et al. 2015, 2019). We adopted a few different values for the scatter σ_α. We explain in Sect. 4.1.4 how these were chosen.

4.1.2. Stellar mass density distribution

We described the surface mass density distribution of each galaxy as an elliptical de Vaucouleurs profile (a Sérsic profile with n = 4). Given the observed stellar mass of a galaxy, we assigned a half-mass radius by drawing it from the following distribution in log R_e:

(24)

Then, we assigned an axis ratio q by drawing it from the following beta distribution:

(25)

with α = 6.28 and β = 2.05. The choice for these distributions was motivated by observations of a sample of early-type galaxies. In Appendix we explain how this sample was defined and how the coefficients of Eqs. (24) and (25) were determined.

4.1.3. Dark matter distribution

We modelled the dark matter distribution of each lens galaxy with an elliptical gNFW halo (see Sect. 3.1 for its definition). The parameters of the radial density profile were assigned as follows. Given the stellar mass of a lens, we first determined the value of its halo virial mass, M_h. We defined this as the mass enclosed within a spherical shell with average density equal to 200 times the critical density of the Universe. Using current weak lensing constraints on the halo mass of elliptical galaxies as a reference (Sonnenfeld et al. 2022), we drew halo masses from the following distribution:

(26)

which is a Gaussian in log M_h with a mean that scales with stellar mass and scatter σ_h. The intrinsic scatter in halo mass is not well constrained observationally; therefore, we ran simulations with different values of σ_h, as for the stellar population synthesis mismatch parameter.

To determine the density profile of the halo we relied on a theoretically motivated model that takes into account the effect of baryons on the dark matter. We assumed the dark matter distribution to be initially described by an NFW profile with a concentration⁴ of five, then used the prescription of Cautun et al. (2020) to model the response of the halo to the infall of baryons. In the Cautun et al. (2020) model the halo response is approximated with an analytical function that depends on the present stellar mass distribution, and typically results in a more concentrated and steeper density profile compared to the original NFW model. Finally, we fitted a gNFW profile to the surface mass density of the contracted halo. By doing this, we were able to fully describe the dark matter density profile with three parameters: M_h, γ_DM and r_s. The dark matter density profile defined in this way is determined uniquely by the halo mass, the stellar mass and the half-light radius (the more concentrated the stellar distribution, the stronger the halo response and the steeper the dark matter density profile). In principle we could have allowed for additional degrees of freedom, for instance by relaxing the assumption of a fixed initial halo concentration. In practice, as we explain in Sect. 6.5, our main results are not affected by this choice.

Given the radial profile of the dark matter halo, we obtained an elliptical version of it by applying a transformation of the kind of Eq. (12) to its projected surface mass density. We assumed that the axis ratio and orientation of the halo is the same as that of the stellar component.

4.1.4. Intrinsic scatter

The distribution in stellar mass, halo mass and dark matter inner slope of our sample of simulated foreground galaxies depends on parameters describing the intrinsic scatter in these properties, namely σ_sps and σ_h. Direct observational constraints on these quantities are poor. However, we can derive upper limits on them on the basis of observed scaling relations.

Early-type galaxies lie on the stellar mass fundamental plane, a scaling relation between stellar mass, half-light radius and central velocity dispersion (Hyde & Bernardi 2009; de Graaff et al. 2021):

(27)

The existence of a fundamental plane relation is a consequence of the virial theorem: the velocity dispersion of a galaxy in dynamical equilibrium is directly related to the three-dimensional mass distribution in its inner regions, which is typically dominated by the stellar component. At fixed observed stellar mass distribution, however, the central velocity dispersion can vary depending on the stellar population synthesis mismatch parameter, on the dark matter mass and on the dark matter density profile. This gives rise to a spread in the values of the velocity dispersion given and R_e. Therefore, we can use the observed scatter in velocity dispersion to put an upper limit on the intrinsic scatter in the stellar population synthesis mismatch parameter, halo mass and dark matter slope.

We used Jeans modelling for this purpose. We generated samples of early-type galaxies with the recipes described above and with different values of σ_sps and σ_h. Using the spherical Jeans equation under the assumption of isotropic orbits, we predicted the central velocity dispersion of each galaxy in the sample. Then, we fitted a fundamental plane relation to this mock sample and measured the predicted scatter in velocity dispersion at fixed and R_e. We then varied σ_sps and σ_h to match the observed and the predicted scatter. We did not attempt to match the other parameters of the fundamental plane (i.e. the constant of proportionality of Eq. (27) and the power-law indices β_σ and ξ_σ), because these are sensitive to the orbital anisotropy of the galaxies, which we are asserting to be zero. The details of this procedure are given in Appendix B.

We settled on three different sets of intrinsic scatter parameters, as indicated in Table 1. We label them the fiducial, the low-scatter and the high-scatter scenarios. The high-scatter scenario is ruled out both by the fundamental plane and by weak lensing constraints (Sonnenfeld et al. 2022). Moreover, our dynamical model is quite simplistic, as it neglects the effects of orbital anisotropy and departures from spherical symmetry, which are additional sources of scatter. Nevertheless, we use it in our experiment in order to obtain a more conservative upper limit on the amplitude of the strong lensing bias.

4.2. Background sources

4.2.1. Galaxies

Our background galaxy population is taken from the SURFS-based KiDS-Legacy-Like Simulation (SKiLLS) input catalogue, a hybrid simulation catalogue integrating cosmological simulation with high-quality imaging observations (Li et al. 2023). The cosmological simulation is obtained from the Synthetic UniveRses For Surveys (SURFS) simulations, a set of N-body simulations from Elahi et al. (2018). The galaxy properties, including the star formation history and the metallicity history, are from an open-source semi-analytic model named SHARK⁵ (Lagos et al. 2018). The original photometry is drawn from a stellar population synthesis technique using stellar synthesis libraries with physically motivated dust attenuation and re-emission models (Robotham et al. 2020). Li et al. (2023) further applied an empirical correction to the original synthetic photometry to better agree with the COSMOS2015 observations (Laigle et al. 2016). The galaxy morphology is described by a Sérsic profile with three parameters: the half-light radius in angular units θ_e, s, the Sérsic index n_s, and the axis ratio q_s. These structural parameters are learned from the imaging data obtained with the Advanced Camera for Surveys (ACS) instrument on the Hubble Space Telescope (Griffith et al. 2012). We refer to Li et al. (2023) for details on the learning algorithm and validation.

The complete SKiLLS catalogue contains ∼108 deg² of galaxies with redshift up to 2.5 and r-band apparent magnitude down to 27. We applied a lower limit to the source redshift, by selecting only sources with z_s > 0.8. This ensures that all of the sources lie behind all of the lenses. A similar cut could be applied to a real survey using photometric redshifts, to reduce the incidence of false positives (e.g. arc-like features physically associated with the lens galaxy) in the lens finding phase. The resulting number density of sources is 70 arcmin⁻². We approximated their spatial distribution as uniform in the sky, that is, we neglected clustering of the sources.

4.2.2. Quasars

We described the population of background quasars with the following double power-law luminosity function in the rest-frame UV absolute magnitude, M:

(28)

Following Manti et al. (2017), we set α_Q = −1.35 (faint-end slope), β_Q = −3.23 (bright-end slope), and adopted a redshift-evolving normalisation

(29)

and characteristic magnitude

(30)

Given the redshift and rest-frame UV luminosity of a quasar, we then computed the apparent magnitude in the observed i-band, m_qso, using a quasar spectral template⁶ built from optical and near-infrared spectra obtained by Vanden Berk et al. (2001), Glikman et al. (2006).

For the sake of consistency with the population of background extended sources, we limited the redshift distribution of quasars to the range 0.8 < z_qso < 2.5. We then truncated the distribution in m_qso at two magnitudes fainter than the detection limit (which is specified in Sect. 4.3). Finally, we randomly placed quasars in the source plane with a projected number density of 70 arcmin⁻². This is a much larger number density than observed in the real universe, but we are allowed to do so because we are not interested in predicting the absolute number of strong lenses, so this is a legitimate choice. The advantage of boosting the number density of quasars is that it allows us to produce a large number of lenses without the need for generating too big a population of foreground galaxies.

4.3. Observations

For each of the three intrinsic scatter scenarios, we drew a population of foreground galaxies covering 1000 sq. deg. The expectation value of the number of galaxies in the corresponding volume, given the redshift and stellar mass cuts described in Sect. 4.1.1, is around 300 000. We approximated the foreground galaxies as isolated: when determining whether a galaxy acts as a strong lens, we only modelled the contribution to the lensing signal from the galaxy itself, and neglected that of the environment. We discuss the possible implications of this approximation in Sect. 6. For each lens, we determined its caustics relative to the highest source redshift, z_s = 2.5. We then placed sources randomly behind the lens. If at least one source fell within a circular region enclosing the caustics, we proceeded to compute its lensed images.

For the simulation with extended sources, we produced images with properties similar to those expected for the Euclid Wide survey (Euclid Collaboration 2022). We used a pixel size of 0.1″ and we applied a Moffat PSF with an FWHM of 0.2″ and a β parameter of 5.0. Finally, we assumed a background noise level such that an extended source with half-light radius 0.5″ and an apparent magnitude in the absence of lensing of m_s = 25 is detected with S/N = 10. We then applied the peak detection-based lens selection criterion introduced in Sect. 2.4 to find the strong lenses.

The fiducial scatter simulation with extended background sources produced a sample of 2113 lenses, corresponding to a number density of 2.1 deg⁻². This is about a factor of five smaller than the number density predicted by Collett (2015) for the Euclid survey. This is a result of differences in the description of the source population, in the criteria used to define a strong lensing event, and in the redshift and stellar mass cuts that we applied to define the foreground galaxy population.

For the lensed quasars we did not simulate pixel-level data, but simply computed the observed magnitudes of the multiple images. Following the definition of Sect. 2.3, we included in the sample of strong lenses only systems with at least two images brighter than a limiting magnitude m_lim. We set m_lim = 23.3, which corresponds to the 10σ detection limit of the Legacy Survey of Space and Time (LSST) in a single visit (Oguri & Marshall 2010). This is motivated by the fact that the quasar lenses are meant to simulate a sample assembled for the purpose of carrying out time-delay measurements, which in turn require combining single-visit detections over many epochs. The scenario that we are simulating, then, is that of a lens search in a Euclid-like survey, followed-up with LSST time-domain observations.

The resulting number of quasar lenses in our simulation with the fiducial scatter is 1621. This number is meaningless, given that the simulation was created with an unrealistically large number density of quasars. More interesting is the relative number of quad lenses with respect to the total, which is about 9%. This is a slightly smaller value than the fraction of quads predicted by Oguri & Marshall (2010). The reason for this discrepancy lies in the differences between the lens mass models in the two simulations.

4.4. Lens finding probability

Figure 11 shows the distribution in Einstein radius of the simulated lens samples. In all cases, the distribution peaks at θ_Ein ≈ 0.7″. Although all of the lenses in these samples are detected, it does not necessarily follow that they would all be included in a strong lensing study. There can be a few reasons for excluding certain lenses from a sample. One is the low accuracy of lens finders: current automated lens finding algorithms tend to produce lens candidates samples with low purity (see e.g. Sonnenfeld et al. 2018; Petrillo et al. 2019; Savary et al. 2022). Such samples are then visually inspected, and only those candidates that can be clearly distinguished from false-positives are kept. This visual inspection step tends to disfavour lenses with a small image separation, because of the contamination from the light of the lens galaxy.

Fig. 11. Einstein radius distribution of the simulated lens samples. These are: galaxy-galaxy lenses in the fiducial, high scatter, and low scatter scenarios and galaxy-quasar lenses with fiducial scatter, considering all lenses or only lenses that produce four images.

Another possible reason for refining a sample of lens candidates is the availability of redshift measurements. Redshifts of both the lens and the source are needed in order to convert a lens model into a measurement of mass. When working with large samples of lenses, obtaining spectroscopic measurements is not a viable option, and photometric redshifts are a necessity. Measuring the photometric redshift of a strongly lensed source, however, is challenging, especially when the Einstein radius is small and the source light is blended with the light from the lens (Langeroodi et al. 2023).

Both of these scenarios can result in samples that are incomplete below a certain value of the Einstein radius. We simulate this situation via the following Einstein radius-dependent lens finding probability:

(31)

In words, all lenses with Einstein radii larger than θ_Ein, min are included in the sample, while all those with smaller Einstein radius are excluded. We refer to θ_Ein,min as the completeness limit: our simulated lens samples are complete down to θ_Ein = θ_Ein,min. We explored scenarios with different values of θ_Ein,min. As the next section shows, the larger the minimum Einstein radius, the higher the strong lensing bias.

5. Results

In this section we present the results of the lens population simulations. Section 5.1 shows the results of the galaxy-galaxy lens experiment, while Sect. 5.2 focuses on the population of galaxy-quasar lenses. Given the number of parameters that are needed to describe our model, providing a complete characterisation of the strong lensing bias is a problem with relatively high dimensionality, and is beyond the scope of this paper. For the sake of conciseness, we focus instead on the quantities that we consider most important. Nevertheless, the output of our simulations is available online⁷. We encourage readers who are interested in studying aspects of the strong lensing bias that are not covered in this section to download our data and analyse them directly.

5.1. Galaxy-galaxy lenses

In this section we show the results of the experiments with populations of galaxy-galaxy lenses. We first present the results in a qualitative way, then proceed to quantify the amplitude of the strong lensing bias in various quantities of interest.

Figure 12 shows the distribution in the parameters of the foreground galaxies of the lens systems, compared to those of the parent population, for the fiducial scatter scenario with two different values of the minimum Einstein radius: 0.5″ and 1.0″. A completeness limit of 0.5″ is close to what can currently be achieved via visual inspection of high-resolution space-based images (Garvin et al. 2022), while the 1.0″ limit can be seen as a more conservative case.

Fig. 12. Comparison between the properties of lens samples and the parent population: distribution in foreground galaxy parameters. Filled contours represent the distribution of the parent sample, solid green lines the distribution of the lenses with Einstein radii larger than 0.5″, and solid red lines the distribution of the lenses with Einstein radii larger than 1.0″.

The most striking difference between the samples is in the stellar mass: strong lensing selects preferentially galaxies with larger values of . Lenses tend to also have a larger halo mass, a smaller half-light radius at fixed stellar mass, and a larger stellar population synthesis mismatch parameter. The distribution in ellipticity and inner dark matter slope of the lenses instead look very similar to that of the parent population. Additionally, we can see that the amplitude of the strong lensing bias appears to be always larger in the lens sample with the more restrictive selection on Einstein radius. We quantify the amplitude of these biases later in this section.

Figure 13 shows the distribution in parameters describing the background source population, for the same simulations of Fig. 12. Additionally, Fig. 13 shows the subset of the parent population that consists of detectable sources. These are background galaxies that, in the absence of lensing, can be detected according to the same criterion used for the lensed sources (i.e. the S/N over their 2σ footprint is larger than ten). The detection limit of the survey is at m_s ≈ 25 (the actual limit varies depending on the surface brightness distribution parameters). Because our simulated background source population extends to much fainter magnitudes, the distribution of detectable sources differs substantially from that of the parent population. We can then consider two different strong lensing bias definitions: one that quantifies the difference in lensed source properties with respect to the parent population, and one that describes the difference with respect to the detectable source population. We are mostly interested in the second definition.

Fig. 13. Comparison between the properties of lens samples and the parent population: distribution in background source parameters. Filled contours show the distribution of the parent sample, solid black lines the distribution of the detectable sources, solid green lines the distribution of the lenses with Einstein radii larger than 0.5″, and solid red lines the distribution of the lenses with Einstein radii larger than 1.0″.

Strong lensing tends to preferentially select sources at higher redshift, especially in the more restrictive case with θ_Ein > 1.0″. This is because, at fixed lens properties, increasing the source redshift lowers the critical surface mass density and, consequently, increases the size of the caustics and the Einstein radius. The distribution in the magnitude of the lensed sources is also very different from that of the parent distribution, as it drops rapidly for values larger than m_s ≈ 25. Interestingly, however, there does not seem to be a large difference with respect to the distribution of detectable sources in the absence of lensing. This result can appear to be somewhat counter-intuitive: lensing magnification should allow the detection of sources that are intrinsically fainter than the detection limit. To some extent, this is the case: the distribution of lensed sources shows a slight excess of fainter galaxies compared to the unlensed case. However, as we quantify later in this section, the difference is far from large. The reason for this behaviour lies in the fact that, in the detection of both lensed and unlensed sources, the most important quantity is surface brightness, which is preserved by gravitational lensing when the source is larger than the PSF size. As we showed in Sect. 3.4, the lensing cross-section drops to zero once the surface brightness of the source reaches a value that would make it undetectable in the absence of lensing. For this reason, also the distribution in half-light radius is very similar between the lensed sources and the detectable source population.

Qualitatively, the results shown in Figs. 12 and 13 match our expectations from Sect. 3. In the rest of this section we quantify the lensing bias. We present the results in three different parts. First, we focus on the properties of the lens galaxies that can be observed directly. These are quantities that can be derived from photometry and spectroscopy with minimal assumptions: the lens redshift, the observed stellar mass, the half-light radius and the axis ratio⁸. Any bias in these quantities can be determined relatively easily in a real strong lens survey. In the second part, we consider the parameters related to the mass distribution: the stellar population synthesis mismatch parameter and the dark matter distribution parameters. Determining the bias on these parameters is much more difficult, but these are quantities of great interest from a galaxy science point of view. In the third part, we focus on source distribution parameters.

5.1.1. Bias in observable lens parameters

Figure 14 shows the median redshift, median , median size for a given , and median axis ratio of various lens population simulations, as a function of the minimum Einstein radius. We defined the median size at fixed by fitting the following mass-size relation to the R_e − distribution:

(32)

Fig. 14. Bias on lens observable properties as a function of the minimum Einstein radius. First panel: median redshift. Second panel: median observed stellar mass. Third panel: median half-light radius at an observed mass of log = 11.4. Fourth panel: median axis ratio. In each panel, the dashed black line indicates the value of the parent population. Error bars indicate the standard deviation of the mean of the lens sample.

The quantity shown in the third panel of Fig. 14 is the parameter μ_R, 0, which is the average log R_e at an observed stellar mass of log = 11.4. The horizontal dashed line in each panel shows the value of the parent population: the larger the distance between the curve of a simulation and this line, the higher the strong lensing bias.

As was already visible in Fig. 12, there is a clear bias towards lower redshift, higher stellar mass and smaller sizes, with the bias becoming stronger for more restrictive cuts in Einstein radius. Biases in stellar mass and size are stronger for the simulation with low scatter in α_sps and M_h. This can be explained as follows. For a given value of the observed stellar mass, strong lensing selection favours galaxies with a larger α_sps or M_h, or with a smaller half-light radius. Since the galaxy stellar mass function of the parent population is steep, lenses tend to have a relatively small and large values of α_sps or M_h for their observed stellar mass (this is shown more clearly in the next section). If the intrinsic scatter in α_sps and M_h is low, however, the number of galaxies with a small and a large stellar or halo mass is greatly reduced. Only galaxies with a large or a small size can therefore act as strong lenses. This is an interesting result, because it suggests that it is in principle possible to use the strong lensing bias on or R_e, which is observable, as a way to constrain the amplitude of the intrinsic scatter in the mass parameters, which is poorly known.

5.1.2. Bias in lens mass parameters

Figure 15 shows the median of the distribution in various mass-related quantities, as a function of the minimum Einstein radius. The top panel is the median stellar population synthesis mismatch parameter. In all simulations, strong lenses are biased towards larger values than the parent population. The bias is larger the higher the intrinsic scatter in α_sps and M_h, and increases with increasing θ_Ein, min. For the fiducial model, the bias on α_sps can be as small as 0.03 dex (7%), if no cut on Einstein radius is applied. However, it can rise up to 0.09 dex (23%) when considering only lenses with θ_Ein > 2″. For comparison, we also show two reference values of α_sps, corresponding to a Chabrier IMF (Chabrier 2003) and a Salpeter IMF (Salpeter 1955). These values roughly bracket the current systematic uncertainty on α_sps.

Fig. 15. Bias on lens mass properties as a function of the minimum Einstein radius. First panel: median log αsps. Dotted lines indicate values of αsps corresponding to a Chabrier and a Salpeter IMF (αsps = 1 corresponds to a Kroupa IMF). Second panel: median Mh at log = 11.4. Third panel: median MDM, 5 at log = 11.4 and log Re = 1.2. Fourth panel: median γDM at log = 11.4 and log Re = 1.2. In each panel, the dashed line indicates the value of the parent population. Error bars indicate the standard deviation of the mean of the lens sample.

The second panel of Fig. 15 shows the median halo mass of galaxies with an observed stellar mass of log = 11.4. We measured this quantity by fitting the following relation to the M_h − distribution,

(33)

and taking the resulting value of the parameter μ_h, 0. Similarly to the α_sps case, the strong lenses are biased towards larger values of the halo mass, with the bias being larger for higher-scatter simulations and more restrictive cuts on θ_Ein. In the fiducial scatter scenario with a completeness limit of θ_Ein, min = 1.0″, the halo masses of lenses are on average 0.16 dex larger than those of their parent population.

The third panel of Fig. 15 shows the median projected dark matter mass enclosed within an aperture of 5 kpc, M_DM, 5, at fixed observed stellar mass and half-light radius. We obtained this quantity by fitting the following relation to the M_DM, 5 − − R_e distribution:

(34)

Figure 15 shows the value of μ_DM, 0. The bias on this quantity is qualitatively similar to that on the total halo mass, but much smaller in amplitude. There are two reasons for this. First, the dependence of M_DM, 5 on M_h is shallower than linear. This is because, as the virial mass increases, the virial radius increases as well: the extra mass is spread over a larger volume, and therefore the mass within the inner region does not increase proportionally. Second, while μ_h, 0 is the halo mass at fixed stellar mass, μ_DM, 0 is measured at fixed half-light radius as well. The density profile of the dark matter halos in our simulation have a dependence on galaxy size: the response of dark matter to baryons is stronger for more concentrated stellar distributions (see Sect. 4.1.3). By capturing this dependence in the model of Eq. (34), the residual scatter in M_DM, 5 around the mean is reduced, and so is the strong lensing bias. This, however, is a minor effect: we repeated the analysis while setting the dependence on size to zero and found minimal differences on the derived values of μ_DM, 0. In the fiducial scenario, the bias on M_DM, 5 is as small as 0.02 dex for θ_Ein, min < 0.5″, and is negligible in the low-scatter scenario.

In the fourth panel of Fig. 15 we show the average inner dark matter slope at fixed observed stellar mass and half-light radius, which we measured by fitting the following model

(35)

and taking the resulting value of μ_γ, 0. Different simulations show different trends in the average γ_DM. These are the results of underlying correlations between the dark matter density profile and the stellar and halo mass. In our model, galaxies of a given size with a larger stellar mass have a steeper dark matter slope, because the response of the dark matter to the infall of baryons is stronger. Vice versa, galaxies with a larger halo mass have a shallower density slope. At fixed and R_e, strong lenses have both a larger stellar mass (because their α_sps is larger) and a larger halo mass. In the simulation with low scatter, the bias in stellar mass is more important; therefore, γ_DM has a positive bias. In the simulation with large scatter, the bias in halo mass dominates, and therefore γ_DM is negatively biased.

5.1.3. Bias in source parameters

Figure 16 shows the lensing bias in source-related parameters, with respect to the population of sources that are detectable without lensing. As previously seen in Fig. 13, lensed sources are biased towards higher redshift (top panel of Fig. 16), with the trend being larger for more restrictive cuts on the Einstein radius. The second panel shows the bias in the source magnitude. This bias is negative for small values of θ_Ein, min, that is, lensed sources tend to be brighter than their field counterparts. Naively one might have expected the opposite trend, as strong lensing magnification allows the detection of sources that are intrinsically fainter than the detection limit. However, as the analysis of Sect. 3 clearly shows, the lensing cross-section is always larger for brighter sources, and that explains the sign of the bias. At the same time, the median does not capture the whole picture of the bias in m_s. For instance, in the fiducial simulation with no cut on the Einstein radius, the 90th percentile of the m_s distribution is 26.05, while that of the population of detectable sources is 0.14 mag brighter: indeed, strong lensing allows the detection of fainter sources.

Fig. 16. Bias on source properties as a function of the minimum Einstein radius. First panel: median source redshift. Second panel: median source magnitude. Third panel: median half-light radius in the magnitude bin 24.8 < ms < 25.2. Fourth panel: median Sérsic index in the magnitude bin 24.8 < ms < 25.2. Fifth panel: median axis ratio. In each panel, the dashed line indicates the value of the population of detectable sources. Error bars indicate the standard deviation of the mean of either the full sample (first, second, and fifth panel) or the bin (third and fourth panel).

The third and fourth panels of Fig. 13 show the bias in half-light radius and Sérsic index, respectively, for sources within a magnitude bin centred on m_s = 25 and with a 0.4 mag width. We detected no clear sign of bias. Finally, the fifth panel shows the bias in source axis ratio. Also in this, case no obvious sign of bias was detected, except in the largest values of θ_Ein, min.

5.2. Galaxy-quasar lenses

5.2.1. Bias in lens parameters

The main goal of the experiment with lensed quasars is to check whether there are any differences in the strong lensing bias with respect to the extended source case, at fixed properties of the foreground galaxy population. For this reason, we only ran simulations with lensed quasars in the fiducial scatter scenario and compared the results with those from the extended source simulation. We are also interested in understanding how the subset of quad lenses differs from the entire population of quasar lenses; therefore, we also analysed that subsample on its own.

The first four rows of Fig. 17 show the strong lensing bias as a function of minimum Einstein radius in the following quantities: stellar population synthesis mismatch parameter, halo mass at fixed stellar mass, enclosed dark matter mass at fixed stellar mass and half-light radius, and axis ratio. The biases of the population of lensed quasars (red curves) are very similar to that of lensed galaxies (blue curves), especially for values of θ_Ein, min < 1″, When considering only quad lenses, however, there are some differences, the most remarkable of which is the bias in the axis ratio: quad lenses tend to be on average galaxies with a higher ellipticity. This is a well-known effect (see e.g. Keeton et al. 1997) that can be explained in the context of the results of Sect. 3.2: the higher the ellipticity, the larger the area enclosed within the inner caustic, which is where a source needs to lie in order to produce four or more images. Quad lenses tend to also have a slightly larger halo mass at fixed stellar mass, for θ_Ein, min < 1″.

Fig. 17. Bias in the population of quasar lenses (all and quads only), compared to the extended source simulation. First panel: median log αsps. Second panel: median Mh at log M* = 11.5. Third panel: median MDM, 5 at log M* = 11.5 and log Re = 1.2. Fourth panel: median lens axis ratio. Fifth panel: median mass-sheet transformation parameter, defined in Appendix C. In each panel, the dashed line indicates the value of the parent population. Error bars indicate the standard deviation of the mean of the lens sample.

5.2.2. Impact on H₀

A measurement of the time delay between the images of a lensed quasar can be converted into an estimate of the expansion rate of the Universe, H₀. In general, though, the value of H₀ is degenerate with the mass distribution of the lens (and of matter along the line of sight), which can be difficult to constrain on an individual lens basis. One possible strategy to break this degeneracy is to use information gathered from a larger set of lenses as a prior for the time-delay lenses. This was the strategy adopted by Birrer et al. (2020) in the most recent measurement of H₀ by the TDCOSMO Collaboration: they combined a larger sample of extended source lenses with a smaller set of time-delay lenses. In this section we want to quantify the possible bias on H₀ that can result from this approach.

The lens mass model used by Birrer et al. (2020) consists in a power-law density profile modified by a mass-sheet transformation. The corresponding projected density profile is

(36)

where λ_mst is called the mass-sheet transformation parameter. Strong lensing data can constrain the power-law index γ very well, but not at all λ_mst. Given the time delay, the inferred value of the Hubble constant H₀ scales linearly with 1/λ_mst. Birrer et al. (2020) measured the population distribution of λ_mst on the sample of extended source lenses, by combining strong lensing and stellar kinematics data, and then used it as a prior in the time-delay analysis. Although the mass model used in our work is different from that of Eq. (36), we can estimate the bias associated with this approach by computing the value of λ_mst for each lens and then comparing the strong lensing bias on λ_mst of extended source and quasar lenses. Appendix C explains how we obtained λ_mst for each lens. The strong lensing bias on the median value of λ_mst is shown in the fifth panel of Fig. 17⁹.

At fixed minimum Einstein radius, there is a small difference between the median λ_mst in extended source and quasar lenses, typically of order 1%. The associated bias on H₀ is, therefore, equally small. The value of λ_mst, however, increases with increasing minimum θ_Ein. This implies that if the two samples probe different ranges in Einstein radius, the bias on H₀ can be larger than that. The bias can also be larger if the properties of the general population of galaxies (as opposed to the strong lenses) are used to inform the structure of time-delay lenses (see Collett & Cunnington 2016). This is simply because quasar lenses are on average more similar to galaxy-galaxy lenses than they are to the underlying galaxy population (the blue and red curves in Fig. 17 are closer to each other than they are to the dashed horizontal lines).

6. Discussion

6.1. Key results

In order for a galaxy-source pair to be included as a strong lens in a survey, three conditions must be met. First of all, at least part of the source must be multiply imaged. Second, the multiple images must be detectable. Third, the lens must be recognised as such. Each one of these conditions introduces a bias with respect to the parent population of foreground galaxies and background sources. Together, they define the lens selection probability term P_sel of Eq. (1). The first two points are intrinsic to a strong lensing survey and constitute an unavoidable source of bias. The best case scenario occurs when the efficiency of including a detected strong lens in a survey is always one, and the sample is 100% complete. In this case, the third condition does not introduce any further selection and the strong lensing bias is minimised. We explore this scenario in Sect. 5 when setting the minimum Einstein radius to zero. If only lenses with Einstein radii larger than a given threshold are selected, however, the bias generally increases.

The strong lensing bias affects all quantities that are related to the mass distribution of the lens, as well as the redshifts of lens and source and the source surface brightness parameters. Some of these quantities, such as the observed stellar mass and half-light radius of the lenses, can be directly measured, and it is straightforward to quantify their lensing bias. Other quantities, however, such as the stellar population synthesis mismatch parameter or the dark matter content, are difficult to obtain via traditional, non-lensing, observations. Our simulations are particularly useful to quantify the bias on these properties.

In the fiducial scatter scenario, the one consistent with the observed scatter around the fundamental plane, the bias on α_sps varies from 0.03 dex, when no restrictions on the lens Einstein radius are applied, to 0.09 dex, corresponding to the extreme case in which only lenses with θ_Ein > 2.0″ are selected (see Fig. 15). A reasonable value for the minimum Einstein radius in a space-based survey like Euclid or the Chinese Space Station Telescope is θ_Ein, min = 0.5″. In this case, the bias on α_sps is slightly smaller than 0.04 dex. The current systematic uncertainty on α_sps is 0.2 − 0.3 dex: this is roughly the difference in measurements of the stellar mass of a galaxy obtained with a Chabrier or a Salpeter IMF. Compared with this uncertainty, the amplitude of the strong lensing bias on α_sps is small: strong lensing observations can be used directly to discriminate between these two alternative choices of IMF, without the need to correct for selection effects. Such a goal could be reached with a sample size of a thousand lenses and a statistical study of the kind proposed by Sonnenfeld & Cautun (2021).

Strong lenses are also biased towards larger halo masses. Nevertheless, when focusing on the dark matter content in the inner regions, the amplitude of the strong lensing bias is relatively small, especially when controlling for the stellar distribution of the lens galaxies. For instance, at fixed observed stellar mass and half-light radius, the bias on the projected dark matter mass enclosed within 5 kpc, M_DM, 5, is only a few percent in the fiducial scatter scenario with θ_Ein, min < 1.0″. This means that strong lenses can indeed be used to understand the inner dark matter distribution of galaxies, as long as the dependence of the dark matter distribution on the properties of the baryonic component is accurately modelled (for example, by following the approach of Sonnenfeld & Cautun 2021).

We also looked at the bias on source-related parameters. Strong lensing causes the luminosity function of background sources to be broadened, compared to the distribution of detectable sources. On the one hand, it preferentially selects brighter sources, because the lensing cross-section increases with increasing source brightness. On the other hand, it allows for the detection of sources that are intrinsically fainter than the detection limit in the absence of lensing. Interestingly, we did not find any significant bias on the source size, at fixed magnitude. This result follows from the fact that our lens detection criterion relies on a surface brightness threshold, and surface brightness is preserved by lensing. However, it appears to be in contradiction with the work of Oldham et al. (2017), who argued that their strong lens sample selected preferentially compact sources. The origin of this discrepancy probably lies in the differences between the criteria used to define a strong lens in the two studies. The Oldham et al. (2017) sample was selected primarily via spectroscopy, by looking for signatures of two galaxies at different redshifts in the Sloan Digital Sky Survey (SDSS; York et al. 2000) data. SDSS spectra were taken in fibres with a 1.5″ radius. Lensed galaxies that are comparable in size to this scale, or larger, are less likely to be detected, because part of their flux extends outside of the fibre. Compact galaxies with a large overall magnification, instead, are more likely to be detected (see the discussion in Sect. 5.3 of Oldham et al. 2017). We refer the reader to Arneson et al. (2012) for a thorough study of selection effects associated with spectroscopy.

There is also an indirect way in which strong lensing could preferentially probe more compact sources than possible in the field. Sources that are smaller than the size of the PSF can be easily confused as stars in our own Galaxy. Since these stars cannot be strongly lensed into multiple images, the detection of a strong lens automatically confirms the extragalactic nature of a lensed source. This type of selection effect could explain, for example, the detection of an apparent outlier in the magnitude-size relation of Lyman-break galaxy by means of strong lensing, by Jaelani et al. (2020). It is, however, improper to refer to this effect as strong lensing bias, since the star-galaxy separation is a bias that primarily affects field observations. Our analysis showed that the tendency to select compact sources is not a general feature of strong lens samples.

6.2. Combining lensed quasars with lensed galaxies

Our experiments were also useful for understanding the possible biases that might incur when combining information from samples of galaxy-galaxy lenses with samples of galaxy-quasar lenses. On the one hand the biases in the mass-related quantities (α_sps and M_DM, 5) are very similar in simulations with the same foreground galaxy population and different background sources (see Fig. 17). On the other hand, these biases are a function of the completeness limit. In order to use a sample of strong lenses as a prior for another sample, then, it is important to make sure that (1) the parent population of foreground galaxies among which lenses are searched for is the same in both surveys, and (2) the two surveys can probe the same distribution in Einstein radius.

The above argument applies to samples of quasars selected regardless of the number of multiple images. When dealing with quad lenses, however, the situation is more complicated. First of all, because our experiment revealed differences in the bias on the dark matter distribution between quads and the entire sample of lensed quasars. Second, because quad lenses tend to have a preferentially higher ellipticity. The ellipticity is important in the context of stellar dynamical analyses, which are often used in combination with strong lensing to constrain lens mass parameters (see e.g. Yıldırım et al. 2020). In order to correctly interpret stellar dynamics data, assumptions on the three-dimensional structure of a lens must be made. This is directly related to the projected ellipticity: galaxies with an axis ratio close to one tend to be preferentially elongated in the line-of-sight direction, and vice versa. When using stellar dynamics-based mass measurements to inform the properties of a sample of quad lenses, it is therefore important to take into account possible biases due to the different three-dimensional structure of the two samples.

6.3. The importance of the source population properties

Strictly speaking, the results shown in Sect. 5 apply only to samples of lenses and background sources with the same properties as our simulations. Real lens samples can have more complex source distributions. For example, many existing lens-finding algorithms are tuned to identify sources with different colours from the lenses, introducing a selection that can modify the source distribution in redshift-magnitude space. However, the comparison of Sect. 5.2 shows that, at fixed foreground galaxy population, the strong lensing bias for the sample of extended sources is indistinguishable from that of the lensed quasar population, at least for values of θ_Ein, min < 1.0″. This result suggests that the dependence of the strong lensing bias on the details of the background source population is very weak, as long as no additional selection on the image configuration is applied (e.g. by selecting only quad lenses). This is not a coincidence, but is a consequence of the fact that the dependence of the strong lensing cross-section on the lens parameters is a weak function of source magnitude¹⁰, especially for magnitudes close to the detection limit (see Sect. 3.2). We can then conclude that the details of the properties of the background source population play a secondary role in determining the strong lensing bias. A selection aimed at specific kinds of sources, such as blue star-forming galaxies or quasars, results primarily in different numbers of lenses being discovered, but does not change how biased the lenses are with respect to the general galaxy population.

6.4. Mitigation strategies

As we discussed throughout this paper, the amplitude of the strong lensing bias depends on the intrinsic scatter in the mass parameters of the foreground galaxy population. One way to minimise the bias, therefore, is to identify scaling relations between observable quantities and mass-related properties that can account for part of the scatter. Describing the inner dark matter distribution as a function of stellar mass and half-light radius, as we did in Sect. 5.1, is the first step in this direction: the bias in M_DM, 5 at fixed M_* and R_e is smaller than the overall shift in the median M_DM, 5 of the lens population (albeit by a marginal amount, as we pointed out earlier). This description could be extended by including the central velocity dispersion as an additional control parameter. The velocity dispersion is directly related to the mass distribution of the lens galaxy, which means that, for example, the distribution in α_sps at fixed velocity dispersion should be narrower than its global distribution marginalised over the whole population. The central velocity dispersion, however, is very sensitive to the orbital anisotropy, to the three-dimensional structure, and to gradients in stellar mass-to-light ratio, which are not well known. Therefore, it is difficult to quantitatively estimate the benefit of including it in the description of the lensing bias.

If one wishes to directly account for the strong lensing bias, the formally correct procedure is to explicitly model all of the selection steps in a Bayesian hierarchical formalism, as explained by Sonnenfeld (2022). Although this can be computationally challenging, machine learning can offer an efficient alternative (Legin et al. 2023). In order for either of these approaches to work, however, it is essential that the lens selection procedure can be simulated. This, in turn, requires having an objective definition of a strong lensing event. A peak-based definition such as that introduced in Sect. 2.4 could be used for this purpose.

Nevertheless, it can still be difficult to fully forward model strong lensing selection if visual inspection by humans is involved in the definition of the lens sample. In that case, a possible alternative is to approximate the lens finding probability P_find. For example, we can make the assumption that P_find depends purely on the Einstein radius. This is a reasonable assumption, as long as the lens finding procedure does not selectively pick lenses with different image configurations depending on their Einstein radius. The dependence of P_find on θ_Ein could then be described empirically and inferred during the analysis. This is essentially the approach adopted by Sonnenfeld et al. (2019) in the analysis of strong and weak lensing data from the Hyper Suprime-Cam Survey.

6.5. Limitations of our analysis

The simulations on which our analysis is based are as complex as required by the goal of the analysis itself, which is to estimate the amplitude of the strong lensing bias in a few key quantities. We did, however, make some simplifying assumptions. One such assumption consisted in neglecting the contribution from line-of-sight structure and the environment to the lensing signal, which typically introduce an external shear and convergence. The effect of including external shear in a lens model is similar to that of changing the ellipticity of the lens. External convergence mimics the effect of adding or removing a constant sheet of invisible mass, producing an effect similar to varying the distribution of dark matter. Typical values of the external shear and external convergence are |γ|< 0.1 and |κ_ext|< 0.1 (Millon et al. 2020): these values are small compared to the typical ellipticities and dark matter fractions of our simulated lenses. Therefore, while including them might modify slightly the amplitude of the strong lensing bias in the axis ratio and dark matter distribution parameters, the conclusions of our analysis would not be affected.

We also assumed that the dark matter density profile of the galaxies is completely determined by the halo mass and by the stellar mass distribution (see Sect. 4.1.3). In reality we expect there to be a range of density profiles, for example as a result of halos having a non-zero scatter in their initial (i.e. before baryonic infall) concentration parameter. We could in principle carry out an experiment with such an additional source of variation in the dark matter profile. Adding a scatter in concentration while keeping the halo mass scatter parameter σ_h fixed results in a larger spread in the inner dark matter distribution, which would increase the amplitude of the lensing bias. However, the value of σ_h that we chose for our fiducial experiment was tuned on the basis of the predicted scatter in velocity dispersion around the fundamental plane. That scatter would also be increased by adding flexibility to the dark matter density profile. Then, in order to be consistent with the analysis carried out in this paper, the value of σ_h would need to be lowered accordingly. This, in turn, would reduce the amplitude of the strong lensing bias and introduce a bias on the halo concentration. Since lensing and dynamics are sensitive to the mass distribution at comparable scales (the half-light radius and Einstein radius are similar for most of the lenses), we expect these two effects to cancel out to first approximation. The lensing bias on quantities directly related to the inner mass distribution, such as α_sps and M_DM, 5, would be left roughly unchanged. The main effect of adding scatter to the concentration is to reduce the correlation between the total halo mass and the dark matter mass enclosed within the Einstein radius. This would then reduce the amplitude of the strong lensing bias on the halo mass.

Another simplifying assumption that we made was that of adopting a constant mass-to-light ratio for the stellar component. Massive galaxies are known to have colour gradients, the main effect of which is to cause the stellar half-mass radius to be smaller than the half-light radius (Szomoru et al. 2013; Suess et al. 2019). Since, as shown in Sect. 5, strong lensing bias tends to select galaxies with a more compact stellar distribution, it would also preferentially select galaxies with a steeper mass-to-light ratio gradient. However, testing for the amplitude of the strong lensing bias on the mass-to-light ratio gradient is beyond the scope of this paper.

6.6. Lens versus source selection

The simulations described in Sect. 4, on which our experiment is based, are examples of so-called lens-selected surveys: the lens samples are built by first defining a population of possible lens galaxies and then searching for strong lenses among them. Photometric lens searches are typically carried out in this way (see e.g. Gavazzi et al. 2012; Sonnenfeld et al. 2018; Petrillo et al. 2019). An alternative approach is to build a source-selected sample, by first searching for signatures of strong lensing within a population of possible background sources and then obtaining information on the lens galaxies (see e.g. Hezaveh et al. 2013). In principle, these two approaches lead to different strong lensing biases. For example, from Fig. 12 we can see that the stellar mass distribution of the lenses extends below the minimum mass cut that we applied to the population of foreground galaxies: a source-based search that is able detect those lenses would result in a different sample and a different associated strong lensing bias. However, when considering only lenses with θ_Ein > 1″, lens-selection and source-selection yield the same sample (see the red histogram in Fig. 12): in the limit of large Einstein radius, the two operations (identifying suitable lens galaxies and detecting a strongly lensed source) commute.

In general, there is no such thing as a purely lens- or source-selected sample. A lens-based search is always also source-based, because a detection of the lensed source is required to identify a lens. A source-based search is also lens-based, because the lens properties determine important observational features of the lensed source, such as the number of images, the image separation and the total flux.

6.7. Alternative lens definitions

Our extended source analysis was based on a peak detection criterion to define strong lensing events, as explained in Sect. 2.4. However, many previous studies aimed at simulating populations of lenses adopted a magnification-based definition, in which a lens-source system is classified as a strong lens if the magnification of the source is larger than a minimum value (see e.g. Collett 2015; Robertson et al. 2020). In order to check whether our results are robust to the choice of definition of a strong lensing event, we re-analysed the fiducial sample with this alternative definition. In particular, we defined as strong lenses all systems with a detected background source and μ_tot > 3, where μ_tot is the total magnification (obtained by integrating the flux of all of the lensed images and dividing it by the intrinsic source flux). The total number of strong lenses increased by 34% with this new definition, mostly due to the inclusion of fainter sources with no detected multiple images. Figures 18 and 19 show the strong lensing bias on a series of lens- and source-related properties. Besides the change in the source redshift and magnitude (first two panels in Fig. 19), we see very little differences in the bias with respect to the peak-based definition of the fiducial analysis. We then conclude that the exact definition of a strong lensing event does not have a big impact on the strong lensing bias of lens-related properties. We stress, however, that magnification is not a directly observable quantity, and therefore we advise against applying magnification-based lens definitions to define lens samples, for the reasons discussed in Sect. 2.4.

Fig. 18. Lens strong lensing bias in the case of a magnification-based lens definition (brown curves), compared to that of the fiducial analysis (blue curves). First panel: median log αsps. Second panel: median MDM, 5 at log M* = 11.5 and log Re = 1.2. Third panel: median γDM at log = 11.4 and log Re = 1.2. Fourth panel: median axis ratio. Fifth panel: median mass-sheet transformation parameter, defined in Appendix C. In each panel, the dashed line indicates the value of the parent population.

Fig. 19. Source strong lensing bias in the case of a magnification-based lens definition (brown curves), compared to that of the fiducial analysis (blue curves). First panel: median source redshift. Second panel: median source magnitude. Third panel: median half-light radius in the magnitude bin 24.8 < ms < 25.2. The dashed line indicates the median value of the population of sources detectable in the absence of lensing.

7. Conclusions

Strong lensing is a very active line of research, but a comprehensive understanding of the selection effects associated with it has so far been lacking. This work takes a major step towards filling that knowledge gap. After a thorough investigation, we have learned several lessons regarding the strong lensing bias in photometrically selected lens samples. The following are the most important ones.

1. The strong lensing cross-section increases primarily with increasing lens mass, with decreasing half-mass radius, and with increasing source brightness. At fixed stellar distribution and fixed dark matter mass enclosed within a given aperture, varying the inner dark matter slope has little impact on the strong lensing cross-section (as long as the aperture within which the dark matter mass is normalised is comparable to the Einstein radius).

2. The strong lensing cross-section has little dependence on the size of the source if the former is smaller than the Einstein radius and is detectable in the absence of lensing. Sources with a surface brightness that is too low to be detected are still undetected when strongly lensed.

3. Lens galaxies tend to be more massive and more compact than their non-lens counterparts. Their redshift distribution also differs from that of the general galaxy population. At fixed observed stellar mass (i.e. inferred by means of stellar population synthesis), lens galaxies have a higher intrinsic stellar mass (i.e. a larger stellar population synthesis mismatch parameter, α_sps) and a higher dark matter halo mass. At fixed stellar mass and size, lens galaxies are still biased towards a larger dark matter content.

4. The amplitude of the strong lensing bias depends on how broad the distribution of the parameters describing the lens population is. Important quantities in this regard are the intrinsic scatter in the stellar population synthesis mismatch parameter, σ_sps, and in the dark matter halo mass, σ_h. Increasing the values of σ_sps and σ_h results in a stronger bias on α_sps and on the dark matter mass, and a weaker bias on the observed stellar mass and half-light radius. This implies that, in principle, we could constrain the intrinsic scatter parameters, which are currently poorly known, by measuring the amplitude of the strong lensing bias on and R_e, which is easily observable.

5. The strong lensing bias varies depending on the completeness of the lens sample as a function of the Einstein radius. Surveys that can discover lenses with a smaller Einstein radius have a smaller associated strong lensing bias in all quantities.

6. Under reasonable assumptions regarding the intrinsic scatter parameters, for a Euclid-like survey that is complete down to θ_Ein = 0.5″ the bias on α_sps is smaller than 0.04 dex (10%). This bias is much smaller than the current systematic uncertainty on the stellar population synthesis-based stellar masses. Therefore, strong lensing measurements could be used directly to calibrate stellar mass measurements of massive galaxies to 10% accuracy, without the need to correct for selection effects. Under the same assumptions, the strong lensing bias on the average halo mass at fixed stellar mass is 0.07 dex, while that on the inner dark matter distribution at fixed stellar mass and size is 0.02 dex.

7. Strong lensing selection broadens the magnitude distribution of background sources, compared to the population of objects that are detectable without lensing. At the same time, we did not find any evidence for a bias in the size distribution of background sources at fixed magnitude.

8. Simulations with lensed quasars in place of extended sources showed that the amplitude of the strong lensing bias in the lens-related parameters is not very sensitive to the details of the source population. This result has positive implications for time-delay lensing studies: it means that information from a sample of galaxy-galaxy lenses can be used as a prior on the properties of a set of galaxy-quasar lenses, as long as the two samples probe the same range in Einstein radius and lens observable properties. The associated bias on the inference of H₀ is of the order of 1%.

9. Samples of quad lenses are biased towards galaxies with larger ellipticities, which implies that their three-dimensional structure is also biased. This means that particular care must be taken when stellar dynamics measurements obtained on galaxy-galaxy lenses are used to inform the properties of quads.

In conclusion, strong lensing selection introduces unavoidable biases in the properties of the lens galaxy and background source populations. Biases that affect observable properties, such as the redshift and the light distribution of the lens, can be easily quantified. Biases on mass-related quantities, such as the stellar mass-to-light ratio or the dark matter distribution, are more difficult to measure directly and must be modelled taking selection effects into account. Designing strong lensing surveys with clearly defined and easily modellable selection criteria would help greatly in this task.

Acknowledgments

The collaboration leading to this work was initiated at the 2022 Lorentz Center workshop “Bridging gaps between dynamical probes of galaxies”. A.S. and S.S.L. are supported by NOVA, the Netherlands Research School for Astronomy. Support for this work was provided by NASA through the NASA Hubble Fellowship grant HST-HF2-51492 awarded to AJS by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., for NASA, under contract NAS5-26555.

References

Appendix A: Lens galaxy surface brightness distribution

To assign half-light radii and ellipticities to the simulated lenses we relied on observations of a sample of early-type galaxies selected from the SDSS (York et al. 2000). This sample was selected as follows. Starting from the SDSS spectroscopic sample, we defined a narrow redshift slice around z = 0.2. Then we applied a selection in colour, by choosing objects with g − r > 1.2, and on Sérsic index, by selecting only galaxies with n > 2.5. We used the Sérsic fit measurements by Meert et al. (2015) for this purpose. These cuts produced a sample of 8078 galaxies. We then focused on the r−band de Vaucouleurs model-based photometric measurements of Meert et al. (2015). Using the r−band total flux from the de Vaucouleurs model and the stellar mass-to-light ratio estimates of Mendel et al. (2014), we obtained measurements of . Finally, we fitted the stellar mass-size relation with the following model:

(A.1)

We obtained μ_R, 0 = 1.20, β_R = 0.63 and σ_R = 0.14.

We then proceeded to fit for the axis ratio distribution of the same sample of galaxies. Figure A.1 shows a histogram of the observed distribution. We fitted this with a beta distribution:

(A.2)

Fig. A.1. Distribution in axis ratio of a sample of early-type galaxies. The sample consists of 8078 galaxies from the SDSS, selected by means of cuts in redshift, colour, and Sérsic index, as explained in the text. Measurements of the axis ratio are taken from the de Vaucouleurs model fits of Meert et al. (2015). The model curve is a beta distribution (Eq. A.2) with α = 6.28 and β = 2.05.

We obtained α = 6.28 and β = 2.05.

Appendix B: Upper limits on the intrinsic scatter parameters

We used the fundamental plane of early-type galaxies to set an upper limit on the intrinsic scatter parameters of the simulation: σ_sps, σ_h and σ_γ. First, we measured the fundamental plane of the sample of SDSS early-type galaxies introduced in Appendix A. We took measurements of the line-of-sight stellar velocity dispersion within the SDSS spectroscopic aperture, σ_ap, from the SDSS data release 16 catalogue (Ahumada et al. 2020). Then, we fitted the following model to the distribution of σ_ap of the sample:

(B.1)

We accounted for observational uncertainties on σ_ap when doing the fit; therefore, the parameter σ_σ describes the intrinsic scatter in the logarithm of the velocity dispersion, deconvolved from the observational scatter. In principle we should also account for observational uncertainties on the stellar mass measurement, as they too contribute to the inferred scatter in σ_ap. In practice, however, it is difficult to estimate observational uncertainties on stellar population synthesis measurements (but see Dogruel et al. 2023). For this reason we chose not to propagate uncertainties on . As a result, the inferred scatter parameter σ_σ is slightly overestimated. We obtained μ_σ = 2.36, β_σ = 0.33, ξ_σ = −0.17, and σ_σ = 0.035.

We then generated samples of z = 0.2 galaxies from the model of Sect. 4.1 and used the spherical Jeans equation to predict their central velocity dispersion. The spherical Jeans equation is (Binney & Tremaine 1987)

(B.2)

where ρ_*(r) is the three-dimensional distribution of dynamical tracers (i.e. the stars), σ_r the radial component of the velocity dispersion, β(r) the orbital anisotropy parameter, and M(r) the mass enclosed within a spherical shell of radius r. We assumed isotropic orbits (β = 0), then integrated Eq. B.2 to obtain the seeing-convolved, surface brightness-weighted line-of-sight velocity dispersion within the SDSS spectroscopic aperture¹¹. Finally, we fitted the fundamental plane relation of Eq. B.1 to the observed stellar mass, half-light radius and central velocity dispersion of the mock sample. We repeated this procedure for different values of the intrinsic scatter parameters σ_sps, σ_h and σ_γ, then settled on the three scenarios listed in Table 1.

Appendix C: Mass-sheet transformation parameter

In the time-delay study of Birrer et al. (2020), lenses are described with a mass density profile given by Eq. 36. Birrer et al. (2020) defined the profile by first fitting a pure power-law mass model to the strong lensing data, thus constraining the Einstein radius and power-law index γ, and then applying a mass-sheet transformation of the following kind:

(C.1)

The mass-sheet transformation parameter λ_mst can be interpreted as a quantity describing a departure from a pure power-law profile. One possible way of estimating λ_mst for our lenses is to fit a power-law to the simulated images and then optimise for the value of λ_mst that best matches the true density profile. That, however, would be too time-consuming. Instead, we emulate this process as follows.

Sonnenfeld (2018) showed that, when the main arc and the counter-image of an extended source are well resolved, their relative widths constrain the following combination of derivatives of the lens potential at the Einstein radius:

(C.2)

In other words, the above quantity can be measured robustly (i.e. in a model-independent way). When fitting a pure power-law lens model to lensing data, then, the value of γ that we obtain is the one that reproduces the true ψ″′/(1 − ψ″)¹². This is given by

(C.3)

Hence, given a lens in our sample, we first define a pure power-law model with the same Einstein radius as the true one and with γ given by the equation above. Then, given a pure power-law model, we find λ_mst such that the first three derivatives of the lens potential match the true values. In particular, since a mass-sheet transformation changes the second derivative of the potential as

(C.4)

we set

(C.5)

where is the value of ψ″ of the pure power-law lens model.

All Tables

All Figures

Fig. 1. Lens equation of axisymmetric lenses. Solid lines represent the right-hand side of Eq. (7) for five lenses with different density profiles. Given a source at position β, its lensed images form at the values of θ where the solid line intersects a horizontal line located at β on the vertical axis. Dotted lines represent positions of the radial caustics. Sources located within the radial caustic of a given lens produce three images. The lens models used in this simulation consist of a stellar component and a dark matter halo, as described in Sect. 3.1. Their density profile is plotted in Fig. 4.

In the text

Fig. 2. Caustics of lenses with fixed radial structure and different ellipticity. Source-plane angular coordinates are in units of the Einstein radius. The outer curve is the radial caustic, while the inner diamond is the tangential caustic. Point sources located outside the caustic are not strongly lensed. Sources that lie in the region enclosed by the two caustics produce three images, while sources inside the tangential caustic produce five images. The lens model adopted for this experiment consists of a stellar component and a dark matter halo, as described in Sect. 3.1. The two components have the same ellipticity.

In the text

Fig. 3. Criterion used to classify lensed images of extended sources. Six examples are shown. First column: caustics (red curves) and source position (blue circle). The radius of the circle corresponds to the radius at which the surface brightness is equal to 2σ the sky background rms fluctuation for a single pixel. In other words, the blue circle delimits the area of the source that can be detected. Second column: mock image of the lensed source, with added noise. Third column: footprint of the source. The purple footprint is obtained with a 2σ detection criterion. The coral region is the largest footprint with the highest number of detected images.

In the text

Fig. 4. Dimensionless surface mass density profile of Sérsic + gNFW composite models. The Sérsic index of the baryonic component is fixed to n = 4, and the scale radius of the dark matter component is fixed to ten times the half-light radius. All of the profiles are normalised in such a way that the Einstein radius is equal to the half-light radius. Solid lines show the total density profile, dotted lines the dark matter density profile, and dashed lines the baryonic density profile. The dashed blue, orange, and green lines are identical, as they correspond to profiles with the same fraction of baryonic mass within the half-light radius.

In the text

	Fig. 5. Magnification of the secondary image as a function of the source position. The lenses are axisymmetric composite models. Their density profile is plotted in Fig. 4. Vertical dotted lines mark the position of the radial caustic.
In the text

Fig. 6. Strong lensing cross-section of an axisymmetric lens and a point source. The lens galaxy is a composite model, introduced in Sect. 3.1, with an angular half-light radius equal to the Einstein radius. The left panel shows the effect of varying the slope γDM for fixed fDM; conversely, the right panel shows variations as a function of fDM with γDM fixed. The system is defined as a strong lens if at least two images are detected. Different lines correspond to the difference, Δm, between the source magnitude and the survey detection limit for a point source.

In the text

Fig. 7. Absolute value of the strong lensing cross-section as a function of various lens properties. The reference lens is a galaxy at z = 0.3 with log M* = 11.5, Re = 7 kpc, γDM = 1.5, log MDM, 5 = 11.0, rs = 100 kpc, and a source at z = 1.5 In each panel, only one property of the lens is varied, as indicated on the label of the horizontal axis. Each curve corresponds to a different value of Δm, in accordance with Fig. 6. The dashed line in each panel shows the quantity .

In the text

Fig. 8. Point source strong lensing cross-section as a function of lens axis ratio. Solid lines show the cross-section based on the lens event definition of Eq. (13) and dashed lines the cross-section based on the detection of four images (quad cross-section). Lines of different colours correspond to sources of different intrinsic magnitudes. The dashed blue line, corresponding to the brightest source magnitude, overlaps completely with the dashed orange line. The dashed purple line is zero: very faint sources cannot produce any detectable quads. The parameters of the lens density profile are fDM = 0.5, γDM = 1.5, rs = 10Re, and θe = θEin. The baryonic and dark matter components have the same ellipticity and direction of the major axis.

In the text

Fig. 9. Strong lensing cross-section of an extended source: small source size regime (θe, s < θEin). The cross-section is plotted as a function of the ratio between the half-light radius of the source and the Einstein radius of the lens. The lens model is fixed to be the same as in Fig. 8, with axis ratio q = 0.7. The source is a circular exponential profile. Each solid line corresponds to a different value of the total unlensed flux of the source. The flux, f, is expressed in terms of the background noise rms fluctuation measured over an area equal to . The vertical dashed lines correspond to the maximum size for which a galaxy with a given flux can be detected in the absence of lensing. The horizontal dotted line indicates the cross-section for a very bright point source (i.e. the area enclosed within the caustics).

In the text

Fig. 10. Strong lensing cross-section of an extended source: large source size regime (θe, s ≳ θEin). The cross-section is plotted as a function of the ratio between the Einstein radius of the lens and the half-light radius of the source. The source model is fixed to a circular Sérsic profile with n = 1, detected out to the half-light radius. The source half-light radius measures 20 pixels, and the PSF has an FWHM of two pixels. The lens model is varied: starting from the model of Fig. 9, the critical surface mass density is increased to produce progressively smaller values of θEin. The solid line shows the lensing cross-section, in units of the footprint area of the source in the absence of lensing, as a function of the ratio between the lens Einstein radius and the source half-light radius. The dashed line delimits the area enclosed within the lens caustics, and the dotted line shows the value of the cross-section corresponding to the area of the source. The dash-dotted line is the FWHM of the PSF.

In the text

	Fig. 11. Einstein radius distribution of the simulated lens samples. These are: galaxy-galaxy lenses in the fiducial, high scatter, and low scatter scenarios and galaxy-quasar lenses with fiducial scatter, considering all lenses or only lenses that produce four images.
In the text

	Fig. 12. Comparison between the properties of lens samples and the parent population: distribution in foreground galaxy parameters. Filled contours represent the distribution of the parent sample, solid green lines the distribution of the lenses with Einstein radii larger than 0.5″, and solid red lines the distribution of the lenses with Einstein radii larger than 1.0″.
In the text

Fig. 13. Comparison between the properties of lens samples and the parent population: distribution in background source parameters. Filled contours show the distribution of the parent sample, solid black lines the distribution of the detectable sources, solid green lines the distribution of the lenses with Einstein radii larger than 0.5″, and solid red lines the distribution of the lenses with Einstein radii larger than 1.0″.

In the text

Fig. 14. Bias on lens observable properties as a function of the minimum Einstein radius. First panel: median redshift. Second panel: median observed stellar mass. Third panel: median half-light radius at an observed mass of log = 11.4. Fourth panel: median axis ratio. In each panel, the dashed black line indicates the value of the parent population. Error bars indicate the standard deviation of the mean of the lens sample.

In the text

Fig. 15. Bias on lens mass properties as a function of the minimum Einstein radius. First panel: median log αsps. Dotted lines indicate values of αsps corresponding to a Chabrier and a Salpeter IMF (αsps = 1 corresponds to a Kroupa IMF). Second panel: median Mh at log = 11.4. Third panel: median MDM, 5 at log = 11.4 and log Re = 1.2. Fourth panel: median γDM at log = 11.4 and log Re = 1.2. In each panel, the dashed line indicates the value of the parent population. Error bars indicate the standard deviation of the mean of the lens sample.

In the text

Fig. 16. Bias on source properties as a function of the minimum Einstein radius. First panel: median source redshift. Second panel: median source magnitude. Third panel: median half-light radius in the magnitude bin 24.8 < ms < 25.2. Fourth panel: median Sérsic index in the magnitude bin 24.8 < ms < 25.2. Fifth panel: median axis ratio. In each panel, the dashed line indicates the value of the population of detectable sources. Error bars indicate the standard deviation of the mean of either the full sample (first, second, and fifth panel) or the bin (third and fourth panel).

In the text

Fig. 17. Bias in the population of quasar lenses (all and quads only), compared to the extended source simulation. First panel: median log αsps. Second panel: median Mh at log M* = 11.5. Third panel: median MDM, 5 at log M* = 11.5 and log Re = 1.2. Fourth panel: median lens axis ratio. Fifth panel: median mass-sheet transformation parameter, defined in Appendix C. In each panel, the dashed line indicates the value of the parent population. Error bars indicate the standard deviation of the mean of the lens sample.

In the text

Fig. 18. Lens strong lensing bias in the case of a magnification-based lens definition (brown curves), compared to that of the fiducial analysis (blue curves). First panel: median log αsps. Second panel: median MDM, 5 at log M* = 11.5 and log Re = 1.2. Third panel: median γDM at log = 11.4 and log Re = 1.2. Fourth panel: median axis ratio. Fifth panel: median mass-sheet transformation parameter, defined in Appendix C. In each panel, the dashed line indicates the value of the parent population.

In the text

Fig. 19. Source strong lensing bias in the case of a magnification-based lens definition (brown curves), compared to that of the fiducial analysis (blue curves). First panel: median source redshift. Second panel: median source magnitude. Third panel: median half-light radius in the magnitude bin 24.8 < ms < 25.2. The dashed line indicates the median value of the population of sources detectable in the absence of lensing.

In the text

	Fig. A.1. Distribution in axis ratio of a sample of early-type galaxies. The sample consists of 8078 galaxies from the SDSS, selected by means of cuts in redshift, colour, and Sérsic index, as explained in the text. Measurements of the axis ratio are taken from the de Vaucouleurs model fits of Meert et al. (2015). The model curve is a beta distribution (Eq. A.2) with α = 6.28 and β = 2.05.
In the text