© The Authors 2024

1. Introduction

Strong gravitational lensing is a unique tool for measuring galaxy masses at cosmological distances. Strong lensing constraints have enabled us to gain insight on the mass structure of the most massive galaxies, as well as on cosmology (see Shajib et al. 2022, and references therein). Strong lenses, however, are a special subset of the general galaxy population, for two reasons. First, the probability of a galaxy to act as a strong lens varies as a function of its properties: for example, the more massive and more concentrated a galaxy is, the more likely it is to be a lens. Second, lenses of different properties have different probabilities of being discovered by a strong lensing survey: in general, lenses with a larger Einstein radius (and hence a larger mass, everything else being fixed) are more likely to be identified (see e.g. Herle et al. 2024). The ratio between the probability distribution of galaxy properties that is sampled by a strong lens survey and that of the parent population from which the lenses are drawn is referred to as the strong lensing bias. This bias, if not accounted for, makes it difficult to compare strong lensing measurements to other measurements or models.

One possible approach to mitigate the strong lensing bias is to control for observed galaxy properties when comparing lenses with non-lens galaxies. For example, the stellar mass distribution of strong lens galaxies is very different from that of the rest of the population. However, if a given property X correlates with stellar mass, differences in X between the lens and non-lens sample are reduced when comparing them at a fixed stellar mass.

In practice, even at a fixed stellar mass strong lens galaxies are noticeably different from non-lens galaxies: for example, they are on average more compact (see e.g. Auger et al. 2010). Controlling for central velocity dispersion appears to be a better choice: lens galaxies from the Sloan Lens ACS Survey (SLACS, Bolton et al. 2006) have been shown to be indistinguishable from non-lens galaxies in terms of their fundamental plane¹ (Treu et al. 2006) and environment (Treu et al. 2009), at fixed velocity dispersion. This can be explained by noting that, under the assumption of an isothermal density profile, the probability of a galaxy to create multiple lensed images of a background source is a function of velocity dispersion alone. Thus, to first approximation, strong lens samples can be viewed as being velocity dispersion-selected. However, galaxies are not exactly isothermal systems, and the selection of lenses involves more criteria than the simple production of multiple images (as we discuss below). Therefore we expect that not all of the strong lensing bias can be captured by velocity dispersion alone.

In general, finding analogues to galaxies acting as strong lenses is an ill-posed problem, because the strong lensing bias depends on the mass structure of galaxies, which is difficult to know exactly (and usually is the very property that one seeks to constrain with strong lensing). Instead, a more productive approach to dealing with the strong lensing bias is to forward model the selection process of lenses in a survey. Sonnenfeld (2022) showed an example of how this can be done when working with a complete sample of lenses and well-defined selection criteria. For most of the strong lensing surveys carried out so far, however, it has been very challenging to model the selection function from first principles. This is because existing lens samples are partly the result of human decisions on the classification of lens candidates, which are difficult to reproduce and include in a model.

In this paper we propose a practical solution to the problem of correcting for the strong lensing bias in the analysis of a strong lens sample. We tackled the problem by making some simplifying assumptions on the lens detection probability, and by using empirical models to describe those aspects of the problem that are most difficult to predict. These empirical models are to be inferred from the data along with the mass properties of the population of galaxies. We first explain our methodology in general terms, then show its application to the SLACS strong lens sample.

We used our technique to revisit measurements of the total density profile of massive galaxies. Strong lensing and stellar dynamics observations based on the SLACS sample have revealed how, under the assumption of a power-law density profile, ρ(x)∝x^−γ, the average density slope of early-type galaxies is ⟨γ⟩ = 2.09, corresponding to a distribution slightly steeper than isothermal, and scales positively with stellar surface mass density (Koopmans et al. 2006; Auger et al. 2010; Sonnenfeld et al. 2013). The trend with stellar surface mass density can be easily explained with the dark matter halo contributing to a smaller fraction of the total mass within the half-light radius in more compact galaxies (Shankar et al. 2017). However, when combining SLACS with other lens samples at higher redshift, such as the BOSS Emission-Line Lens Survey (BELLS, Bolton et al. 2012) or the Strong Lens Legacy Survey (SL2S, Gavazzi et al. 2012), lensing and dynamics constraints indicate that the average slope γ increases from z = 1 to the present (Ruff et al. 2011; Bolton et al. 2012; Sonnenfeld et al. 2013), in a way that is difficult to reconcile with theory. Hydrodynamical simulations and empirical models tend to predict opposite trends (Remus et al. 2013; Xu et al. 2017; Shankar et al. 2018), with the exception of HORIZON-AGN (Peirani et al. 2019). The issue of the evolution in γ is related to the late growth history of massive galaxies: these objects are generally believed to grow as a result of dissipationless mergers, which tend to make the density slope shallower (Sonnenfeld et al. 2014). One possible way to reconcile theory and observations is then to allow for some amount of dissipation and related gas infall associated with mergers (Sonnenfeld et al. 2014). As an alternative to the dissipational merger scenario, it has been suggested that strong lensing selection effects are responsible for the observed trend (e.g. Remus et al. 2013). With this work we aim to bring new insight to the problem, by re-analysing the SLACS sample while accounting for the strong lensing bias.

The structure of this paper is as follows. In Section 2 we describe the lens sample that we use for our case study, with a focus on the lens selection process. In Section 3 we introduce a formal definition of the strong lensing bias and describe the general methodology that we adopt to account for it. In Section 4 we describe the model that we fit to the lensing data. In Section 5 we show the results of our analysis. 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. Throughout this work we assumed a flat Λ cold dark matter model, with H₀ = 70 km s⁻¹ Mpc⁻¹, and Ω_m = 0.3.

2. The strong lens sample

Our study focuses on lenses drawn from the SLACS survey. In this section we first describe the selection process of these lenses, and then their associated data.

2.1. The selection process

SLACS lenses were discovered as follows. Foreground galaxies were selected among the main and the luminous red galaxy (LRG) spectroscopic samples (Eisenstein et al. 2001; Strauss et al. 2002) of the Sloan Digital Sky Survey (SDSS, York et al. 2000). The SDSS spectra of the foreground galaxies were scanned for the presence of emission lines associated with objects at higher redshift. Among the systems with a spectroscopic detection of a background galaxy, those expected to have the largest strong lensing cross-section (we explain below how it was estimated) were followed-up with the Hubble Space Telescope (HST) for high-resolution imaging (Bolton et al. 2006). Lenses with a clear detection of a set of strongly lensed arcs were then modelled for the purpose of measuring their Einstein radius.

The resulting selection function of the SLACS sample is rather complex. In order for the source redshift to be measured, multiple emission lines had to be detected. These lines are typically [OIII], Hβ and the [OII] doublet. Since the wavelength coverage of the SDSS spectrograph extends only up to ∼ 9200 Å, only sources up to redshift ∼0.9 could effectively be discovered. At fixed source redshift, the lens detection probability also depends on the source emission line strength and surface brightness distribution. Additionally, the finite size of the SDSS spectroscopic fibre, 1.5″ in radius, put constrains on the Einstein radius and image configuration of the lens. For instance, lenses with a very large Einstein radius have a main arc that falls outside of the fibre, and could only be detected if one of the counter-images is sufficiently close to the centre, or if sufficient flux is brought into the fibre by the point spread function (PSF). We refer to the simulation study of Arneson et al. (2012) for more details on fibre size effects.

Finally, the prioritisation in terms of lensing cross-section applied by Bolton et al. (2006) introduced another non-trivial selection effect. Bolton et al. (2006) estimated the lensing cross-section of each lens candidate as follows. First, they assumed a singular isothermal density profile for each lens. Then, they used the observed stellar velocity dispersion measured by the SDSS, σ_ap^(obs), as a proxy for the velocity dispersion of the singular isothermal mass distribution. Under this assumption, they estimated the Einstein radius of a lens as

$$\begin{array}{c}\hfill {\theta }_{\mathrm{E}}^{(\mathrm{est})}=4\pi {\left(\frac{{\sigma }_{\mathrm{ap}}^{(\mathrm{obs})}}{c}\right)}^2\frac{D_{\mathrm{ds}}}{D_{\mathrm{s}}},\end{array}$$(1)

where c is the speed of light, D_ds is the angular diameter distance between the lens and the source, and D_s is the angular diameter distance between the observer and the source. Finally, they defined the strong lensing cross-section as the area in the source plane that is mapped into sets of multiple images², which for a singular isothermal sphere is equal to πθ_E². We indicate this quantity as σ_SL^(B06). In summary, σ_SL^(B06) is a combination of observed stellar velocity dispersion, lens redshift and source redshift. The way that σ_SL^(B06) enters the SLACS selection function, however, is rather subtle: it was used to rank candidates identified via spectroscopy for the purpose of scheduling follow-up HST observations. This means that the probability of a lens being included in the SLACS sample depended on the size of the parent sample of lens candidates and on the number of systems with HST data: in the limit of infinite HST time, all candidates would have been followed-up. Explicitly incorporating these aspects in a model of the SLACS selection function can be cumbersome. Moreover, a selection based on the observed velocity dispersion, as opposed to its true value (which is not directly accessible), can lead to the velocity dispersion of the sample to be systematically overestimated.

In principle, the selection function of the SLACS sample is also determined by the finite resolution of the photometric data: lenses with an image separation that is close to the size of the PSF may not be recognised as such. In the case of the SLACS sample, however, we argue that resolution does not play a primary role. This is because the smallest Einstein radius in the SLACS sample is θ_E = 0.79″, much larger than the FWHM of HST photometric data. Indeed, the same search methodology at the base of the SLACS selection has been used by Shu et al. (2017) to discover lenses with significantly smaller Einstein radii, down to θ_E ≈ 0.5″. Shu et al. (2017) used a different criterion to prioritise spectroscopically detected candidates for HST follow-up, compared to SLACS. We then concluded that the main source of selection on the Einstein radius distribution of the SLACS sample is the lensing cross-section prioritisation, which depends on the observed velocity dispersion. In Section 3 we explain how we can model the overall selection function of the SLACS sample.

2.2. Sample definition and data

We took 59 lenses from the SLACS survey, for our study. These are the systems for which Auger et al. (2010) carried out a joint lensing and dynamical analysis, measuring their mass density slope. All of the lenses were selected to be early-type galaxies: although the SLACS campaign also led to the discovery of late-type lenses (Treu et al. 2011), these were not included in the Auger et al. (2010) study and are not considered here either. The lenses span the redshift range 0.06 < z < 0.36, with a median redshift of 0.19. For the purposes of our analysis, we needed measurements of the lens redshift, Einstein radius, stellar mass, half-light radius and velocity dispersion, as well as the redshift of the source. We took these data from Auger et al. (2009). Einstein radii were measured by fitting singular isothermal ellipsoid models to HST images of the lenses. Stellar masses were obtained by fitting stellar population synthesis models to broad band photometric measurements. Auger et al. (2009) provided measurements of the stellar mass under the assumption of a Chabrier or a Salpeter stellar initial mass function (IMF). We used Chabrier stellar masses for our analysis (though all of our results can be converted to a Salpeter IMF by rescaling stellar masses upwards by 0.25 dex). Rest-frame V-band half-light radii were obtained by Auger et al. (2009) by fitting a de Vaucouleurs profile to the HST images. Figure 1 (grey dots) show the distribution in redshift, observed stellar mass, half-light radius and stellar velocity dispersion of the SLACS lenses.

Fig. 1. SLACS lenses and parent population. Distribution in redshift, observed stellar mass, half-light radius, and stellar velocity dispersion. The green contours and histograms are obtained by weighting the parent sample distribution by σap4. Vertical dashed or dotted lines indicate the median of the distributions of the corresponding colour. Contours enclose 68% and 95% of the distribution. Half-light radius measurements of the parent population were not used in our analysis. Instead, we assumed a stellar mass-size relation from Hyde & Bernardi (2009), which is shown as a dashed red line. The black line in the redshift-stellar mass panel is the inferred truncation stellar mass mt, introduced in section 4.2.

2.3. The parent population

As we explained in Section 2.1, the first step in the selection of SLACS lenses consisted in scanning SDSS spectra of a large sample of early-type galaxies³. We refer to that sample as the parent population of foreground galaxies, or simply the parent sample. The goals of our are analysis are to (1) understand the strong lensing bias of SLACS with respect to the parent population, and (2) recover the mass structure of the parent population galaxies. To achieve both goals, it is necessary to obtain an accurate description of the parent sample. Obtaining the exact sample that was used originally for the SLACS search is nearly impossible, because Bolton et al. (2006) carried out the search in internal data releases of the SDSS, which may not correspond to any of the public data releases available today. For this reason, we built our fiducial SLACS parent sample on the SDSS data release 18 (DR18, Almeida et al. 2023) instead, under the assumption that it samples the same galaxy population studied by Bolton et al. (2006). Starting from the main and LRG spectroscopic samples, we first applied a redshift cut, limiting the sample to galaxies in the range 0.05 < z < 0.40. Then, in order to select early-type galaxies, we removed objects with Hα emission, by applying an upper limit of 0 Å to the Hα equivalent width (positive values corresponding to emission). Finally, we removed galaxies with blue colour. For this purpose, we used the rest-frame absolute magnitudes provided by the Granada Flexible Stellar Population Synthesis (FSPS, Conroy et al. 2009) value-added catalogue of DR18 (Montero-Dorta et al. 2016), which are K- and E-corrected to z = 0.55, and imposed a minimum u − r value of 2.0. The resulting sample consists of 401576 galaxies.

Given this sample, our inference method required us to obtain a model for their distribution in redshift, stellar mass and half-light radius space. We obtained stellar mass measurements from the Granada FSPS catalogue, using the masses based on the early star-formation with dust model. The Granada FSPS masses for the SLACS lenses are on average 0.23 dex larger than the measurements of Auger et al. (2009), most likely as a result of a different assumption on the stellar IMF. In order to make these measurements consistent, we applied a −0.23 dex shift to the Granada FSPS masses of the entire parent sample.

Finally, we relied on a mass-size relation from the literature to describe the distribution in half-light radius of the parent sample. In particular, we used the quadratic fit of Hyde & Bernardi (2009), which describes the mass-size relation of early-type galaxies over the full stellar mass range of our parent sample. We also took half-light radius measurements directly from the SDSS, using the r-band de Vaucouleurs model from DR18. However, these can suffer from contamination from neighbouring objects and deblending issues (see e.g. Meert et al. 2015), and therefore cannot be used directly for precision measurements without significant additional work. We used these half-light radius measurements for illustrative purposes only.

Figure 1 shows the distribution of the parent sample in redshift, stellar mass, half-light radius and stellar velocity dispersion, along with that of the SLACS lenses. Stellar velocity dispersion measurements were taken from DR18. In Section 1 we argued that, to first approximation, strong lens samples are selected in velocity dispersion. To verify the validity of this assertion, we also plot in Figure 1 the distribution of the parent sample weighted by σ_ap⁴, which is proportional to the estimated lensing cross-section of a singular isothermal sphere with velocity dispersion equal to σ_ap. The medians of the four quantities for the σ_ap⁴-weighted distribution are remarkably close to those of the SLACS lenses. However, the detailed shapes of the distributions are significantly different. This suggests that a σ_ap⁴-weighting can account for large part of the strong lensing bias, but not its entirety.

3. Methods

In this section we develop a formalism to correct for the strong lensing bias in a strong lensing inference. First we examine the problem in its full complexity, then we apply a few approximations to make it more tractable in practice.

3.1. The problem

Sonnenfeld et al. (2023) developed a formalism describing how strong lenses differ from their parent galaxy population. We follow the same notation in this paper. We considered a population of foreground galaxies and a population of background sources. Each foreground galaxy is described by a set of parameters ψ_g, drawn from a probability distribution P_g(ψ_g). These parameters describe all aspects that are relevant for predicting the lensing properties (i.e. the mass structure) and the selection of lens candidates (e.g. the stellar mass or the redshift). Each background source is described by parameters ψ_s, drawn from a probability distribution P_s(ψ_s). These parameters describe all aspects necessary to predict its observed features, when lensed by a foreground object (i.e. its redshift, position, surface brightness distribution and spectrum). Then the joint probability distribution in (ψ_g, ψ_s) of a sample of strong lenses from a given survey is

$$\begin{array}{c}\hfill {\mathrm{P}}_{\mathrm{SL}}({\mathit{\psi }}_{\mathrm{g}},{\mathit{\psi }}_{\mathrm{s}})\propto {\mathrm{P}}_{\mathrm{g}}({\mathit{\psi }}_{\mathrm{g}}){\mathrm{P}}_{\mathrm{s}}({\mathit{\psi }}_{\mathrm{s}}){\mathrm{P}}_{\mathrm{sel}}({\mathit{\psi }}_{\mathrm{g}},{\mathit{\psi }}_{\mathrm{s}}),\end{array}$$(2)

where P_sel(ψ_g, ψ_s) is the probability for a galaxy-source pair to produce a strong lens that is identified as such and included in the survey. The term P_sel, in other words, describes the selection function of a strong lens sample.

Equation (2) describes how the distribution in the properties of lens galaxies and lensed sources (the left-hand side of the equation) differ from that of their parent populations (P_g and P_s on the right-hand side). If we are exclusively interested in the foreground galaxies, we can marginalise over the source distribution to obtain

$$\begin{array}{c}\hfill {\mathrm{P}}_{\mathrm{SL}}({\mathit{\psi }}_{\mathrm{g}})\propto {\mathrm{P}}_{\mathrm{g}}({\mathit{\psi }}_{\mathrm{g}})b_{\mathrm{SL}}({\mathit{\psi }}_{\mathrm{g}}),\end{array}$$(3)

where b_SL is the formal definition of the strong lensing bias and is given by

$$\begin{array}{c}\hfill b_{\mathrm{SL}}({\mathit{\psi }}_{\mathrm{g}})\equiv B{\displaystyle \int d}{\mathit{\psi }}_{\mathrm{s}}{\mathrm{P}}_{\mathrm{s}}({\mathit{\psi }}_{\mathrm{s}}){\mathrm{P}}_{\mathrm{sel}}({\mathit{\psi }}_{\mathrm{g}},{\mathit{\psi }}_{\mathrm{s}}).\end{array}$$(4)

The constant B is set by the condition that P_SL(ψ_g) is normalised to unity. In order to compare a strong lens sample with non-lens galaxies in an accurate way, it is necessary to either invert Equation (3) and obtain P_g, or to multiply the distribution of non-lenses by b_SL, each approach requiring good knowledge of the strong lensing bias. However, b_SL depends on survey-specific details, such as the distribution of background sources (in particular their redshift, as we discuss later) and the lens finding probability. In other words, two different strong lensing surveys targeting the same population of foreground galaxies have in general different strong lensing biases.

We can address the problem by modelling the properties of the foreground population of galaxies, of the background sources, and of the lens finding process. We can assume that all of these aspects of the problem can be fully described by a set of parameters η, and try to infer the values of these parameters from the data of a strong lens survey, d. The posterior probability distribution of the parameters η is

$$\begin{array}{c}\hfill \mathrm{P}(\mathit{\eta }|\mathit{d})\propto \mathrm{P}(\mathit{\eta })\mathrm{P}(\mathit{d}|\mathit{\eta }),\end{array}$$(5)

where P(η) is the prior probability and P(d|η) the likelihood of observing the data given the model. Assuming that measurements on individual lenses are independent of each other, the likelihood can be factorised as

$$\begin{array}{c}\hfill \mathrm{P}(\mathit{d}|\mathit{\eta })={\displaystyle \prod _i^{N_{\mathrm{Lens}}}}\mathrm{P}({\mathit{d}}_i|\mathit{\eta }),\end{array}$$(6)

where N_Lens is the total number of lenses.

In order to evaluate individual factors in Equation (6), it is necessary to introduce parameters describing each system and average over them. That allows us to expand the likelihood of an individual lens as

$$\begin{array}{c}\hfill \mathrm{P}({\mathit{d}}_i|\mathit{\eta })={\displaystyle \int d}{\mathit{\psi }}_{\mathrm{g},i}d{\mathit{\psi }}_{\mathrm{s},i}\mathrm{P}({\mathit{d}}_i|{\mathit{\psi }}_{\mathrm{g},i},{\mathit{\psi }}_{\mathrm{s},i}){\mathrm{P}}_{\mathrm{SL}}({\mathit{\psi }}_{\mathrm{g},i},{\mathit{\psi }}_{\mathrm{s},i}|\mathit{\eta }).\end{array}$$(7)

Assuming that the lens selection process can be forward modelled completely, then it is possible to compute integrals of the kind of Equation (7), and thus evaluate the likelihood. Unfortunately, that is not the case for virtually all existing lens surveys, including SLACS. In the next section we provide a strategy for substantially simplifying Equation (7).

3.2. A solution

The first step consists of simplifying the data vector d_i. The data can be though of as the sum of source- and lens-related observables. The former are the redshift, emission line flux, and surface brightness distribution of the lensed images; the latter are the redshift, stellar mass distribution and velocity dispersion of the lens galaxy. First, we dropped the source emission line flux from d_i. This does not mean that we ignored the spectroscopic detection process: as we explain below, we described it in an approximate way. Rather, we gave up modelling explicitly the spectroscopic selection from first principles. Second, we compressed the information from the source surface brightness into a model-independent summary observable: the Einstein radius. This means that, similarly to the spectroscopic selection, also the photometric detection is not fully modelled, but is described in an approximate way. The data related to a single lens then becomes

$$\begin{array}{c}\hfill {\mathit{d}}_i\equiv \{z_{\mathrm{s}},{\theta }_{\mathrm{E}}^{(\mathrm{obs})},{\mathit{d}}_{\mathrm{g}}\},\end{array}$$(8)

where d_g includes all of the observables related to the lens galaxy. From here onward we omit the subscript i on individual observables and parameters for ease of notation. These assumptions allowed us to simplify the first factor in the integrand of Equation (7), P(d_i|ψ_g, ψ_s), as now the only source-related parameter on which the data depends explicitly is the redshift.

As the next step, we proceeded to simplify the strong lens population distribution term P_SL. We began by expanding P_sel as follows:

$$\begin{array}{c}\hfill {\mathrm{P}}_{\mathrm{sel}}({\mathit{\psi }}_{\mathrm{g}},{\mathit{\psi }}_{\mathrm{s}})={\mathrm{P}}_{\det }({\mathit{\psi }}_{\mathrm{g}},{\mathit{\psi }}_{\mathrm{s}}){\mathrm{P}}_{\mathrm{find}}({\mathit{\psi }}_{\mathrm{g}},{\mathit{\psi }}_{\mathrm{s}}).\end{array}$$(9)

In the above equation, P_det is the probability of a galaxy-source pair to produce a detectable lens in the survey considered. The term P_find is instead the probability of the lens being found in the survey, given that it is detectable. The concept of detectability is rather subtle and deserves further explanation. First of all, we defined a detected strong lens as a system with at least two images of the same source that are detected and spatially resolved in photometric data, and a spectroscopic measurement of the lens and source redshift. Therefore, our definition of detection covers both photometry and spectroscopy. Second, we defined a lens as detectable if the survey data and hypothetical follow-up observations would yield a detection according to our definition. Crucially, P_det is different from the probability of a lens being actually detected: the latter depends on the availability of data, which in turn is determined by human decisions on which candidates are followed-up with high-resolution imaging. These human factors are captured by the factor P_find in Equation (9).

The SLACS sample constitutes only a subset of the strong lenses with a spectroscopic detection in SDSS data. Some of those lenses that are not in SLACS would have been detected photometrically if followed-up with HST (many, in fact, have been detected by Shu et al. 2017). The difference between the pool of detectable lenses and the actual SLACS sample is due to the criteria used for the prioritisation of targets for follow-up, which can be described with P_find.

The use of the concept of detectability might appear like a nuisance, but is made necessary by the fact that the SLACS lens finding process cannot be split cleanly into a detection and a search phase. Equation (9) instead allowed us to separate the problem into a part that depends purely on data, P_det, and a part that describes human decision processes, P_find. In the rest of this section, we provide a practical way to describe the factors entering Equation (9) explicitly.

We began by writing the source parameters as

$$\begin{array}{c}\hfill {\mathit{\psi }}_{\mathrm{s}}=(z_{\mathrm{s}},{\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}},\mathit{\beta }),\end{array}$$(10)

where z_s is the redshift, β the angular position in the source plane and ${\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}}$ encodes all information on the light distribution (flux, surface brightness profile and spectrum). We then simplified the probability of finding a detectable lens, P_find, by making the following approximation:

$$\begin{array}{c}\hfill {\mathrm{P}}_{\mathrm{find}}({\mathit{\psi }}_{\mathrm{g}},{\mathit{\psi }}_{\mathrm{s}})\approx {\mathrm{P}}_{\mathrm{find}}({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}}).\end{array}$$(11)

That is, we assumed that the lens finding probability does not depend on any source property other than the redshift. In the case of SLACS, this is a very good approximation: the decision process that led to the follow-up and subsequent discovery of the SLACS lenses was based purely on the measured values of the lens redshift, lens velocity dispersion, and source redshift (source flux and spectrum only enter P_det).

With the approximations introduced so far, Equation (7) could be rewritten as

$$\begin{array}{cc}\hfill \mathrm{P}({\mathit{d}}_i|\mathit{\eta })={\displaystyle \int }& d{\mathit{\psi }}_{\mathrm{g}}d{\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}}\mathrm{P}({\mathit{d}}_{\mathrm{g}}|{\mathit{\psi }}_{\mathrm{g}})\mathrm{P}({\theta }_{\mathrm{E}}^{(\mathrm{obs})}|{\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}}){\mathrm{P}}_{\mathrm{g}}({\mathit{\psi }}_{\mathrm{g}}|\mathit{\eta })\times \hfill \\ \hfill & {\mathrm{P}}_{\mathrm{s}}(z_{\mathrm{s}},{\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}}|\mathit{\eta }){\mathrm{P}}_{\mathrm{find}}({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}}|\mathit{\eta }){\sigma }_{\mathrm{SL}}({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}},{\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}}).\hfill \end{array}$$(12)

In writing Equation (12) we used Equation (9), we introduced the strong lensing cross-section,

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{SL}}({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}},{\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}})\equiv {\displaystyle \int d}\mathit{\beta }{\mathrm{P}}_{\det }({\mathit{\psi }}_{\mathrm{g}},{\mathit{\psi }}_{\mathrm{s}}),\end{array}$$(13)

we carried out the integral over the source redshift (because it is known spectroscopically), and we implicitly assumed that the sources are distributed uniformly in the sky.

To further simplify the problem, we assumed that the strong lensing cross-section can be separated into a geometric factor g that depends purely on the lens properties and source redshift, and a factor l that depends only on the source light distribution:

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{SL}}({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}},{\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}})\propto l({\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}})g({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}}).\end{array}$$(14)

The geometric factor g can be thought of as the strong lensing cross-section computed for a reference value of the source light distribution. The ansatz of Equation (14) is, in general, not true. We discuss the limitations of this approximation in Section 4.4.

Nevertheless, assuming Equations (11) and (14) to hold, we could compute the integral over the source light parameters in Equation (12) to obtain

$$\begin{array}{cc}\hfill \mathrm{P}({\mathit{d}}_i|\mathit{\eta })={\displaystyle \int }& d{\mathit{\psi }}_{\mathrm{g}}\mathrm{P}({\mathit{d}}_{\mathrm{g}}|{\mathit{\psi }}_{\mathrm{g}})\mathrm{P}({\theta }_{\mathrm{E}}^{(\mathrm{obs})}|{\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}}){\mathrm{P}}_{\mathrm{g}}({\mathit{\psi }}_{\mathrm{g}}|\mathit{\eta })\times \hfill \\ \hfill & {\mathrm{P}}_{\mathrm{s}}^{(\mathrm{eff})}(z_{\mathrm{s}}|\mathit{\eta }){\mathrm{P}}_{\mathrm{find}}({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}}|\mathit{\eta })g({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}}),\hfill \end{array}$$(15)

where

$$\begin{array}{c}\hfill {\mathrm{P}}_{\mathrm{s}}^{(\mathrm{eff})}(z_{\mathrm{s}}|\mathit{\eta })={\displaystyle \int d}{\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}}{\mathrm{P}}_{\mathrm{s}}(z_{\mathrm{s}},{\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}})l({\mathit{\psi }}_{\mathrm{s}}^{\mathrm{l}})\end{array}$$(16)

is an effective probability distribution of the redshift of sources that can be strongly lensed. One subtlety with this last definition is that P_s^(eff)(z_s) is neither the redshift distribution of the background source population, which describes all sources regardless of whether they are detectable or not, nor the redshift distribution of the lensed sources. So it cannot be trivially inferred from the data without taking selection effects into account.

The factors that do not depend on the data in the integrand of Equation (15) can be written in compact form, to indicate the effective distribution of strong lens parameters (averaged over the source light distribution):

$$\begin{array}{c}\hfill {\mathrm{P}}_{\mathrm{SL}}^{(\mathrm{eff})}({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}})={\mathrm{P}}_{\mathrm{g}}({\mathit{\psi }}_{\mathrm{g}}|\mathit{\eta }){\mathrm{P}}_{\mathrm{s}}^{(\mathrm{eff})}(z_{\mathrm{s}}|\mathit{\eta }){\mathrm{P}}_{\mathrm{find}}({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}}|\mathit{\eta })g({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}}).\end{array}$$(17)

An important point is that the above probability must be properly normalised:

$$\begin{array}{c}\hfill {\displaystyle \int \mathrm{d}}{\mathit{\psi }}_{\mathrm{g}}\mathrm{d}{z}_{\mathrm{s}}{\mathrm{P}}_{\mathrm{SL}}^{(\mathrm{eff})}({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}})=1.\end{array}$$(18)

The normalisation can be absorbed by the geometric cross-section factor g (hence the proportionality sign in Equation (14)). With these approximations the strong lensing bias was also simplified to

$$\begin{array}{c}\hfill b_{\mathrm{SL}}({\mathit{\psi }}_{\mathrm{g}})\propto {\displaystyle \int \mathrm{d}}{z}_{\mathrm{s}}{\mathrm{P}}_{\mathrm{s}}^{(\mathrm{eff})}(z_{\mathrm{s}}){\mathrm{P}}_{\mathrm{find}}({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}})g({\mathit{\psi }}_{\mathrm{g}},z_{\mathrm{s}}).\end{array}$$(19)

In summary, we have approximated the expression for the lens finding probability and we have eliminated the explicit dependence of the model on the source light distribution. The price to pay for this choice is losing sensitivity to the total number of lenses: had we known the number of spectroscopic detections and photometric non-detections from the SLACS search, together with the distribution of sources in redshift, surface brightness and emission line flux space, we could have used that information to further constrain the model (in a similar vein to the approach of Zhou & Sonnenfeld 2024). We chose not to treat that aspect of the problem explicitly, for the sake of simplicity.

To compute the likelihood, one still has to evaluate the integrand of Equation (15). The geometric factor g can be computed by means of simulation. The factors P_s^(eff)(z_s) and P_find can be inferred directly from the data. In Section 4 we show a practical example of how this can be done.

4. The model

4.1. Mass density profile and individual lens parameters

Our goal was to infer the mass density profile of the population of massive galaxies, using strong lensing and stellar kinematics constraints. We did so under the assumption of a spherical power-law density profile. Although the true density profile of galaxies might not be a perfect power law, deviations larger than ∼10% are inconsistent with stellar kinematics and strong lensing constraints (Cappellari et al. 2015; Tan et al. 2024).

A power-law profile has two degrees of freedom in the radial direction. We chose to parameterise it in terms of the total mass within a cylinder of radius 5 kpc, M₅, and the three-dimensional density slope γ. The choice of M₅ was motivated by the fact that lensing directly constraints the projected mass, and the value of 5 kpc is very close to the median Einstein radius of the SLACS lenses (which is 4.2 kpc). Therefore we expected it to be relatively well constrained by the data, without the need for major extrapolations. Expressed in terms of M₅ and γ, the three-dimensional density profile of a galaxy reads

$$\begin{array}{c}\hfill \rho (x)=\frac{3-\gamma }{2{\pi }^{3/2}}\frac{\mathrm{\Gamma }[\gamma /2]}{\mathrm{\Gamma }[(\gamma -1)/2]}\frac{M_5}{{(5\mathrm{kpc})}^{3-\gamma }}{x}^{-\gamma },\end{array}$$(20)

where x is the three-dimensional distance from the lens centre and Γ is the Gamma function.

In addition to the two parameters of the power-law, we described each foreground galaxy by means of its redshift z_g, its stellar mass M_*, and its half-light radius R_e (for simplicity also referred to as size). We needed the redshift, first of all, to compute the lensing properties. We also needed redshift and size (more generally, the stellar distribution) to predict the stellar velocity dispersion, to be compared with stellar kinematics measurements. Finally, as we explain in the next section, our model for the population distribution allows for a scaling of M₅ and γ with M_* and R_e. We used the stellar mass purely as a label: in our model, the value of the stellar mass of a galaxy does not directly translate into a mass density profile. In summary, the set of individual galaxy parameters ψ_g is

$$\begin{array}{c}\hfill {\mathit{\psi }}_{\mathrm{g}}\equiv \{z,m_{\ast },r_{\mathrm{e}},m_5,\gamma \},\end{array}$$(21)

where

$$\begin{array}{cc}\hfill m_{\ast }& \equiv log\frac{M_{\ast }}{M_{\odot }},\hfill \end{array}$$(22)

$$\begin{array}{cc}\hfill r_{\mathrm{e}}& \equiv log\frac{R_{\mathrm{e}}}{\mathrm{kpc}},\hfill \end{array}$$(23)

$$\begin{array}{cc}\hfill m_5& \equiv log\frac{M_5}{M_{\odot }}.\hfill \end{array}$$(24)

4.2. Foreground galaxies distribution

In this section we introduce the probability distribution P_g describing the population of foreground galaxies.

We began by writing P_g(ψ_g) as follows:

$$\begin{array}{cc}\hfill {\mathrm{P}}_{\mathrm{g}}({\mathit{\psi }}_{\mathrm{g}})=& \mathcal{M}(z,m_{\ast })\mathcal{R}(r_{\mathrm{e}}|m_{\ast })\mathcal{P}(m_5,\gamma |m_{\ast },r_{\mathrm{e}}).\hfill \end{array}$$(25)

The first factor on the right-hand side of the above equation is the distribution in redshift and stellar mass of the parent sample. This sample is the sum of galaxies from the main and the LRG spectroscopic samples, and modelling the corresponding distribution from first principles is challenging. Nevertheless we expected, at any given redshift, the distribution to be volume-limited at the high-mass end, and to drop to zero at lower masses. Therefore we assumed the following functional form:

$$\begin{array}{c}\hfill \mathcal{M}(z,m_{\ast })\propto {\left({10}^{m_{\ast }-\overline{m}}\right)}^{\alpha +1}exp\{-{10}^{m_{\ast }-\overline{m}}\}\frac{\mathrm{d}V}{\mathrm{d}z}{f}_t(z,m_{\ast }),\end{array}$$(26)

where dV/dz is the derivative of the comoving volume with respect to the redshift, and f_t is a function that truncates the distribution at low values of m_*. We chose the following form for f_t:

$$\begin{array}{c}\hfill f_t(z,m_{\ast })=\frac{1}{\pi }arctan\left[\frac{m_{\ast }-m_t(z)}{{\sigma }_t}\right]+\frac{1}{2},\end{array}$$(27)

where m_t is a redshift-dependent truncation mass, which we modelled with a fifth-order polynomial function. Because the strong lenses occupy the high-mass end of the stellar mass distribution, the exact details of the stellar mass distribution at low M_* are not critical for the accuracy of the model.

The factor ℛ in Equation (25) describes the distribution in half-light radius as a function of redshift and stellar mass. Following Sonnenfeld et al. (2019), we modelled it as

$$\begin{array}{c}\hfill \mathcal{R}(r_{\mathrm{e}})=\mathcal{N}({\mu }_R(m_{\ast }),{\sigma }_R^2),\end{array}$$(28)

where the notation 𝒩(μ, σ²) indicates a Gaussian distribution with mean μ and variance σ². We set the mean μ_R to scale quadratically with log M_* as

$$\begin{array}{c}\hfill {\mu }_R={\mu }_{R,0}+{\mu }_{R,1}{m}_{\ast }+{\mu }_{R,2}{m}_{\ast }^2.\end{array}$$(29)

We fixed the values of the coefficients in Equation (29) to those measured by Hyde & Bernardi (2009): μ_R, 0 = 7.55, μ_R, 1 = −1.84, μ_R, 2 = 0.11. Then, we set the intrinsic scatter in size to σ_R = 0.112, which is the value measured by Sonnenfeld et al. (2019) on a large sample of elliptical galaxies at z = 0.55 (Hyde & Bernardi 2009, did not provide an estimate of the intrinsic scatter around the mass-size relation).

The last factor in Equation (25) describes the distribution in the power-law parameters, which is the main quantity of interest for our measurement. We assumed a Gaussian functional form in both m₅ and γ:

$$\begin{array}{c}\hfill \mathrm{P}(m_5,\gamma )=\mathcal{N}({\mu }_5(m_{\ast },r_{\mathrm{e}}),{\sigma }_5^2)\mathcal{N}({\mu }_{\gamma }(m_{\ast },r_{\mathrm{e}}),{\sigma }_{\gamma }^2).\end{array}$$(30)

We let the average of m₅ and γ to scale with stellar mass and size as

$$\begin{array}{cc}\hfill {\mu }_5(m_{\ast },r_{\mathrm{e}})& ={\mu }_{5,0}+{\beta }_5(m_{\ast }-11.3)+{\xi }_5(r_{\mathrm{e}}-{\mu }_R(m_{\ast }))\hfill \end{array}$$(31)

$$\begin{array}{cc}\hfill {\mu }_{\gamma }(m_{\ast },r_{\mathrm{e}})& ={\mu }_{\gamma ,0}+{\beta }_{\gamma }(m_{\ast }-11.3)+{\xi }_{\gamma }(r_{\mathrm{e}}-{\mu }_R(m_{\ast })).\hfill \end{array}$$(32)

This model is similar to that adopted by Sonnenfeld et al. (2013), but differs from it in a few aspects. First, here we modelled explicitly the normalisation of the total density profile, described by means of m₅, while Sonnenfeld et al. (2013) focused only on γ. Second, the coefficients ξ₅ and ξ_γ describe a scaling with excess size, that is the ratio between the half-light radius of a galaxy and the average size of a galaxy of the same stellar mass (the definition of ξ_γ in Sonnenfeld et al. 2013, is different). Third, we did not allow for a redshift evolution, since our sample spans a much smaller redshift range than that of Sonnenfeld et al. (2013). As we show in Section 5, this model is sufficiently complex to provide a good fit to the available data.

4.3. Source galaxy distribution

The next ingredient needed for our analysis was the redshift probability distribution of the sources that can be strongly lensed, P_s^(eff)(z_s), introduced in Section 3.2. According to Equation (16), this quantity is the redshift distribution of the background source population, marginalised over the source light parameters, and weighted by the light-dependent component of the lens finding probability, l. Our approach was to directly model P_s^(eff)(z_s) empirically, so a full characterisation of l was not needed. Our model is a Gaussian distribution, with mean μ_zs and dispersion σ_zs:

$$\begin{array}{c}\hfill {\mathrm{P}}_{\mathrm{s}}^{(\mathrm{eff})}(z_{\mathrm{s}})=\frac{1}{\sqrt{2\pi }{\sigma }_{z_{\mathrm{s}}}}exp\{-\frac{{(z_{\mathrm{s}}-{\mu }_{z_{\mathrm{s}}})}^2}{2{\sigma }_{z_{\mathrm{s}}}^2}\}.\end{array}$$(33)

4.4. Lensing cross-section

A key ingredient for correcting our inference for selection effects is the strong lensing cross-section, σ_SL. In Section 3.2 we split σ_SL into a part that depends on the source brightness, and a geometric factor g. This factor is the strong lensing cross-section of a source with a reference light distribution. In order to compute g, then, we needed to define exactly what a detection of a strong lens consists of, and what the light distribution of this reference source is. The detection of SLACS lenses is both spectroscopic and photometric. We then defined the spectroscopic detection on the basis of the total emission line flux within a circular aperture of 1.5″ (the spectroscopic fibre), after convolution with the PSF of the SDSS, and the photometric detection based on whether there are multiple images in the HST data. This entailed specifying both the spectrum and surface brightness distribution of the reference source.

To simplify the problem, we approximated the mass of the lenses as circular and the sources as point-like. The first assumption is justified by the fact that the strong lensing cross-section is a weak function of ellipticity; the second is a good approximation of cases in which the source size is smaller than the Einstein radius of the lens (see Figs. 8 and 9 of Sonnenfeld et al. 2023). Under these assumptions, we defined the emission line flux of the reference source to be 1/3 of the total fibre flux needed for a spectroscopic detection, and the broad-band flux to be the same as that needed for a photometric detection of an image. In other words, our reference source would be classified as strongly lensed if (1) the total emission line flux within the fibre, from all images and smeared by the PSF, is magnified by at least a factor of three with respect to the intrinsic flux; and (2) at least two images with magnification larger than one are produced.

Figure 2 shows the strong lensing cross-section of a power-law lens with θ_E = 1.2″ with respect to the reference source (green curve), as a function of the density slope γ. We computed this quantity from Equation (13) by simulating lensed sources, solving the lens equation, predicting image fluxes, convolving the emission line flux by the SDSS PSF, which we took to have an FWHM of 1.5″, and averaging over the source position. Along with it, we also computed σ_SL for sources with both emission line and broadband flux rescaled up to a factor of two with respect to the reference. This allowed us to assess the accuracy of the separability assumption of Equation (14). For values of γ > 1.8, σ_SL increases with the brightness of the source in a way that is weakly dependent on the density slope, roughly confirming our ansatz. For smaller values, however, σ_SL becomes independent of the source brightness: this is because lenses with a shallow density profile produce sets of highly magnified images, which can be detected regardless of the intrinsic brightness of the source (down to some limiting flux). The approximation of Equation (14) therefore breaks down for very bright sources or for lenses with a very shallow density profile. More generally, the accuracy of our approximation is higher the more similar our reference source is to the typical strongly lensed source of the SLACS sample. Since the choice of the reference source is arbitrary, we also carried out the analysis with a different definition of it, by setting its emission line flux to half of that necessary for a spectroscopic detection (instead of a third). The results did not change appreciably.

Fig. 2. Strong lensing cross-section of a circular power-law lens and a point source. A lensing event is defined as the detection of an emission line within an aperture of 1.5″, and of at least two images in photometric data. Top: σSL of a lens with θE = 1.2″ as a function of the slope γ. The thick green line corresponds to a source with intrinsic emission line flux equal to one third of the flux required for a spectroscopic detection, and a broadband flux equal to that needed for a photometric detection. Different lines were obtained by rescaling both the emission line and broadband flux by the amount indicated in the legend. Bottom: geometric factor of the lensing cross-section, that is the value of σSL for the reference source (corresponding to the green line in the top panel), as a function of both θE and γ.

The bottom panel of Figure 2 shows the geometric factor of the lensing cross-section, as a function of θ_E and γ, that is σ_SL of our reference source. For γ > 2, the lensing cross-section peaks at values of the Einstein radius close to the fibre size, then drops to zero for larger values of θ_E. This is a distinct signature of the spectroscopic detection requirement: for lenses with a steep density profile and large Einstein radius, any images that form within the spectroscopic fibre are not sufficiently magnified to lead to a detection (see also Arneson et al. 2012).

4.5. Lens finding probability

Finally, we needed to prescribe a model for the probability of finding a detectable lens. As we explained in Section 2.1, the main quantity that determined the discovery of a SLACS lens was the estimated lensing cross-section, which is the square of Equation (1). Therefore, we asserted that

$$\begin{array}{c}\hfill {\mathrm{P}}_{\mathrm{find}}({\mathit{\psi }}_{\mathrm{g}},{\mathit{\psi }}_{\mathrm{s}})={\mathrm{P}}_{\mathrm{find}}\left({\theta }_{\mathrm{E}}^{(\mathrm{est})}\right),\end{array}$$(34)

which is a special case of Equation (11), since ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ depends only on lens properties and on the source redshift (see Equation (1)). Our choice for the functional form of P_find is the following sigmoid, described by two parameters:

$$\begin{array}{c}\hfill {\mathrm{P}}_{\mathrm{find}}({\theta }_{\mathrm{E}}^{(\mathrm{est})})=\frac{1}{1+exp[-a({\theta }_{\mathrm{E}}^{(\mathrm{est})}-{\theta }_0)]}.\end{array}$$(35)

According to our model, the lens finding probability is zero for small values of the estimated Einstein radius, and reaches one for values of ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ larger than θ₀.

Predicting ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ for the parent population of lenses requires accounting for the effects of observational noise on the velocity dispersion measurements. We adopted a 6.25% uncertainty on all σ_ap measurements, which is the median relative uncertainty on σ_ap of the SLACS lenses. By comparison, Bolton et al. (2008) quoted a value of 7% as representative of the uncertainty on σ_ap measurements of the SLACS lenses. Birrer et al. (2020) argued, on the basis of a mismatch between the inferred properties of the SLACS and a different set of lenses, that the velocity dispersion measurements of SLACS lenses have been underestimated. They allowed for a systematic uncertainty, which they constrained to be approximately 6%, to be added in quadrature to the SLACS values. We repeated the analysis with these added uncertainties, but found that the goodness-of-fit worsened. Therefore we stuck with the uncertainties provided by Auger et al. (2009) for our analysis.

4.6. Fundamental plane constraint

Our goal was to infer the properties of the parent population of galaxies, and we could use information on the non-lens galaxies for this purpose. The model that we introduced so far makes use of the stellar mass-redshift distribution of the parent population of the lenses, and of the mass-size relation of early-type galaxies. A natural extension would have been to also consider the velocity dispersion distribution. Our model makes predictions for the velocity dispersion of a galaxy, by means of Jeans modelling. Therefore, it is in principle possible to fit the model directly to the observed velocity dispersion measurements of the entire parent population. We chose not to, for the following reasons. First, with a sample of ∼4 × 10⁵ objects, the constraints from the dynamical model of the parent population would dominate over the lensing data. Second, the simplified assumptions on which our dynamical analysis is based would end up completely dominating the error budget: we would need to adopt a much more complex model of the stellar orbits and of the measurement systematics to properly fit such a large dataset. Third, the inference is sensitive to the input values of the stellar mass, size and velocity dispersion of the galaxies, and to their uncertainties. Since the stellar mass and half-light radius of the lenses and non-lenses have been measured with different procedures, separate dynamical analyses of the two samples can easily produce inconsistent posterior probabilities, which therefore cannot be combined.

In light of these difficulties, our compromise solution was to adopt a prior on the velocity dispersion distribution of the parent population that matches the observations, but allowing for some degree of flexibility. In particular, we set a prior on the coefficients of the stellar mass fundamental plane, the observed correlation between stellar mass, half-light radius and central velocity dispersion of early-type galaxies (see e.g. Bernardi et al. 2020). We defined the stellar mass fundamental plane as follows:

$$\begin{array}{c}\hfill log({\sigma }_{\mathrm{ap}})\sim \mathcal{N}({\mu }_{\mathrm{FP},0}+{\beta }_{\mathrm{FP}}(m_{\ast }-11.3)+{\xi }_{\mathrm{FP}}(r_{\mathrm{e}}-{\mu }_R(m_{\ast })),{\sigma }_{\mathrm{FP}}^2).\end{array}$$(36)

We could not fit directly Equation (36) to the parent sample galaxies, because we lacked measurements of their half-light radius (or rather, we lacked measurements that are consistent with those of the SLACS lenses). So, we fitted the stellar mass-velocity dispersion relation instead, to constrain the first two coefficients in Equation (36). We assumed a quadratic dependence of log σ_ap on m_* for the fit. We obtained μ_FP, 0 = 2.342 and β_FP = 0.258 for galaxies with stellar mass close to the pivot value of m_* = 11.3, and adopted a Gaussian prior centred on these values. The statistical errors on these coefficients are very small, but we allowed for a 0.03 dispersion to accommodate possible inconsistencies between the measurement on the SLACS sample and the parent population. Finally, we put a prior on the intrinsic scatter parameter σ_FP. This is an important quantity that defines the tightness of the fundamental plane relation. We defined our prior on σ_FP in a way to accommodate both measurements from the literature and from the SLACS sample. Sonnenfeld (2020) measured σ_FP = 0.047 ± 0.02 on a sample of early-type galaxies at z = 0.2, and, when fitting Equation (36) to the SLACS data used in this paper, we obtained σ_FP = 0.039 ± 0.008. On the basis of these results, we chose the following prior:

$$\begin{array}{cc}\hfill \mathrm{P}({\mu }_{\mathrm{FP},0})& =\mathcal{N}(2.342,0.{030}^2)\hfill \\ \hfill \mathrm{P}({\beta }_{\mathrm{FP}})& =\mathcal{N}(0.258,0.{030}^2)\hfill \end{array}$$(37)

$$\begin{array}{cc}\hfill \mathrm{P}({\sigma }_{\mathrm{FP}})& =\mathcal{N}(0.047,0.{008}^2).\hfill \end{array}$$(38)

In practice, for a given value of the model parameters η we fitted Equation (36) to a large mock population of galaxies, then weighted the corresponding posterior probability by Equation (37). Since we were mostly interested in obtaining an accurate model for the high-mass end of the galaxy distribution, where the bulk of the SLACS lenses lie, we restricted this fit to the stellar mass range m_* > 11.0.

We found that the prior of Equation (37) had a significant impact on our inference. Therefore, in Section 5, we show results obtained when both including or excluding it.

4.7. Inference procedure

Our fiducial model is described by the following set of free parameters:

$$\begin{array}{cc}\hfill \mathit{\eta }\equiv \{& \overline{m},\alpha ,m_t^{(0)},\dots ,m_t^{(5)},{\sigma }_t,{\mu }_{5,0},{\beta }_5.{\xi }_5,{\mu }_{\gamma ,0},{\beta }_{\gamma },{\xi }_{\gamma },{\sigma }_{\gamma },{\mu }_{z_{\mathrm{s}}},\hfill \\ \hfill & {\sigma }_{z_{\mathrm{s}}},{\theta }_0,a\},\hfill \end{array}$$(39)

where m_t⁽⁰⁾, …, m_t⁽⁵⁾ are the coefficients of the polynomial function describing the truncation mass m_t in Equation (27). We first constrained the stellar mass function parameters, by fitting the model of Equation (26) to the stellar mass measurements of the parent population. The inferred values are reported in Table 1. In virtue of the very large sample size, formal uncertainties are very small. Therefore we fixed these parameters to their best-fit value and proceeded to infer the remaining parameters by fitting the model to the SLACS lenses.

The data vector of each SLACS lens consists of the lens redshift, source redshift, lens stellar mass, half-light radius, stellar velocity dispersion, and Einstein radius:

$$\begin{array}{c}\hfill {\mathit{d}}_i\equiv \{z_{\mathrm{g}},m_{\ast }^{(\mathrm{obs})},r_{\mathrm{e}},{\sigma }_{\mathrm{ap}}^{(\mathrm{obs})},{\theta }_{\mathrm{E}}^{(\mathrm{obs})},z_{\mathrm{s}}\}.\end{array}$$(40)

Stellar mass and velocity dispersion have associated uncertainties, on the order of 25% and 6% respectively. Redshift and half-light radius are instead measured very precisely, and we neglected their uncertainties. Einstein radii are measured with a precision of a few percent (Etherington et al. 2022). To further simplify the inference procedure, we neglected their uncertainties as well. We justified this approximation after the inference, by noting that our measurement on the average mass within 5 kpc has a larger relative uncertainty, which is therefore dominated by other aspects of the data (most likely the velocity dispersion measurements).

We assumed a Gaussian form for the likelihood of the stellar mass and velocity dispersion measurements. We calculated the model stellar velocity dispersion, σ_ap, using the spherical Jeans equation under the assumption of isotropic orbits. We discuss the possible implications of these assumptions in Section 6.

From Equations (15) and (17) we could write the likelihood term for a SLACS lens as

$$\begin{array}{cc}\hfill \mathrm{P}({\mathit{d}}_i|\mathit{\eta })={\displaystyle \int }& \mathrm{d}{m}_{\ast }\mathrm{d}\gamma \mathrm{d}{m}_5\mathrm{P}(m_{\ast }^{(\mathrm{obs})}|m_{\ast })\mathrm{P}({\sigma }_{\mathrm{ap}}^{(\mathrm{obs})}|{\sigma }_{\mathrm{ap}}(m_5,\gamma ))\times \hfill \\ \hfill & \mathrm{P}({\theta }_{\mathrm{E}}^{(\mathrm{obs})}|{\theta }_{\mathrm{E}}(m_5,\gamma )){\mathrm{P}}_{\mathrm{SL}}^{(\mathrm{eff})}(m_{\ast },m_5,\gamma ,z_{\mathrm{s}}|\mathit{\eta }).\hfill \end{array}$$(41)

In the factor P(σ_ap^(obs)|σ_ap(m₅, γ)), which describes the likelihood of the observed velocity dispersion, the model-predicted velocity dispersion is obtained via the spherical Jeans equation with isotropic orbits. Under the assumption of negligible uncertainty on the Einstein radius, the likelihood term in ${\theta }_{\mathrm{E}}^{(\mathrm{obs})}$ became a Dirac delta function, which we could integrate over after a variable change:

$$\begin{array}{c}\hfill (\gamma ,m_5)\to (\gamma ,{\theta }_{\mathrm{E}}(\gamma ,m_5)).\end{array}$$(42)

The result is

$$\begin{array}{cc}\hfill \mathrm{P}({\mathit{d}}_i|\mathit{\eta })& ={\displaystyle \int \mathrm{d}}{m}_{\ast }\mathrm{d}\gamma \mathrm{P}(m_{\ast }^{(\mathrm{obs})}|m_{\ast })\mathrm{P}({\sigma }_{\mathrm{ap}}^{(\mathrm{obs})}|{\sigma }_{\mathrm{ap}}(m_5^{(\mathrm{obs})}(\gamma ),\gamma ))\times \hfill \\ \hfill & \times {\mathrm{P}}_{\mathrm{SL}}^{(\mathrm{eff})}(m_{\ast },m_5^{(\mathrm{obs})}(\gamma ),\gamma ,z_{\mathrm{s}}|\mathit{\eta }){\left|\frac{\mathrm{d}{m}_5}{\mathrm{d}{\theta }_{\mathrm{E}}}\right|}_{{\theta }_{\mathrm{E}}={\theta }_{\mathrm{E}}^{(\mathrm{obs})}},\hfill \end{array}$$(43)

where $m_5^{(\mathrm{obs})}(\gamma )$ is the value of m₅ that reproduces the observed Einstein radius, for a given γ.

We assumed a uniform prior in all model hyper-parameters, multiplied by the fundamental plane prior of Equation (37). We sampled the posterior probability distribution with a Markov chain Monte Carlo (MCMC). We carried out the integrals of Equation (40), as well as the integral of Equation (18) to normalise P_SL^(eff), via Monte Carlo integration.

4.8. Lens-only inference

One of the goals of this work was to quantify the impact of selection effects on a strong lensing measurement of the mass-structure of galaxies. To do this, we also carried out the inference while ignoring all lensing selection terms. We did this by taking the prior factor P_SL^(eff) in Equation (42) to simply be the product between the distribution in (m₅, γ) of Equation (30) and a Gaussian distribution in m_*, meant to describe solely the SLACS lenses.

5. Results

5.1. Individual lens constraints

As a first thing we examined the constraining power of each individual lens. A power-law model has two free parameters: m₅ and γ. Each lens provided two constraints: the Einstein radius, which we assumed to be measured exactly, and the stellar velocity dispersion within the spectroscopic aperture, σ_ap. Therefore, the range of parameter values allowed by the data is a line in m₅ − γ space. These constraints are shown in Figure 3, along with the values of the velocity dispersion corresponding to the model. Each line spans the 68% credible region of the posterior probability, assuming a flat prior on both m₅ and γ.

Fig. 3. Constraints on m5 and γ, and predicted values of the velocity dispersion, for individual lenses. Each line covers the 68% credible region of the posterior probability. Lines are colour-coded by the physical size of the Einstein radius. Dots mark the maximum posterior probability values.

The mass enclosed within 5 kpc is generally very well constrained: the median uncertainty on m₅ is as small as 0.013. The uncertainty is largest for the lenses with Einstein radius that differs the most from 5 kpc: for these lenses, obtaining m₅ requires an extrapolation of the power-law profile, which introduces a correlation between γ and m₅. In all cases, γ is strongly correlated with the velocity dispersion: given the Einstein radius of a lens, the steeper the density profile the larger the value of σ_ap. Together with the velocity dispersion selection criterion that we described in Section 2.1, this has important implications for the selection function of the SLACS sample.

5.2. Parameter inference

In Table 2 we report the median and 68% credible region of the marginal posterior probability of each model parameter. We list the values obtained with and without the fundamental plane prior of Equation (37), as well as the values obtained when fitting the SLACS lenses while not accounting for selection effects. In Figures 4, 5 and 6 we show projections of the posterior probability on selected pairs of model parameters. Given the large number of free parameters, we avoid showing all possible pairwise correlations. The MCMC samples that we used to make these figures is available online for further analysis.

Fig. 4. Posterior probability distribution of the model. Parameters shown: average log M5 and γ at fixed stellar mass and half-light radius, and intrinsic scatter around the average. Filled contours: model describing SLACS lenses, without correction for selection effects (see Section 4.8). Thick solid curves: fiducial model. Thin solid curves: model with no fundamental plane prior of Equation (37). Contour levels correspond to 68% and 95% enclosed probability.

Fig. 5. Same as Figure 4, but showing parameters describing the scaling of m5 and γ with stellar mass and size. Dashed lines indicate values of the parameters corresponding to the case in which galaxies are homologous systems with an average density slope of γ = 2.00.

Fig. 6. Same as Figure 4, but showing parameters describing the lens finding probability, defined in Equation (35), and the average density slope parameter.

Figure 4 focuses on the mean m₅ and γ at fixed stellar mass and size, and scatter around the average. According to the lens-only inference (cyan filled contours), SLACS lenses of m_* = 11.3 and average size for their mass have an average density slope of μ_γ, 0 = 2.091 ± 0.024. This value is in good agreement with the measurement of Auger et al. (2010) on the same dataset. Once corrected for selection effects, the same quantity for the parent sample galaxies (magenta solid curves) is μ_γ, 0 = 1.99 ± 0.03. The average projected mass within 5 kpc, at the pivot stellar mass and average size, is inferred to be μ_5, 0 = 11.358 ± 0.010 in the lens-only analysis and μ_5, 0 = 11.332 ± 0.013 after the selection effects correction. Based on these results, we can say that SLACS lenses are on average slightly more massive than parent sample galaxies (they have a larger m₅ at fixed stellar mass), and have an overall steeper density profile. The larger mass can be easily explained by the increase in lensing cross-section with m₅. The bias on γ is less intuitive, because, as the bottom panel of Figure 2 shows, the lensing cross-section is generally larger for values of γ < 2: this alone would lead to the lenses having a shallower profile than the non-lenses. Instead, as we show in Section 5.4, the main reason for the bias in γ lies in the lens finding probability term, P_find.

When comparing the inference based on the fiducial model (magenta curves) with that without the fundamental plane prior (black curves), we can see that the strong lensing bias (the difference with respect to the SLACS-only posterior probability) is larger in the latter case. To understand this fact, it is useful to consider the intrinsic scatter parameters. We can see that σ_γ and σ₅ are anti-correlated with the respective average parameters, μ_γ, 0 and μ_5, 0. This means that the strong lensing data can be reproduced both with an underlying population that has a relatively large bias and intrinsic scatter, or with a smaller bias and scatter. But a larger scatter in m₅ and γ translates into a larger scatter around the fundamental plane. Therefore the solution with the larger bias is disfavoured by the prior of Equation (37).

Figure 5 shows the posterior probability distribution in the parameters describing how m₅ and γ scale with stellar mass and size. In order to gain intuition on how to interpret these measurements, it is useful to first predict the values of these parameters in an ideal case. For instance, we can make the assumption that early-type galaxies are homologous systems, which is to say that galaxies with different values of the stellar mass and half-light radius are scaled-up versions of one another. Although strict homology is already ruled out by dynamical observations, such as the tilt of the fundamental plane (see e.g. Bertin et al. 2002), we can still use it as a reference limiting case. Under the homology assumption, the average density slope γ is a constant, the total projected mass within the half-light radius is also a constant, and, at fixed half-light radius, the total mass within any aperture scales linearly with stellar mass. The values of the parameters of our model under the homology assumption are β₅ = 0.35, ξ₅ = −1.00, β_γ = 0 and ξ_γ = 0.

Our measured values are inconsistent with the homology assumption. In particular, β₅ = 0.59 ± 0.05, which means that the mass enclosed within 5 kpc grows more rapidly than it would in the case of perfect homology. It could be explained with the dark matter mass within the half-light radius growing more than linearly with stellar mass, or also with the ratio between the true stellar mass and our estimate m_* being an increasing function of m_*. This is not a new result, but was previously noted by Auger et al. (2010) in their analysis of SLACS lenses. Our work demonstrates that it holds also when accounting for selection effects.

The dependence of m₅ on excess size is also very different from that of homologous systems: ξ₅ = −0.11 ± 0.13, instead of −1. This means that m₅ responds relatively mildly to a change in half-light radius. One possible explanation is that the dark matter distribution is roughly independent of galaxy size. Weak lensing observations have shown that halo mass is a weak function of size, at fixed stellar mass (Sonnenfeld et al. 2018, 2019, 2022), so they support this hypothesis (though see Charlton et al. 2017, for a counter-example). Finally, we put an upper limit on the correlation between the density slope γ and stellar mass, β_γ = 0.03 ± 0.08, and found a negative correlation with galaxy size at fixed stellar mass: ξ_γ = −0.67 ± 0.31. This is consistent with the measurements of Sonnenfeld et al. (2013).

Among the model parameters are those describing the effective source redshift distribution, P_s^(eff), introduced in Section 4.3. The mean and dispersion of the Gaussian of Equation (16) are constrained to be μ_zs = 0.48 ± 0.04 and σ_zs = 0.22 ± 0.02, and the two parameters are anti-correlated. By comparison, the mean and standard deviation of the redshift of the SLACS sources is 0.61 and 0.18, respectively. The fact that we inferred a broader distribution, centred at a smaller redshift, was expected: as we explained in section 3.2, P_s^(eff) describes the redshift distribution of unlensed sources that can be detected. Because the lensing cross-section increases for higher-redshift sources, the distribution of SLACS sources is shifted to higher redshifts. We found no signs of correlation between μ_zs or σ_zs and any of the other parameters of the model. From this we concluded that the exact choice of parameterisation of P_s^(eff) is not a critical element in our model.

Finally, we constrained the lens finding probability parameters, defined in Equation (35). Figure 6 shows the posterior probability in the cutoff radius θ₀ and log a, defined in Equation (35), as well as in the average density slope μ_γ, 0. The cutoff on ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ is θ₀ = 0.93 ± 0.08, meaning that lenses with an estimated Einstein radius smaller than this value had a small probability of being selected for photometric follow-up in the SLACS survey. The parameters θ₀ and a are anti-correlated. This indicates that the data cannot constrain the detailed shape of the lens finding probability function, but can allow either a P_find with a sharp cutoff at a smaller ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$, or one with a more progressive increase centred at a larger value of ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$. Neither θ₀ nor log a are strongly correlated with the average density slope parameter μ_γ, 0, as can be seen in Figure 6, or with any other parameter of the model.

5.3. Goodness of fit

We assessed the goodness of fit of the model with posterior predictive tests, as follows. We sampled sets of values of the model parameters η from the posterior probability distribution. For each set of values, we randomly drew a set of 59 lenses (the same as the SLACS sample) from the strong lens distribution P_SL^(eff). For each lens, we predicted the set of observables that we used in the analysis, adding observational noise to the velocity dispersion. We adopted a Gaussian noise with a 6.25% relative uncertainty for this purpose. We then defined an ensemble of test quantities {T_i} based on observations, and quantified the probability of drawing more extreme values of each test quantity from a posterior predicted sample. Very small (P < 5%) or very large (P > 95%) probabilities would mean that it is unlikely that the observed sample has been drawn from the model. The choice of test quantities is not unique, but must be made on the basis of what aspects of the data we want the model to successfully reproduce. We chose them to be the mean, standard deviation, 10%- and 90%-ile of the distribution in observed Einstein radius, velocity dispersion and density slope γ.

Figure 7 shows the posterior predicted distributions of these twelve test quantities. The probability of producing samples of lenses with more extreme values of each test quantity is shown in the corresponding panel of Figure 7. None of the observed test quantities correspond to very small or very large percentiles of the posterior predicted distribution. Therefore our model is able to reproduce these data.

Fig. 7. Posterior predictive tests of goodness of fit. Each panel shows the posterior predicted (histogram) and observed value (dashed line) of the test quantity indicated at the bottom. Units of Einstein radius are in arcsec and of the velocity dispersion are in km s−1.

We also checked whether the model can recover the velocity dispersion distribution of the parent sample. The bottom panel of Figure 8 shows the posterior predicted average σ_ap (magenta band) as a function of stellar mass, and our quadratic fit to the m_* − σ_ap distribution of the parent sample (green curve). The model matches the data over a broad range in stellar mass, but departs from it at the high-mass end. This is not a goodness-of-fit issue: as we explained in Section 4.6, we only used observations around the pivot stellar mass m_* = 11.3 to constrain the model. Nevertheless, it could mean that the model does not extrapolate well at large stellar masses, where it tends to overpredict the stellar velocity dispersion. However, the source of the discrepancy could more simply be an inconsistency between the way stellar mass is defined for the SLACS lenses and the parent sample galaxies. Although we rescaled the m_* measurements of the parent sample galaxies so that the stellar masses of the SLACS lenses would be consistent on average, we cannot exclude a mass-dependent offset, for instance as a result of different choices in the photometric measurement or stellar population synthesis phase. It could also be that a quadratic model is not sufficiently flexible to describe the high-mass end of the m_* − σ_ap relation. Our main goal is to understand the properties of galaxies at m_* ≈ 11.3, where the bulk of the SLACS lenses lie. Therefore we did not investigate this issue further.

Fig. 8. Posterior predicted trends with stellar mass, for galaxies of average size. Top: average projected mass enclosed within 5 kpc. Middle: average density slope. Bottom: average stellar velocity dispersion. Each band marks the 68% credible region. The magenta band shows the prediction for the parent population. The cyan band shows the prediction from the lens-only analysis, without accounting for selection effects. Dashed blue lines show the prediction for the population of SLACS lenses. Solid black lines show the prediction for the population of detectable lenses (obtained by removing the lens finding probability Pfind from the posterior prediction). Error bars are the values measured for the SLACS lenses, assuming a flat prior on both m5 and γ. Uncertainties on the stellar mass are uncorrelated with those on m5, γ and σap. In the bottom panel, the green line shows our quadratic fit to the m* − σap relation of the parent sample.

5.4. Insight from the model

Figure 8 shows the inferred trend of m₅, γ and velocity dispersion with stellar mass, as predicted from the posterior probability of the model parameters, along with the data. In the middle panel we can see a significant difference between the average γ of the parent population (magenta band), and that inferred from the lens-only analysis (cyan band). According to these results, the average density slope of the SLACS lenses appears to be ∼0.1 larger than the corresponding value of their parent population, at the pivot stellar mass m_* = 11.3. This raises the question of whether this bias is intrinsic to strong lenses in general, or is a peculiarity of SLACS.

One aspect that makes the SLACS sample unique is the pre-selection of targets for photometric follow-up based on the estimated Einstein radius, which we described empirically with the lens finding probability term P_find. We can assess the impact of this choice with posterior prediction. For a given galaxy population drawn from the posterior probability, we considered the population of detectable lenses: these are strong lens systems that are detected spectroscopically and that would be detected photometrically if followed-up with HST. We obtained such a sample by setting P_sel = P_det, dropping P_find in Equation (9). We measured the average m₅ and γ of detectable lenses in bins of stellar mass. The resulting 68% credible regions are shown in Figure 8, enclosed within the black lines. At stellar masses around m_* ∼ 11.3, where the bulk of the SLACS data lies, the posterior predicted average γ of detectable lenses is indistinguishable from the value inferred for the lens population (magenta band). This means that, according to our model, the strong lensing bias on the density slope is entirely attributable to the ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$-based target prioritisation criterion. At fixed lens and source redshift, selecting lenses with larger ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ means selecting preferentially lens galaxies with larger velocity dispersion. Since velocity dispersion correlates positively with γ, this results in lenses with steeper density profiles on average.

The above argument explains part of the difference in the average density slope that we inferred with or without taking selection effects into account, but not its entirety. From our model, we predicted the average γ of SLACS lenses, this time using the full expression for P_sel. The corresponding 68% confidence region is show by the dashed blue lines in Figure 8. This is indeed larger than the corresponding value for the lens population, but is also smaller than the value that we inferred from the lens-only analysis (cyan band in Figure 8). The difference in γ is about 0.04 at the pivot stellar mass m_* = 11.3. This result might seem counterintuitive, especially given that the posterior predicted tests showed no sign of inconsistency with the observed values of γ (bottom panels of Figure 7).

The root of this residual discrepancy is observational errors. Not only the selection in ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ tends to favour galaxies with intrinsically larger velocity dispersion, but also those for which σ_ap has been systematically overestimated. The probability of a detected lens being included in the SLACS sample increases with ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$. This means that, at fixed lens and source redshift and velocity dispersion, lenses with a larger σ_ap^(obs) are more likely to be selected. Due to the correlation between σ_ap and γ, such lenses have also a larger observed value of the density slope. Figure 9 illustrates this effect. The top panel shows the values of ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ and γ of the SLACS lenses that are allowed by the velocity dispersion measurements. The bottom panel shows the inferred P_find distribution, which describes the probability of a detected lens to be included in the SLACS sample, and which depends on ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$. The uncertainties on the velocity dispersion translate into relatively large errors on ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$, compared to the scale where P_find varies. Therefore, observational errors can play an important role in shaping the lens sample, particularly at low values of ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$.

Fig. 9. Estimated Einstein radius, density slope and lens finding probability. Top: 68% probable values of ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ and density slope γ of the SLACS lenses, where ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ is defined in Equation (2). Bottom: 68% credible region of the lens finding probability, defined in Equation (35). The correlation between ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ and γ makes so that lenses with overestimated γ are more likely to be included in the sample.

To verify this conjecture, we carried out the following test. Given a posterior predicted sample of SLACS-like lenses, we considered the estimated Einstein radius under the assumption of a singular isothermal profile and perfect knowledge of the stellar velocity dispersion:

$$\begin{array}{c}\hfill {\theta }_{\mathrm{E}}^{(\mathrm{SIS})}({\sigma }_{\mathrm{ap}})=4\pi {\left(\frac{{\sigma }_{\mathrm{ap}}}{c}\right)}^2\frac{D_{\mathrm{ds}}}{D_{\mathrm{s}}}\end{array}$$(44)

(the difference with Equation (2) is in the value of the velocity dispersion that is used). Then, we looked at the relative error on σ_ap and the error on γ as a function of ${\theta }_{\mathrm{E}}^{(\mathrm{SIS})}$. These are shown in Figure 10, marginalised over the values of the parameters η. As expected, both quantities are positively biased at small values of ${\theta }_{\mathrm{E}}^{(\mathrm{SIS})}$. The average posterior predicted observed γ is, on average (i.e. when marginalising over the model parameters), 2.10, which is larger than the true value by 0.04. This bias disappears for the population of detectable lenses. In conclusion, SLACS lenses have on average a steeper density slope than their parent sample, and the value of γ of some of their lenses has likely been overestimated.

Fig. 10. Bias on velocity dispersion and density slope of SLACS-like lenses. Top: error on γ as a function of ${\theta }_{\mathrm{E}}^{(\mathrm{SIS})}$, defined in Equation (43), of posterior predicted samples of lenses, marginalised over the values of the model parameters η. Bottom: relative error on σap. The solid line indicates the average in bins of ${\theta }_{\mathrm{E}}^{(\mathrm{SIS})}$.

Observational errors aside, SLACS lenses have on average a larger velocity dispersion than parent sample galaxies of the same stellar mass, as can be seen in the bottom panel of Figure 8. The difference is most obvious at stellar masses m_* < 11.5. This is again due to the ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$-based selection, because the posterior predicted distribution in σ_ap of all detectable lenses (black lines in Figure 8) is instead very close to that of the parent population. Since the velocity dispersion is related to both stellar mass and size through the fundamental plane relation, it is interesting to determine whether the SLACS lenses lie on the fundamental plane or not. If they do, then their larger σ_ap must be due to them being more compact than the average at fixed stellar mass. If not, then their values of σ_ap are genuinely larger for their stellar mass and size. We found the latter to be the main reason.

To illustrate this point, we show in Figure 11 the posterior predicted distribution of the stellar mass-size relation (top), and the fundamental plane relation (bottom). SLACS lenses have indeed smaller sizes than the parent population, on average, but only at the low-mass end, and only by 0.04 dex at m_* = 11. When fitting the fundamental plane relation of Equation (36) to the model for the parent population, we inferred an average value of the coefficient describing the scaling of σ_ap with excess size of ξ_FP ≈ −0.20. This means that a 0.04 dex offset in size translates into an average difference in velocity dispersion of 0.008 dex (2%), corresponding to 3.5 km s⁻¹ at σ_ap = 180 km s⁻¹. This is much smaller than the difference in σ_ap between SLACS lenses and the parent population seen in the bottom panel of Figure 8. At the bottom of Figure 11 we show the posterior predicted distribution of the stellar velocity dispersion, as a function of the position on the fundamental plane (the mean of the distribution of Equation (36)), marginalised over the model parameters η. SLACS lenses lie on average above the fundamental plane: their velocity dispersion is 0.035 dex (8%) higher at values of stellar mass and size corresponding to σ_ap = 200 km s⁻¹, and the bias increases towards lower velocity dispersions. Thus, this is the main reason for SLACS lenses having higher values of σ_ap at fixed stellar mass. This conclusion depends to some extent on our choice of adopting the mass-size relation of Hyde & Bernardi (2009), which lies very close to that of the SLACS lenses (see also Auger et al. 2010). For a more robust assessment of the relative position of SLACS lenses and their parent population in (M_*, R_e, σ_ap) space we need consistent measurements of the three quantities, which we do not have. Nevertheless, the main point of this investigation is that, although the fundamental plane is a relatively tight relation, it is sufficiently thick to allow for a significant bias in velocity dispersion between the SLACS lenses and the parent population.

Fig. 11. Posterior predicted mass-size and fundamental plane relation. Top: distribution in stellar mass and half-light radius of the parent population (magenta regions) and the SLACS lenses (dashed contours), marginalised over the model parameters η. The dotted black line is the mass-size relation from Hyde & Bernardi (2009). The solid blue line indicates the average size of SLACS lenses as a function of stellar mass. Bottom: distribution in velocity dispersion, as a function of the predicted average given the stellar mass and size (the mean of the distribution in Equation (36)). The solid blue line indicates the average velocity dispersion of the SLACS lenses.

Another interesting question, related to the point above, is to what extent SLACS lenses are similar to galaxies from their parent population, once controlling for velocity dispersion. In Figure 12 we show the posterior predicted average true values of m₅ and γ of parent sample galaxies, SLACS lenses and detectable lenses, measured in bins of velocity dispersion. Strong lenses are remarkably similar to the parent population from this point of view: for instance, the average m₅ of SLACS lenses at σ_ap = 250 km s⁻¹ is smaller than that of the parent sample by only 0.02. Differences in γ are for the most part smaller than 0.01. These results corroborate previous claims that SLACS lenses are representative of galaxies of the same stellar velocity dispersion Treu et al. (2006, 2009), despite occupying a different region in (M_*, R_e, σ_ap) space. However, as shown previously in Figure 10, SLACS lenses have slightly biased velocity dispersion measurements. This bias needs to be taken into account when comparing them to other samples of galaxies.

Fig. 12. Posterior predicted trends with stellar velocity dispersion. Top: average projected dark matter mass within 5 kpc. Bottom: average density slope. The magenta band marks the 68% credible region for the parent population of galaxies. The dashed blue lines enclose the 68% credible region for the SLACS lenses. Grey lines in the top two panels span the values of the SLACS lenses that are allowed with 68% probability.

As a last thing, we computed the strong lensing bias b_SL of the SLACS survey. We introduced b_SL in Section 3.1: it is the function by which the distribution of parent sample galaxies is multiplied to produce the sample of strong lens galaxies observed in a survey. We computed b_SL from Equation (19) for an example galaxy at redshift z = 0.2 and with half-light radius R_e = 7 kpc, as a function of m₅ and γ. Figure 13 shows the result, marginalised over the model parameters η. The strong lensing bias peaks on a relatively narrow band in m₅ − γ space, which can be seen as a window function of the SLACS survey: galaxies that lie on that band are more likely to be selected as SLACS lenses, while those outside of it are missed. In Figure 13 we also overplotted lines of constant velocity dispersion, which are nearly parallel to the lines of constant b_SL. This matches our view that the selection function of SLACS is mostly one in velocity dispersion: at fixed σ_ap, the strong lensing bias does not vary much with either m₅ or γ, hence SLACS lenses have a similar density profile to regular early-type galaxies with the same velocity dispersion.

Fig. 13. SLACS Strong lensing bias bSL for a galaxy at z = 0.2 with Re = 7 kpc, as a function of m5 and γ. White lines connect points of constant velocity dispersion, with σap equal to the value indicated. We computed bSL from Equation (19), while marginalising over the model parameters η.

The expression for the strong lensing bias is independent of the distribution of the parent galaxy population P_g. One could then consider a different parent sample from the one used by SLACS, and use b_SL to predict what kind of lenses would be discovered in a SLACS-like campaign with this different sample. In general this could be a valid use of b_SL, but not in our case. This is because our model for the lens finding probability, which is meant to describe the prioritisation of targets for follow-up observations, depends explicitly on the observed velocity dispersion. But in a hypothetical sample with a different distribution in σ_ap, a SLACS-like prioritisation that selects a fixed number of targets at the top of the ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ distribution would correspond to a different selection in σ_ap, and consequently different P_find. The question of what a SLACS-like survey would look like with a different parent population is then an ill-posed one: more details on the follow-up strategy need to be provided.

6. Discussion

Joint lensing and dynamics studies of the SLACS lenses carried out so far have found that the average density slope of the 59 lenses in the sample is ⟨γ⟩≈2.1. According to our analysis, this value of the average γ is larger than that of regular early-type galaxies of the same stellar mass, by about 0.1. This result can have important implications for studies of the evolution in the density profile of massive early-type galaxies.

Ruff et al. (2011), Bolton et al. (2012), Sonnenfeld et al. (2013) showed how, when combining SLACS with other samples of lenses, the average γ tends to be larger at lower redshift, pointing towards a steepening of the total density profile with time. The common feature among these studies is that they all rely on the SLACS sample as a low-redshift pivot. Thus, if ⟨γ⟩ of SLACS were to be lowered, the conclusions would change dramatically. For instance, Sonnenfeld et al. (2013) measured d⟨γ⟩/dz = −0.31 ± 0.10 when combining SLACS with SL2S lenses. The average redshift of the SL2S sample is ⟨z⟩≈0.5, larger than that of SLACS by ≈0.3. This means that a systematic offset in ⟨γ⟩ of the SLACS lenses by 0.1 could account for the detected trend entirely. A similar argument can be made for the study of Bolton et al. (2012), which is based on the BELLS sample. BELLS lenses were discovered spectroscopically with a very similar method employed for SLACS, including a criterion for the prioritisation of photometric follow-up observations based on the estimated Einstein radius ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$. Crucially, however, BELLS did not rely on the observed velocity dispersion in the estimate of ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$, due to the low signal-to-noise ratio of the available spectroscopic measurements: instead, they assigned a fixed value of σ_ap = 250 km s⁻¹ to all lens candidates (Brownstein et al. 2012). As a result, the BELLS sample is immune from the selection in velocity dispersion that causes the bias in γ of the SLACS lenses.

In light of these findings, the evidence for a steepening with time of the total density profile of massive early-type galaxies needs to be revisited. SLACS lenses cannot be simply compared to lenses at higher redshift from different samples, but their velocity dispersion selection needs to be taken into account. It is possible that the value of d⟨γ⟩/dz will be revised upwards after correcting for selection effects. That would help reconcile lensing and dynamics measurements with some lensing-only studies (Tan et al. 2024; Sahu et al. 2024) and with theoretical predictions (Remus et al. 2013; Sonnenfeld et al. 2014; Xu et al. 2017; Shankar et al. 2018). Nevertheless, in order to know with certainty whether the steepening in γ is real or not, it is necessary to correct for selection effects also the other samples of lenses mentioned above. That, however, is beyond the scope of this work.

Our findings have implications for time-delay cosmography studies as well. Birrer et al. (2020) used the SLACS sample to put a prior on the lens mass model of time-delay lenses from the TDCOSMO sample (Millon et al. 2020), under the assumption that any differences between the SLACS and the TDCOSMO sample can be accounted for by scalings with the ratio between Einstein radius and half-light radius. Given the peculiar selection function of SLACS lenses, which tends to pick objects that are offset from the fundamental plane of early-type galaxies, it seems unlikely that the two samples can be combined in such a way without introducing biases. However, we are unable to quantify such a bias, which might well be within the statistical uncertainties, because the mass model used by Birrer et al. (2020) is more complex than the simple power-law on which our study is based.

All of our analysis was carried out under the assumption of a spherical mass distribution and isotropic orbits for the stars. This assumption allowed us to convert directly the Einstein radius and stellar velocity dispersion of each lens into an inference on the two free parameters of the power-law mass profile, (m₅, γ). Although the amount of flattening and anisotropy of the SLACS lenses has been found to be modest (Barnabè et al. 2011), high precision measurements would require us to allow for a departure from spherical symmetry and from the isotropy assumption. At fixed mass and tracer density profile, radial orbits tend to increase the stellar velocity dispersion. Fitting an isotropic power-law model to a radially anisotropic lens galaxy would then lead to an overestimate of the density slope γ (see Koopmans et al. 2006, for a quantification of this effect). We expect the orbital anisotropy to vary from galaxy to galaxy: this is a source of scatter in the velocity dispersion at fixed mass profile that we are absorbing into the inferred distribution in γ. As a result, we likely overestimated the strength of the correlation between γ and the velocity dispersion. Because, as we established, the nature of the bias in γ between the SLACS lenses and the parent population is entirely due to the selection in velocity dispersion, an intrinsic scatter in orbital anisotropy would mean that what we interpreted exclusively as a bias in γ is partly a bias in anisotropy: SLACS lenses are not as biased in γ with respect to the parent population, but have on average more radially anisotropic orbits. Appendix A explores the degeneracy between density slope and orbital anisotropy in more detail. A similar argument can be made regarding the three-dimensional structure. This does not change the main point of our analysis, which is that the way that SLACS galaxies were selected makes them different from galaxies with similar stellar mass and size, in terms of their internal structure.

One of the limiting factors in our analysis is the lack of consistent measurements of the stellar mass profile of galaxies of the SLACS lenses and the parent population galaxies. This prevented us from taking full advantage of the information available on the parent sample, and forced us to adopt the more relaxed prior on the fundamental plane described in Section 4.6. For future studies, it would be beneficial to have access to the same photometric data for both the lenses and the non-lenses. We expect that strong lens studies based on data from Euclid (Euclid Collaboration 2024) will satisfy this requirement.

7. Conclusions

We carried out a joint lensing and dynamical analysis of the SLACS sample of lenses. Our study differed from previous analyses of the same sample (Auger et al. 2010; Sonnenfeld et al. 2013; Posacki et al. 2015; Shajib et al. 2021) under one crucial aspect: we corrected the inference for selection effects. Our model for the selection function describes both the probability of detecting a lensed source in the spectroscopic and photometric data on which the survey was based, and the probability of a lens candidate being included in the final SLACS sample. This constitutes a step forward in complexity and fidelity with respect to earlier attempts (Sonnenfeld et al. 2015, 2018).

We constrained the distribution in the mass structure of massive early-type galaxies, under the assumption of a spherically symmetric power-law mass density profile and isotropic orbits. We found that galaxies with a stellar mass of log(M_*/M_⊙) = 11.3 and average size for their mass have on average a total projected mass enclosed within 5 kpc of log M₅ = 11.332 ± 0.013, and an average density slope of γ = 1.99 ± 0.03. By contrast, a fit of the same mass model to the SLACS lenses gives a larger average density slope by approximately 0.1. We identified the source of this bias to lie entirely in the way that the SLACS lens candidates were prioritised for photometric follow-up observations, which relied on the observed velocity dispersion. As a result, SLACS lenses have larger velocity dispersion than regular early-type galaxies at fixed stellar mass and size, and their velocity dispersion measurements are likely biased towards larger values. Nevertheless, they are fairly representative of galaxies of similar velocity dispersion (though they have slightly lower masses). We also corroborated evidence towards early-type galaxies not being homologous systems: their total mass grows more strongly with stellar mass and decreases less strongly with half-light radius than a purely homologous structure would suggest.

Because the source of the bias of the SLACS lenses is in the target prioritisation, which is a peculiarity of SLACS, we do not expect other lens samples to be similarly biased. This could have implications for previous claims of an evolution in the total density profile of massive early-type galaxies that relied on the combination of SLACS data with other lens samples.

Our study provides a reference on which to base future investigations that require knowledge of the selection function of a strong lens sample. One immediate use case is testing the predictions of theoretical models on the mass structure of early-type galaxies: this can be done directly by taking the model distribution in (m₅, γ) that we inferred from the data and comparing it with theory. If the goal is building a large sample of galaxies with similar internal structure as the SLACS lenses, for instance to carry out weak lensing measurements (such as in the study of Sonnenfeld et al. 2018), this can be done by matching the two samples in velocity dispersion, provided that differences in mass of 0.02 dex are deemed not important. If one wishes to use information from the SLACS sample as a prior for different samples of lenses, as done in the time-delay lensing analysis of Birrer et al. (2020), additional work is needed: the strong lensing bias of each strong lens sample needs to be determined separately. Finally, by applying the method developed in this paper to more complex models and different samples of lenses, we can remove the strong lensing bias on other properties of the mass structure of galaxies, such as the stellar mass-to-light ratio, the dark matter distribution, and their time evolution.

Data availability

Additional data for this article are available in a public data repository at: https://zenodo.org/records/13366012. The code and accompanying documentation are publicly available on GitHub at: https://github.com/astrosonnen/strong_lensing_tools/tree/main/papers/slacs_selection.

Acknowledgments

This work was initiated at the 2023 KICP workshop “Lensing at different scales: strong, weak, and synergies between the two”. I thank Tommaso Treu, Tom Collett, Tian Li, Anowar Shajib, Natalie Hogg, Birendra Dhanasingham and Raphael Gavazzi for their useful comments and suggestions.

References

Appendix A: Orbital anisotropy

Our measurements of the density profile of early-type galaxies rely mainly on strong lensing and stellar kinematics information. In particular, we used the spherical Jeans equation to fit dynamical models to the observed velocity dispersion data. The entire analysis was carried out under the assumption of isotropic orbits. Here we explore briefly how some of the results would change if we were to relax that assumption.

Orbital anisotropy in spherical systems is defined in terms of the tangential and radial components of the velocity dispersion, as

$$\begin{array}{c}\hfill {\beta }_{\mathrm{ani}}=1-\frac{{\sigma }_t}{{\sigma }_r}.\end{array}$$(A.1)

For isotropic orbits, σ_t = σ_r and β_ani = 0. For radially anisotropic orbits, β_ani > 0, while β_ani < 0 indicates tangential anisotropy.

One effect of fitting an isotropic power-law model to the velocity dispersion of a galaxy with non-zero anisotropy is to bias the inferred value of the density slope γ: at fixed mass profile, the velocity dispersion increases with β_ani, and so does the inferred value of γ. Therefore, to first approximation, a non-zero anisotropy produces an overall shift on the inferred γ, the amplitude of which depends on the average value of β_ani of the population of galaxies.

However, we also expect a bias between the distribution in anisotropy of the SLACS lenses and that of the parent sample, similar to the bias in γ that we uncovered in this paper. This is because anisotropy affects the observed velocity dispersion, which enters the lens selection probability by means of the estimated Einstein radius ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$. Figure A.1 illustrates this effect, by showing how the value of ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ varies as a function of both γ and β_ani, for a lens with fixed true value of the Einstein radius and a spatially constant anisotropy. Qualitatively speaking, ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ increases both with γ and β. Coupled with the fact that lenses with a larger ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ were more likely to be selected, this means that galaxies with more radially anisotropic orbits are likely over-represented in the SLACS sample. For a quantitative assessment of the amplitude of the bias, knowledge of the intrinsic distribution in β_ani of the general population is needed. That quantity, however, is not known very well.

Fig. A.1. Estimated Einstein radius of a power-law lens with θE = 1″, as a function of density slope and orbital anisotropy. The lens is at zg = 0.2, with a half-light radius Re = 7 kpc, and the Einstein radius is computed with respect to a source at zs = 0.6. White lines connect points of constant ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$. These are also points of constant velocity dispersion. In general, ${\theta }_{\mathrm{E}}^{(\mathrm{est})}<{\theta }_{\mathrm{E}}$: this is because the line-of-sight velocity dispersion, which is used to estimate ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$, is smaller than the three-dimensional velocity dispersion of a singular isothermal sphere, due to projection effects.

All Tables

All Figures

Fig. 1. SLACS lenses and parent population. Distribution in redshift, observed stellar mass, half-light radius, and stellar velocity dispersion. The green contours and histograms are obtained by weighting the parent sample distribution by σap4. Vertical dashed or dotted lines indicate the median of the distributions of the corresponding colour. Contours enclose 68% and 95% of the distribution. Half-light radius measurements of the parent population were not used in our analysis. Instead, we assumed a stellar mass-size relation from Hyde & Bernardi (2009), which is shown as a dashed red line. The black line in the redshift-stellar mass panel is the inferred truncation stellar mass mt, introduced in section 4.2.

In the text

Fig. 2. Strong lensing cross-section of a circular power-law lens and a point source. A lensing event is defined as the detection of an emission line within an aperture of 1.5″, and of at least two images in photometric data. Top: σSL of a lens with θE = 1.2″ as a function of the slope γ. The thick green line corresponds to a source with intrinsic emission line flux equal to one third of the flux required for a spectroscopic detection, and a broadband flux equal to that needed for a photometric detection. Different lines were obtained by rescaling both the emission line and broadband flux by the amount indicated in the legend. Bottom: geometric factor of the lensing cross-section, that is the value of σSL for the reference source (corresponding to the green line in the top panel), as a function of both θE and γ.

In the text

	Fig. 3. Constraints on m5 and γ, and predicted values of the velocity dispersion, for individual lenses. Each line covers the 68% credible region of the posterior probability. Lines are colour-coded by the physical size of the Einstein radius. Dots mark the maximum posterior probability values.
In the text

Fig. 4. Posterior probability distribution of the model. Parameters shown: average log M5 and γ at fixed stellar mass and half-light radius, and intrinsic scatter around the average. Filled contours: model describing SLACS lenses, without correction for selection effects (see Section 4.8). Thick solid curves: fiducial model. Thin solid curves: model with no fundamental plane prior of Equation (37). Contour levels correspond to 68% and 95% enclosed probability.

In the text

	Fig. 5. Same as Figure 4, but showing parameters describing the scaling of m5 and γ with stellar mass and size. Dashed lines indicate values of the parameters corresponding to the case in which galaxies are homologous systems with an average density slope of γ = 2.00.
In the text

	Fig. 6. Same as Figure 4, but showing parameters describing the lens finding probability, defined in Equation (35), and the average density slope parameter.
In the text

	Fig. 7. Posterior predictive tests of goodness of fit. Each panel shows the posterior predicted (histogram) and observed value (dashed line) of the test quantity indicated at the bottom. Units of Einstein radius are in arcsec and of the velocity dispersion are in km s−1.
In the text

Fig. 8. Posterior predicted trends with stellar mass, for galaxies of average size. Top: average projected mass enclosed within 5 kpc. Middle: average density slope. Bottom: average stellar velocity dispersion. Each band marks the 68% credible region. The magenta band shows the prediction for the parent population. The cyan band shows the prediction from the lens-only analysis, without accounting for selection effects. Dashed blue lines show the prediction for the population of SLACS lenses. Solid black lines show the prediction for the population of detectable lenses (obtained by removing the lens finding probability Pfind from the posterior prediction). Error bars are the values measured for the SLACS lenses, assuming a flat prior on both m5 and γ. Uncertainties on the stellar mass are uncorrelated with those on m5, γ and σap. In the bottom panel, the green line shows our quadratic fit to the m* − σap relation of the parent sample.

In the text

Fig. 9. Estimated Einstein radius, density slope and lens finding probability. Top: 68% probable values of ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ and density slope γ of the SLACS lenses, where ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ is defined in Equation (2). Bottom: 68% credible region of the lens finding probability, defined in Equation (35). The correlation between ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$ and γ makes so that lenses with overestimated γ are more likely to be included in the sample.

In the text

	Fig. 10. Bias on velocity dispersion and density slope of SLACS-like lenses. Top: error on γ as a function of ${\theta }_{\mathrm{E}}^{(\mathrm{SIS})}$, defined in Equation (43), of posterior predicted samples of lenses, marginalised over the values of the model parameters η. Bottom: relative error on σap. The solid line indicates the average in bins of ${\theta }_{\mathrm{E}}^{(\mathrm{SIS})}$.
In the text

Fig. 11. Posterior predicted mass-size and fundamental plane relation. Top: distribution in stellar mass and half-light radius of the parent population (magenta regions) and the SLACS lenses (dashed contours), marginalised over the model parameters η. The dotted black line is the mass-size relation from Hyde & Bernardi (2009). The solid blue line indicates the average size of SLACS lenses as a function of stellar mass. Bottom: distribution in velocity dispersion, as a function of the predicted average given the stellar mass and size (the mean of the distribution in Equation (36)). The solid blue line indicates the average velocity dispersion of the SLACS lenses.

In the text

Fig. 12. Posterior predicted trends with stellar velocity dispersion. Top: average projected dark matter mass within 5 kpc. Bottom: average density slope. The magenta band marks the 68% credible region for the parent population of galaxies. The dashed blue lines enclose the 68% credible region for the SLACS lenses. Grey lines in the top two panels span the values of the SLACS lenses that are allowed with 68% probability.

In the text

	Fig. 13. SLACS Strong lensing bias bSL for a galaxy at z = 0.2 with Re = 7 kpc, as a function of m5 and γ. White lines connect points of constant velocity dispersion, with σap equal to the value indicated. We computed bSL from Equation (19), while marginalising over the model parameters η.
In the text

Fig. A.1. Estimated Einstein radius of a power-law lens with θE = 1″, as a function of density slope and orbital anisotropy. The lens is at zg = 0.2, with a half-light radius Re = 7 kpc, and the Einstein radius is computed with respect to a source at zs = 0.6. White lines connect points of constant ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$. These are also points of constant velocity dispersion. In general, ${\theta }_{\mathrm{E}}^{(\mathrm{est})}<{\theta }_{\mathrm{E}}$: this is because the line-of-sight velocity dispersion, which is used to estimate ${\theta }_{\mathrm{E}}^{(\mathrm{est})}$, is smaller than the three-dimensional velocity dispersion of a singular isothermal sphere, due to projection effects.

In the text