Press Release
Free Access
Issue
A&A
Volume 643, November 2020
Article Number A165
Number of page(s) 40
Section Cosmology (including clusters of galaxies)
DOI https://doi.org/10.1051/0004-6361/202038861
Published online 20 November 2020

© ESO 2020

1. Introduction

There is a discrepancy in the reported measurements of the Hubble constant from early universe and late universe distance anchors. If confirmed, this discrepancy would have profound consequences and would require new or unaccounted physics to be added to the standard cosmological model. Early universe measurements in this context are primarily calibrated with sound horizon physics. This includes the cosmic microwave background (CMB) observations from Planck with H0 = 67.4 ± 0.5 km s−1 Mpc−1 (Planck Collaboration VI 2020), galaxy clustering and weak lensing measurements of the Dark Energy Survey (DES) data in combination with baryon acoustic oscillations (BAO) and Big Bang nucleosynthesis (BBN) measurements, giving H0 = 67.4 ± 1.2 km s−1 Mpc−1 (Abbott et al. 2018), and using the full-shape BAO analysis in the BOSS survey in combination with BBN, giving H0 = 68.4 ± 1.1 km s−1 Mpc−1 (Philcox et al. 2020). All of these measurements provide a self-consistent picture of the growth and scales of structure in the Universe within the standard cosmological model with a cosmological constant, Λ, and cold dark matter (ΛCDM).

Late universe distance anchors consist of multiple different methods and underlying physical calibrators. The most well established one is the local distance ladder, effectively based on radar observations on the Solar system scale, the parallax method, and a luminous calibrator to reach the Hubble flow scale. The SH0ES team, using the distance ladder method with supernovae (SNe) of type Ia and Cepheids, reports a measurement of H0 = 74.0 ± 1.4 km s−1 Mpc−1 (Riess et al. 2019). The Carnegie–Chicago Hubble Project (CCHP) using the distance ladder method with SNe Ia and the tip of the red giant branch measures H0 = 69.6 ± 1.9 km s−1 Mpc−1 (Freedman et al. 2019, 2020). Huang et al. (2020) used the distance ladder method with SNe Ia and Mira variable stars and measured H0 = 73.3 ± 4.0 km s−1 Mpc−1.

Among the measurements that are independent of the distance ladder are the Megamaser Cosmology Project (MCP), which uses water megamasers to measure H0 = 73.9 ± 3.0 km s−1 Mpc−1 (Pesce et al. 2020), gravitational wave standard sirens with km s−1 Mpc−1 (Abbott et al. 2017) and the TDCOSMO collaboration1 (formed by members of H0LiCOW, STRIDES, COSMOGRAIL and SHARP), using time-delay cosmography with lensed quasars (Wong et al. 2020; Shajib et al. 2020a; Millon et al. 2020). Time-delay cosmography (Refsdal 1964) provides a one-step inference of absolute distances on cosmological scales – and thus the Hubble constant. Over the past two decades, extensive and dedicated efforts have transformed time-delay cosmography from a theoretical idea to a contender for precision cosmology (Vanderriest et al. 1989; Keeton & Kochanek 1997; Schechter et al. 1997; Kochanek 2003; Koopmans et al. 2003; Saha et al. 2006; Read et al. 2007; Oguri 2007; Coles 2008; Vuissoz et al. 2008; Suyu et al. 2010, 2013, 2014; Fadely et al. 2010; Sereno & Paraficz 2014; Rathna et al. 2015; Birrer et al. 2016, 2019; Wong et al. 2017; Rusu et al. 2020; Chen et al. 2019; Shajib et al. 2020a).

The keys to precision time-delay cosmography are: Firstly, precise and accurate measurements of relative arrival time delays of multiple images; Secondly, understanding of the large-scale distortion of the angular diameter distances along the line of sight; and thirdly, accurate model of the mass distribution within the main deflector galaxy. The first problem has been solved by high cadence and high precision photometric monitoring, often with dedicated telescopes (e.g., Fassnacht et al. 2002; Tewes et al. 2013; Courbin et al. 2018). The time delay measurement procedure has been validated via simulations by the Time Delay Challenge (TDC1; Dobler et al. 2015; Liao et al. 2015). The second issue has been addressed by statistically correcting the effect of the line of sights to strong gravitational lenses by comparison with cosmological numerical simulations (e.g., Fassnacht et al. 2011; Suyu et al. 2013; Greene et al. 2013; Collett et al. 2013). Millon et al. (2020) recently showed that residuals from the line of sight correction based on this methodology are smaller than the current overall errors. Progress on the third issue has been achieved by analyzing high quality images of the host galaxy of the lensed quasars with provided spatially resolved information that can be used to constrain lens models (e.g., Suyu et al. 2009). By modeling extended sources with complex and flexible source surface brightness instead of just the quasar images positions and fluxes, modelers have been able to move from extremely simplified models like singular isothermal ellipsoids (Kormann et al. 1994; Schechter et al. 1997) to more flexible ones like power laws or stars plus standard dark matter halos (Navarro et al. 1997, hereafter NFW). The choice of elliptical power-law and stars plus NFW profiles was motivated by their generally good description of stellar kinematics and X-ray data in the local Universe. It was validated post-facto by the small residual corrections found via pixellated models (Suyu et al. 2009), and by the overall goodness of fit they provided to the data.

Building on the advances in the past two decades, the H0LiCOW and SHARP collaborations analyzed six individual lenses (Suyu et al. 2010, 2014; Wong et al. 2017; Birrer et al. 2019; Rusu et al. 2020; Chen et al. 2019) and measured H0 for each lens to a precision in the range 4.3−9.1%. The STRIDES collaboration measured H0 to 3.9% from one single quadruply lensed quasar (Shajib et al. 2020a). The seven measurements follow an approximately standard (although evolving over time) procedure (see e.g. Suyu et al. 2017) and incorporate single-aperture stellar kinematics measurements for each lens. The H0LiCOW collaboration combined their six quasar lenses, of which five had their analysis blinded, assuming uncorrelated individual distance posteriors and arrived at km s−1 Mpc−1, a 2.4% measurement of H0 (Wong et al. 2020). Adding the blind measurement by Shajib et al. (2020a) further increases the precision to ∼2% (Millon et al. 2020).

Given the importance of the Hubble tension, it is crucial, however, to continue to investigate potential causes of systematic errors in time-delay cosmography. After all, extraordinary claims, like physics beyond ΛCDM, require extraordinary evidence.

The first and main source of residual modeling error in time-delay cosmography is due to the mass-sheet transform (MST; Falco et al. 1985). MST is a mathematical degeneracy that leaves the lensing observables unchanged, while rescaling the absolute time delay, and thus the inferred H0. This degeneracy is well known and frequently discussed in the literature (e.g., Gorenstein et al. 1988; Kochanek 2002, 2006, 2020a; Saha & Williams 2006; Read et al. 2007; Schneider & Sluse 2013, 2014; Coles et al. 2014; Xu et al. 2016; Birrer et al. 2016; Unruh et al. 2017; Sonnenfeld 2018; Wertz et al. 2018; Blum et al. 2020). Lensing-independent tracers of the gravitational potential of the deflector galaxy, such as stellar kinematics, can break this inherent degeneracy (e.g., Grogin & Narayan 1996; Romanowsky & Kochanek 1999; Treu & Koopmans 2002). Another way to break the degeneracy is to make assumptions on the mass density profile, which is primarily the strategy adopted by the H0LiCOW/STRIDES collaboration (Millon et al. 2020). Millon et al. (2020) showed that the two classes of radial mass profiles considered by the collaboration, power-law and stars and a Navarro Frenk & White (NFW, Navarro et al. 1997) dark matter halo, yield consistent results2. Sonnenfeld (2018), Kochanek (2020a, b) argued that the error budget of individual lenses obtained under the assumptions of power-law or stars + NFW are underestimated and that, given the MST, the typical uncertainty of the kinematic data does not allow one to constrain the mass profiles to a few percent precision3.

A second potential source of uncertainty in the combined TDCOSMO analysis is the assumption of no correlation between the errors of each individual lens system. The TDCOSMO analysis shows that the scatter between systems is consistent with the estimated errors, and the random measurement errors of the observables are indeed uncorrelated (Wong et al. 2020; Millon et al. 2020). However, correlations could be introduced by the modeling procedure and assumptions made, such as the form and prior on the mass profile and the distribution of stellar anisotropies in elliptical galaxies.

In this paper we address these two dominant sources of potential residual uncertainties by introducing a Bayesian hierarchical framework to analyze and interpret the data. Addressing these uncertainties is a major step forward in the field, however it should be noted that the scope of this framework is broader than just these two issues. Its longer term goal is to take advantage of the expanding quality and quantity of data to trade theoretical assumptions for empirical constraints. Specifically, this framework is designed to meet the following requirements: (1) Theoretical assumptions should be explicit and, whenever possible, verified by data or replaced by empirical constraints; (2) Kinematic assumptions and priors must be justified by the data or the laws of physics; (3) The methodology must be validated with realistic simulations. By using this framework we present an updated measurement of the Hubble constant from time-delay cosmography and we lay out a roadmap for further improvements of the methodology to enable a measurement of the Hubble constant from strong lensing time-delay measurements with 1% precision and accuracy.

In practice, we adopt a parameterization that allows us to quantify the full extent MST in our analysis, addressing point (1) listed above. We discuss the assumptions on the kinematic modeling and the impact of the priors chosen. We deliberately choose an uninformative prior, addressing point (2). We make use of a blind submission to the time-delay lens modeling challenge (TDLMC; Ding et al. 2018, 2020) and validate our approach end to end, including imaging analysis, kinematics analysis and MST mitigation, addressing point (3)4.

In our new analysis scheme, the MST is exclusively constrained by the kinematic information of the deflector galaxies, and thus fully accounted for in the error budget. Under these minimal assumptions, we expect that the data currently available for the individual lenses in our TDCOSMO sample will not constrain H0 to the 2% level. In addition, we take into account covariances between the sample galaxies, by formulating the priors on the stellar anisotropy distribution and the MST at the population level and globally sampling and marginalizing over their uncertainties.

To further improve the constraints on the mass profile and the MST on the population level, we incorporate a sample of 33 lenses from the Sloan Lens ACS (SLACS) survey (Bolton et al. 2006) into our analysis. We make use of the lens model inference results presented by Shajib et al. (2020b), which follow the standards of the TDCOSMO collaboration. We assess the assumptions in the kinematics modeling and incorporate integral field unit (IFU) spectroscopy from VIMOS 2D data of a subset of the SLACS lenses from Czoske et al. (2012) in our analysis. This dataset allows us to improve constraints on the stellar anisotropy distribution in massive elliptical galaxies at the population level and thus reduces uncertainties in the interpretation of the kinematic measurements, hence improving the constraints on the MST and H0. Our joint hierarchical analysis is based on the assumption that the massive elliptical galaxies acting as lenses in the SLACS and the TDCOSMO sample represent the same underlying parent population in regard of their mass profiles and kinematic properties. The final H0 value derived in this work is inferred from the joint hierarchical analysis of the SLACS and TDCOSMO samples.

The paper is structured as follows: Sect. 2 revisits the analysis performed on individual lenses and assesses potential systematics due to MST and mass profile assumptions. Section 3 describes the hierarchical Bayesian analysis framework to mitigate assumptions and priors associated to the MST to a sample of lenses. We first validate this approach in Sect. 4 on the TDLMC data set (Ding et al. 2018) and then move to perform this very same analysis on the TDCOSMO data set in Sect. 5. Next, we perform our hierarchical analysis on the SLACS sample with imaging and kinematics data to further constrain uncertainties in the mass profiles and the kinematic behavior of the stellar anisotropy in Sect. 6. We present the joint analysis and final inference on the Hubble constant in Sect. 7. We discuss the limitations of the current work and lay out the path forward in Sect. 8 and finally conclude in Sect. 9.

All the software used in this analysis is open source and we share the analysis scripts and pipeline with the community5. Numerical tests on the impact of the MST are performed with LENSTRONOMY6 (Birrer & Amara 2018; Birrer et al. 2015). The kinematics is modeled with the LENSTRONOMY.GALKIN module. The reanalysis of the SLACS lenses imaging data is performed with DOLPHIN7, a wrapper around LENSTRONOMY for automated lens modeling (Shajib et al. 2020b) and we introduce HIERARC8 (this work) for the hierarchical sampling in conjunction with LENSTRONOMY. All components of the analysis – including analysis scripts and software – were reviewed internally by people not previously involved in the analysis of the sample before the joint inference was performed. All uncertainties stated are given in 16th, 50th and 84th percentiles. Error contours in plots represent 68th and 95th credible regions.

As in previous work by our team – in order to avoid experimenter bias – we keep our analysis blind by using previously blinded analysis products, and all additional choices made in this analysis, such as considering model parameterization and including or excluding of data, are assessed blindly in regard to H0 or parameters directly related to it. All sections, except Sect. 8.5, of this paper have been written and frozen before the unblinding of the results.

2. Cosmography from individual lenses and the mass-sheet degeneracy

In this section we review the principles of time-delay cosmography and the underlying observables (Sect. 2.1 for lensing and time delays and Sect. 2.2 for the kinematic observables). We emphasize how an MST affects the observables and thus the inference of cosmographic quantities (Sect. 2.3). We separate the physical origin of the MST into the line-of-sight (external MST, Sect. 2.4) and mass-profile contributions (internal MST, Sect. 2.5) and then provide the limits on the internal mass profile constraints from imaging data and plausibility arguments in Sect. 2.6. We provide concluding remarks on the constraining power of individual lenses for time-delay cosmography in Sect. 2.7.

2.1. Cosmography with strong lenses

In this section we state the relevant governing physical principles and observables in terms of imaging, time delays, and stellar kinematics. The phenomena of gravitational lensing can be described by the lens equation, which maps the source plane β to the image plane θ (2D vectors on the plane of the sky)

(1)

where α is the angular shift on the sky between the original unlensed and the lensed observed position of an object.

For a single lensing plane, the lens equation can be expressed in terms of the physical deflection angle as

(2)

with Ds, Dds is the angular diameter distance from the observer to the source and from the deflector to the source, respectively. In the single lens plane regime we can introduce the lensing potential ψ such that

(3)

and the lensing convergence as

(4)

The relative arrival time between two images θA and θB, ΔtAB, originated from the same source is

(5)

where c is the speed of light,

(6)

is the Fermat potential (Schneider 1985; Blandford & Narayan 1986), and

(7)

is the time-delay distance (Refsdal 1964; Schneider et al. 1992; Suyu et al. 2010); Dd, Ds, and Dds are the angular diameter distances from the observer to the deflector, the observer to the source, and from the deflector to the source, respectively.

Provided constraints on the lensing potential, a measured time delay allows us to constrain the time-delay distance DΔt from Eq. (5):

(8)

The Hubble constant is inversely proportional to the absolute scales of the Universe and thus scales with DΔt as

(9)

2.2. Deflector velocity dispersion

The line-of-sight projected stellar velocity dispersion of the deflector galaxy, σP, can provide a dynamical mass estimate of the deflector independent of the lensing observables and joint lensing and dynamical mass estimates have been used to constrain galaxy mass profiles (Grogin & Narayan 1996; Romanowsky & Kochanek 1999; Treu & Koopmans 2002).

The modeling of the kinematic observables in lensing galaxies range in complexity from spherical Jeans modeling to Schwarzschild (Schwarzschild 1979) methods. For example, Barnabè & Koopmans (2007), Barnabè et al. (2009) use axisymmetric modeling of the phase-space distribution function with a two-integral Schwarzschild method by Cretton et al. (1999), Verolme & de Zeeuw (2002). In this work, the kinematics and their interpretation are a key component of the inference scheme and thus we provide the reader with a detailed background and the specific assumptions in the modeling we apply.

The dynamics of stars with the density distribution ρ*(r) in a gravitational potential Φ(r) follows the Jeans equation. In this work, we assume spherical symmetry and no rotation in the Jeans modeling. In the limit of a relaxed (vanishing time derivatives) and spherically symmetric system, with the only distinction between radial, , and tangential, , dispersions, the Jeans equation results in (e.g., Binney & Tremaine 2008)

(10)

with the stellar anisotropy parameterized as

(11)

The solution of Eq. (10) can be formally expressed as (e.g., van der Marel 1994)

(12)

where M(r) is the mass enclosed in a three-dimensional sphere with radius r and

(13)

is the integration factor of the Jeans Equation (Eq. (10)). The modeled luminosity-weighted projected velocity dispersion σs is given by (Binney & Mamon 1982)

(14)

where R is the projected radius and Σ*(R) is the projected stellar density

(15)

The observational conditions have to be taken into account when comparing a model prediction with a data set. In particular, the aperture 𝒜 and the PSF convolution of the seeing, 𝒫, need to be folded in the modeling. The luminosity-weighted line of sight velocity dispersion within an aperture, 𝒜, is then (e.g., Treu & Koopmans 2004; Suyu et al. 2010)

(16)

where is taken from Eq. (14).

The prediction of the stellar kinematics requires a three-dimensional stellar density ρ*(r) and mass M(r) profile. In terms of imaging data, we can extract information about the parameters of the lens mass surface density with parameters ξmass and the surface brightness of the deflector with parameters ξlight. When assuming a constant mass-to-light ratio across the galaxy, the integrals in the Jeans equation can be performed on the light distribution and Σ*(R) can be taken to be the surface brightness I(R). To evaluate the three-dimensional distributions, we rely on assumptions on the de-projection to the three-dimensional mass and light components. In this work, we use spherically symmetric models with analytical projections/de-projections to solve the Jeans equation.

An additional ingredient in the calculation of the velocity dispersion is the anisotropy distribution of the stellar orbits, βani(r). It is impossible to disentangle the anisotropy in the velocity distribution and the gravitational potential from velocity dispersion and rotation measurements alone. This is known as the mass-anisotropy degeneracy (Binney & Mamon 1982).

Finally, the predicted velocity dispersion requires angular diameter distances from a background cosmology. Specifically, the prediction of any σP from any model can be decomposed into a cosmological-dependent and cosmology-independent part, as (Birrer et al. 2016, 2019)

(17)

where J(ξmass, ξlight, βani) is the dimensionless and cosmology-independent term of the Jeans equation only relying on the angular units in the light, mass and anisotropy model. The term ξlight in Eq. (17) includes the deflector light contribution. The deflector light is required for the Jeans modeling (Σ* and deconvolved ρ* terms in the equations above). In practice, the inference of the deflector light profile is jointly fit with other light components, such as source light and quasar flux.

Inverting Eq. (17) illustrates that a measured velocity dispersion, σP, allows us to constrain the distance ratio Ds/Dds, independent of the cosmological model and time delays but while relying on the same lens model, ξlens,

(18)

We note that the distance ratio Ds/Dds can be constrained without time delays being available. If one has kinematic and time-delay data, instead of expressing constraints on Ds/Dds, one can also express the cosmologically independent constraints in terms of Dd (e.g., Paraficz & Hjorth 2009; Jee et al. 2015; Birrer et al. 2019) as

(19)

In this work, we do not transform the kinematics constraints into Ds/Dds or Dd constraints but work directly on the likelihood level of the velocity dispersion when discriminating between different cosmological models.

In Appendix B we illustrate the radial dependence on the model predicted velocity dispersion, σP, for different stellar anisotropy models. Observations at different projected radii can partially break the mass-anisotropy degeneracy provided that we have independent mass profile estimates from lensing observables.

2.3. Mass-sheet transform

The MST is a multiplicative transform of the lens equation (Eq. (1)) (Falco et al. 1985)

(20)

which preserves image positions (and any higher order relative differentials of the lens equation) under a linear source displacement β → λβ. The term (1 − λ)θ in Eq. (20) above describes an infinite sheet of convergence (or mass), and hence the name mass-sheet transform. Only observables related to the absolute source size, intrinsic magnification or to the lensing potential are able to break this degeneracy.

The convergence field transforms according to

(21)

The same relative lensing observables can result if the mass profile is scaled by the factor λ with the addition of a sheet of convergence (or mass) of κ(θ) = (1 − λ).

The different observables described in Sects. 2.1 and 2.2 transform by an MST term λ as follow: The image positions remain invariant

(22)

The source position scales with λ

(23)

The time delay scales with λ

(24)

and the velocity dispersion scales with λ as

(25)

Until now we have only stated how the MST impacts observables directly. However, it is also useful to describe how cosmographic constraints derived from a set of observables and assumptions on the mass profile are transformed when transforming the lens model with an MST (Eq. (8), (18), (19)). The time-delay distance (Eq. (7)) is dependent on the time delay Δt (Eq. (5))

(26)

The distance ratio constrained by the kinematics and the lens model scales as

(27)

Given time-delay and kinematics data the inference on the angular diameter distance to the lens is invariant under the MST

(28)

The Hubble constant, when inferred from the time-delay distance, DΔt, transforms as (from Eq. (9))

(29)

Mathematically, all the MSTs can be equivalently stated as a change in the angular diameter distance to the source

(30)

In other words, if one knows the dependence of any lensing variable upon Ds one can transform it under the MST and scale all other quantities in the same way.

2.4. Line-of-sight contribution

Structure along the line of sight of lenses induce distortions and focusing (or de-focusing) of the light rays. The first-order shear distortions do have an observable imprint on the shape of Einstein rings and can thus be constrained as part of the modeling procedure of strong lensing imaging data. The first order convergence effect alters the angular diameter distances along the specific line of sight of the strong lens. We define Dlens as the specific angular diameter distance along the line of sight of the lens and Dbkg as the angular diameter distance from the homogeneous background metric without any perturbative contributions. Dlens and Dbkg are related through the convergence terms as

(31)

κs is the integrated convergence along the line of sight passing through the strong lens to the source plane and the term 1 − κs corresponds to an MST (Eq. (30))9. To predict the velocity dispersion of the deflector (Eq. (17)), the terms κs and κds are relevant when using background metric predictions from a cosmological model (Dbkg). To predict the time delays (Eq. (5)) from a cosmological model, all three terms are relevant. We can define a single effective convergence, κext, that transforms the time-delay distance (Eq. (7))

(32)

with

(33)

2.5. External vs. internal mass sheet transform

An MST (Eq. (21)) is always linked to a specific choice of lens model and so is its physical interpretation. The MST can be either associated with line-of-sight structure (κs) not affiliated with the main deflector or as a transform of the mass profile of the main deflector itself (e.g., Koopmans 2004; Saha & Williams 2006; Schneider & Sluse 2013; Birrer et al. 2016; Shajib et al. 2020a).

There are different observables and physical priors related to these two distinct physical causes and we use the notation κs to describe the external convergence aspect of the MST and λint to describe the internal profile aspect of the MST. The total transform which affects the time delays and kinematics (see Eqs. (24) and (25)) is the product of the two transforms

(34)

The line-of-sight contribution can be estimated by tracers of the larger scale structure, either using galaxy number counts (e.g., Rusu et al. 2017) or weak lensing of distant galaxies by all the mass along the line of sight (e.g., Tihhonova et al. 2018), and can be estimated with a few per cent precision per lens. The internal MST requires either priors on the form of the deflector profile or exquisite kinematic tracers of the gravitational potential. The λint component is the focus of this work.

2.6. Approximate internal mass-sheet transform

Imposing the physical boundary condition, limr → ∞κ(r) = 0, violates the mathematical form of the MST10. However, approximate MSTs that satisfy the boundary condition of a finite physically enclosed mass may still be possible and encompass the limitations and concerns of strong gravitational lensing in providing precise constraints on the Hubble constant. We specify an approximate MST as a profile without significantly impacting imaging observables around the Einstein radius and resulting in the transforms of the time delays (Eq. (24)) and kinematics (Eq. (25)).

Cored mass components, κc(r), can serve as physically motivated approximations to the MST (Blum et al. 2020). We can write a physically motivated approximate internal MST with a parameter λc as

(35)

where κmodel corresponds to the model used in the reconstruction of the imaging data and λc describes the scaling between the cored and the other model components, in resemblance to λint. Approximating a physical cored transform with the pure MST means that:

(36)

in deriving all the observable scalings in Sect. 2.3.

Blum et al. (2020) showed that several well-chosen cored 3D mass profiles, ρ(r), can lead to approximate MST’s in projection, κc(r), with physical interpretations, such as

(37)

resulting in the projected convergence profile

(38)

where Σcrit is the critical surface density of the lens. The specific functional form of the profile listed above (37) resemble the outer slope of the NFW profile with ρ(r) ∝ r−3.

Figure 1 illustrates a composite profile consisting of a stellar component (Hernquist profile) and a dark matter component (NFW + cored component, Eq. (37)) which transform according to an approximate MST. The stellar component gets rescaled by the MST while the cored component is transforming only the dark matter component.

thumbnail Fig. 1.

Illustration of a composite profile consisting of a stellar component (Hernquist profile, dotted lines) and a dark matter component (NFW + cored component (Eq. (38)), dashed lines) which transform according to an approximate MST (joint as solid lines). The stellar component gets rescaled by the MST while the cored component transforms the dark matter component. Left: profile components in three dimensions. Right: profile components in projection. The transforms presented here cannot be distinguished by imaging data alone and require i.e., stellar kinematics constraints (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_composite_cored.ipynb).

Open with DEXTER

It is of greatest importance to quantify the physical plausibility of those transforms and their impact on other observables in detail. In this section we extend the study of Blum et al. (2020). We perform detailed numerical experiments on mock imaging data to quantify the constraints from imaging data, time delays and kinematics, and we quantify the range of such an approximate transform with physically motivated boundary conditions. Further illustrations and details on the examples given in this section can be found in Appendix A.

2.6.1. Imaging constraints on the internal MST

In this section we investigate the extent to which imaging data is able to distinguish between different lens models with different cored mass components and their impact on the inferred time delay distance in combination with time delay information. We first generate a mock image and time delays without a cored component and then perform the inference with an additional cored component model (Eq. (38)) parameterized with the core radius Rc and the core projected density Σc ≡ (1 − λc) (Eq. (35)). In our specific example, we simulate a quadruply lensed quasar image similar to Millon et al. (2020) (more details in Appendix A and Fig. A.2) with a power-law elliptical mass distribution (PEMD, Kormann et al. 1994; Barkana 1998)

(39)

where γpl is the logarithmic slope of the profile, qm is the axis ratio of the minor and the major axes of the elliptical profile, and θE is the Einstein radius. The coordinate system is defined such that θ1 and θ2 are along the major and minor axis respectively. We also add an external shear model component with distortion amplitude γext and direction ϕext. The PEMD+shear model is one of two lens models considered in the analysis of the TDCOSMO sample. For the source and lens galaxies we use elliptical Sérsic surface brightness profiles. We add a Gaussian point spread function (PSF) with full-width-at-half-maximum (FWHM) of 01, pixel scale of 005 and noise properties consistent with the current TDCOSMO sample of Hubble Space Telescope (HST) images. The time delays between the images between the first arriving image and the subsequent images are 11.7, 27.6, and 94.0 days, respectively. We chose time-delay uncertainties of ±2 days between the three relative delays. The time-delay precision does not impact our conclusions about the MST. The inference is performed on the pixel level of the mock image as with the real data on the TDCOSMO sample.

In the modeling and parameter inference, we add an additional cored mass component (Eq. (38)) and perform the inference on all the lens and source parameters simultaneously, including the core radius Rc and the projected core density Σc. In the limit of a perfect MST there is a mathematical degeneracy if we only use the imaging data as constraints. We thus expect a full covariance in the parameters involved in the MST (Einstein radius of the main deflector, source position, source size etc.) and the posterior inference of our problem to be inefficient in the regime where the cored profile mimics the full MST (κc(θ) acts as Σcrit for Rc → ∞). To improve the sampling, instead of modeling the cored profile κc(θ), we model the difference between the cored component and a perfect MST, Δκc = κc(θ)−Σcrit, with λc (Eq. (35)) instead. Δκc is effectively the component of the model that does not transform under the MST and leads to a physical three-dimensional profile interpretation.

Figure 2 shows the inference on the relevant lens model parameters for the mock image described in Appendix A. The input parameters are marked as orange lines for the model without a cored component. We can clearly see that for small core radii, Rc, the approximate MST parameter λc can be constrained. This is the limit where the additional core profile cannot mimic a pure MST at a level where the data is able to distinguish between them. For core radii Rc = 3θE, the uncertainty on the approximate MST, λc, is 10%. For core radii Rc >  5θE, the approximate MST is very close to the pure MST and the imaging information in our example is not able to constrain λc to better than λc ± 0.4. We make use of the expected constraining power on λc as a function of Rc when we discuss the plausibility of certain transforms. When looking at the inferred time-delay distance λcDΔt, we see that this quantity is constant as a function of Rc and thus the time-delay prediction is accurately being transformed by a pure MST (Eq. (24)). Overall, we find that λc ≈ λint is valid for larger core radii.

thumbnail Fig. 2.

Illustration of the constraining power of imaging data on a cored mass component (Eq. (35)). Shown are the parameter inference of the power-law profile mock quadruply lensed quasar of Fig. A.2 when including a marginalization of an additional cored power law profile (Eq. (38)). Orange lines indicate the input truth of the model without a cored component. λc is the scaled core model parameter (Eq. (35)) resembling the pure MST for large core radii (λc ≈ λint) (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER

Identical tests with a composite profile instead of a PEMD profile result in the same conclusions and are available online11.

2.6.2. Allowed cored mass components from physical boundary conditions

In the previous Sect. 2.6.1 we demonstrated that, for large core radii, there are physical profiles that approximate a pure MST (λc ≈ λint). In this section we take a closer look at the physical interpretation of such large positive and negative cored component transforms with respect to a chosen mass profile. It is possible that the core model itself does not require a physical interpretation as it is overall included in the total mass distribution. The galaxy surface brightness provides constraints on the stellar mass distribution (modulo a mass-to-light conversion factor) and the focus here is a consideration of the distribution of the invisible (dark) matter component of the deflector. Our starting model is a NFW profile and we assess departures from this model by using a cored component.

We apply the following conservative boundary conditions on the distribution of the dark matter component: Firstly, the total mass of the cored component within a three-dimensional radius shall not exceed the total mass of the NFW profile within the same volume, Mcore(<r) ≤ MNFW(<r). This is not a strict bound, but violating this condition would imply changing the mass of the halo itself. Secondly, the density profile shall never drop to negative values, ρNFW + core(r) ≥ 0.

Those two imposed conditions define a physical interpretation of a three-dimensional mass profile as being a redistribution of matter from the dark matter component and a rescaling of the mass-to-light ratio of the luminous component. An independent estimate of the mass-to-light ratio of few per cent is below our current limits of knowledge about the stellar initial mass function, stellar evolution models and dust extinction. Moreover, the mass-to-light ratio can vary with radius. Figure 3 provides the constraints from the two conditions, as well as from the imaging data constraints of Sect. 2.6.1, for an expected NFW mass and concentration profile at a typical lens and source redshift configuration. The remaining white region in Fig. 3 is effectively allowed by the imaging data and simple plausibility considerations. We conclude that the physically allowed parameter space does encompass a pure MST with , with more parameter volume for λint <  1, which corresponds to a positive cored component. We emphasize that the constraining power at small core radii may be due to the angular rather than the radial imprint of the cored profile (see e.g., Kochanek 2020b). However, such a behavior would not alter our conclusions and inference method chosen in the analysis presented in subsequent sections of this work. We also performed this inference for a composite (stellar light + NFW dark matter) model and arrive at the same conclusions.

thumbnail Fig. 3.

Constraints on an approximate internal MST transform with a cored component, λc, of an NFW profile as a function of core radius. In gray are the 1-σ exclusion limits that imaging data can provide. In orange is the region where the total mass of the core within a three-dimensional radius exceeds the mass of the NFW profile in the same sphere. In blue is the region where the transformed profile results in negative convergence at the core radius. The white region is effectively allowed by the imaging data and simple plausibility considerations and where we can use the mathematical MST as an approximation (λc ≈ λint). The halo mass, concentration and the redshift configuration is displayed in the lower left box (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER

2.6.3. Stellar kinematics of an approximate MST

In this section we investigate the kinematics dependence on the approximate MST. To do so, we perform spherical Jeans modeling (Sect. 2.2) and compute the predicted velocity dispersion in an aperture under realistic seeing conditions (Eq. (16)) for models with a cored mass component as an approximation of the MST.

Figure 4 compares the actual predicted kinematics from the modeling of the physical three-dimensional mass distribution κλc (Eq. (35)) and the analytic relation of a perfect MST (Eq. (25)) for the mock lens presented in Appendix A. For this figure, we chose an aperture size of 1″ × 1″ and seeing of FWHM = 07 and an isotropic stellar orbit distribution (βani(r) = 0). For λc in the range [0.8,1.2], the MST approximation in the predicted velocity dispersion is accurate to < 1%. We conclude that, for the λint range considered in this work, the analytic approximation of a perfect MST is valid to reliably compute the predicted velocity dispersion. The precise dependence of the velocity dispersion only marginally depends on the specific core radius Rc and the approximation remains valid for all reasonable and non-excluded core radii and λint. We tested that our conclusions also hold for different anisotropy profiles and observational conditions.

thumbnail Fig. 4.

Comparison of the actual predicted kinematics from the modeling of the physical three-dimensional mass distribution κλint (Eq. (35)) for varying core sizes (solid) and the analytic relation of a perfect MST (Eq. (25), dashed) for the mock lens presented in Fig. A.2. Lower panel: fractional differences between the exact prediction and a perfect MST calculation. The MST prediction matches to < 1% in the considered range. Minor numerical noise is present at the subpercent level (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER

2.7. Constraining power using individual lenses

For each individual strong lens in the TDCOSMO sample, there are four data sets available: (1) imaging data of the strong lensing features and the deflector galaxy, 𝒟img; (2) time-delay measurements between the multiple images, 𝒟td; (3) stellar kinematics measurement of the main deflector galaxy, 𝒟spec; (4) line-of-sight galaxy count and weak lensing statistics, 𝒟los.

These data sets are independent and so are their likelihoods in a joint cosmographic inference. Hence, we can write the likelihood of the joint set of the data 𝒟 = {𝒟img, 𝒟td, 𝒟spec, 𝒟los} given the cosmographic parameters {Dd, Ds, Dds}≡Dd, s, ds as

(40)

(41)

(42)

(43)

In the expression above we only included the relevant model components in the expressions of the individual likelihoods. ξlight formally includes the source and lens light surface brightness. For the time-delay likelihood, we only consider the time-variable source position from the set of ξlight parameters. In Appendix C we provide details on the computation of the combined likelihood, in particular with application in the hierarchical context.

An approximate internal MST of a power law with λint of 10% still leads to physically interpretable mass profiles with the Hubble constant changed by 10% (see Eq. (29)). Imaging data is not sufficiently able to distinguish between models producing H0 value within this 10% range (Kochanek 2020a). The kinematics are changed with good approximation by Eq. (25) through this transform. The kinematic prediction is also cosmology dependent by Eq. (17). The scalings of an MST are analytical in the model-predicted time-delay distance and kinematics and thus its marginalization can be performed in post processing given posteriors for a specific lens model family that breaks the MST, such as a power-law model.

The kinematics information is the decisive factor in discriminating different profile families. The relative uncertainty in the velocity dispersion measurement directly propagates into the relative uncertainty in the MST as

(44)

The current uncertainties on the velocity dispersion measurements, on the order of 5−10% (including the uncertainties due to stellar template mismatch and other systematic errors) limit the precise determination of the mass profile per individual lens. Uncertainties in the interpretation of the stellar anisotropy orbit distribution additionally complicates the problem. Birrer et al. (2016) performed such an analysis and demonstrated that an explicit treatment of the MST (in their approach parameterized as a source scale) leads to uncertainties consistent with the expectations of Kochanek (2020a). Because the kinematic measurement of each lens is not sufficiently precise to constrain the mass profile to the desired level, in this work we marginalize over the uncertainties properly accounting for the priors.

3. Hierarchical Bayesian cosmography

The overarching goal of time-delay cosmography is to provide a robust inference of cosmological parameters, π, and in particular the absolute distance scale, the Hubble constant H0, and possibly other parameters describing the expansion history of the Universe (such as ΩΛ or Ωm), from a sample of gravitational lenses with measured time delays. Based on the conclusions we draw from Sect. 2, it is absolutely necessary to propagate assumptions and priors made on the analysis of an individual lens hierarchically when performing the inference on the cosmological parameters from a population of lenses. In particular, this is relevant for parameters that we cannot sufficiently constrain on a lens-by-lens basis and parameters whose uncertainties significantly propagate to the H0 inference on the population level. In this section, we introduce three specific hierarchical sampling procedures for properties of lensing galaxies and their selection that are relevant for the cosmographic analysis. In particular, these are: (1) an overall internal MST relative to a chosen mass profile, λint, and its distribution among the sample of lenses; (2) stellar anisotropy distribution in the sample of lenses; (3) the line-of-sight structure selection and distribution of the lens sample.

In Sect. 3.1 we formalize the Bayesian problem and define an approximate scheme for the full hierarchical inference that allows us to keep track of key systematic uncertainties while still being able to reuse currently available inference products. In Sect. 3.2 we specify the hyper-parameters we sample on the population level. Section 3.3 details the specific approximations in the likelihood calculation. All hierarchical computations and sampling presented in this work are implemented in the open-source software HIERARC.

3.1. Hierarchical inference problem

In Bayesian language, we want to calculate the probability of the cosmological parameters, π, given the strong lensing data set, p(π|{𝒟i}N), where 𝒟i is the data set of an individual lens (including imaging data, time-delay measurements, kinematic observations and line-of-sight galaxy properties) and N the total number of lenses in the sample.

In addition to π, we introduce ξ that incorporates all the model parameters. Using Bayes rule and considering that the data of each individual lens 𝒟i is independent, we can write:

(45)

In the following, we divide the nuisance parameter, ξ, into a subset of parameters that we constrain independently per lens, ξi, and a set of parameters that require to be sampled across the lens sample population globally, ξpop. The parameters of each individual lens, ξi, include the lens model, source and lens light surface brightness and any other relevant parameter of the model to predict the data. Hence, we can express the hierarchical inference (Eq. (45)) as

(46)

where {ξi}N = {ξ1, ξ2, …, ξN} is the set of the parameters applied to the individual lenses and p(ξi) are the interim priors on the model parameters in the inference of an individual lens. The cosmological parameters π are fully encompassed in the set of angular diameter distances, {Dd, Ds, Dds}≡Dd, s, ds, and thus, instead of stating π in Eq. (46), we now state Dd, s, ds(π). Up to this point, no approximation was applied to the full hierarchical expression (Eq. (45)).

From now on, we assume

(47)

which states that, for the parameters classified as ξ{i}, the interim priors do not propagate into the cosmographic inference and the population prior on those parameters is formally known exactly. The population parameters, ξpop, describe a distribution function such that the values of individual lenses, ξpop, i, follow the distribution likelihood p(ξpop, i|ξpop).

With this approximation and the notation of the sample distribution likelihood, we can simplify expression (46) to

(48)

where

(49)

are the individual likelihoods from an independent sampling of each lens with access to global population parameters, ξpop, and marginalized over the population distribution. The integral in Eq. (49) goes over all individual parameters where a population distribution is applied. Equation (40) is effectively expression (49) without the marginalization over parameters assigned as ξpop.

For parameters in the category ξ{i}, our approximation implies that there is no population prior and that the interim priors do not impact the cosmographic inference. This approximation is valid in the regime where the posterior distribution in ξ{i} is effectively independent of the prior. Although formally this is never true, for many parameters in the modeling of high signal-to-noise imaging data the individual lens modeling parameters are very well constrained relative to the prior imposed.

In the following we highlight some key aspects of the cosmographic analysis and in particular the inference on the Hubble constant where the approximation stated in expression (47) is not valid and thus fall in the category of ξpop. We give explicit parameterizations of these effects and provide specific expressions to allow for an efficient and sufficiently accurate sampling and marginalization, according to Eq. (49), for individual lenses within an ensemble.

3.2. Lens population hyper-parameters

In this section we discuss the choices of population level hyper-parameters we include in our analysis.

3.2.1. Deflector lens model

The deflectors in the quasar lenses with measured time delays of the TDCOSMO sample are massive elliptical galaxies. These galaxies, observationally, follow a tight relation in a luminosity, size and velocity dispersion parameter space (e.g., Faber & Jackson 1976; Auger et al. 2010; Bernardi et al. 2020), exhibiting a high degree of self-similarity among the population.

In Sect. 2.6 we defined λc as the approximate MST relative to a chosen profile of an individual lens and established the close correspondence to a perfect MST (λc ≈ λint). For the inference from a sample of lenses, the sample distribution of deflector profiles is the relevant property to quantify. For the deflector mass profile, we do not want to artificially break the MST based on imaging data and require the kinematics to constrain the mass profile. To do so, we chose as a base-line model a PEMD (Eq. (39)) to be constrained on the lens-by-lens case and we add a global internal MST specified on the population level, λint.

The PEMD lens profile inherently breaks the MST and the parameters of the PEMD profile can be precisely constrained (within few per cent) by exquisite imaging data. In this work, we avoid describing the PEMD parameters at the population level, such as redshift, mass or galaxy environment, and make use of the individual lens inference posterior products derived on flat priors. We note that the power-law slope, γpl, of the PEMD profile inferred from imaging data is a local quantity at the Einstein radius of the deflector. The Einstein radius is a geometrical quantity that depends on the mass of the deflector and lens and source redshift. Thus, the physical location of the measured γpl from imaging data depends on the redshift configuration of the lens system. In a scenario where the mass profiles of massive elliptical galaxies deviate from an MST transformed PEMD resulting in a gradient in the measured slope γpl as a function of physical projected distance, a global joint MST correction on top of the individually inferred PEMD profiles may lead to inaccuracies.

To allow for a radial trend in the applied MST relative to the imaging inferred local quantities, we parameterize the global MST population with a linear relation in reff/θE as

(50)

where λint, 0 is the global MST when the Einstein radius is at the half-light radius of the deflector, reff/θE = 1, and αλ is the linear slope in the expected MST as a function of reff/θE. In this form, we assume self-similarity in the lenses in regard to their half-light radii. In addition to the global MST normalization and trend parameterization, we add a Gaussian distribution scatter with standard deviation σ(λint) at fixed reff/θE.

Wong et al. (2020) and Millon et al. (2020) showed that the TDCOSMO sample results in statistically consistent individual inferences when employing a PEMD lens model. This implies that the global properties of the mass profiles of massive elliptical galaxies in the TDCOSMO sample can be considered to be homogeneous to the level to which the data allows to distinguish differences.

3.2.2. External convergence

The line-of-sight convergence, κext, is a component of the MST (Eq. (34)) and impacts the cosmographic inference. When performing a joint analysis of a sample of lenses, the key quantity to constrain is the sample distribution of the external convergence. We require the global selection function of lenses to be accurately represented to provide a Hubble constant measurement. A bias in the distribution mean of κext on the population level directly leads to a bias of H0.

In this work, we do not explicitly constrain the global external convergence distribution hierarchically but instead constrain p(κext) for each individual lens independently. However, due to the multiplicative nature of internal and external MST (Eq. (34)), the kinematics constrains foremost the total MST, which is the relevant parameter to infer H0. The population distribution of p(κext) only changes the interpretation of the divide into internal vs. external MST and the scatter in each of the two parts.

3.2.3. Stellar anisotropy

The anisotropy distribution of stellar orbits (Eq. (11)) can alter significantly the observed line-of-sight projected stellar velocity dispersion (see Sect. 2.2 and Appendix B). The kinematics can constrain (together with a lens model) the angular diameter distance ratio Ds/Dds (Eqs. (17) and (18)). Having a good quantitative handle on the anisotropy behavior of the lensing galaxies is therefor crucial in allowing for a robust inference of cosmographic quantities. As is the case for an internal MST, the anisotropy cannot be constrained on a lens-by-lens basis with a single aperture velocity dispersion measurement, which impacts the derived cosmographic constraints. It is thus crucial to impose a population prior on the deflectors’ anisotropic stellar orbit distribution and propagate the population uncertainty onto the cosmographic inference.

Observations suggest that typical massive elliptical galaxies are, in their central regions, isotropic or mildly radially anisotropic (e.g., Gerhard et al. 2001; Cappellari et al. 2007); similarly, different theoretical models of galaxy formation predict that elliptical galaxies should have anisotropy varying with radius, from almost isotropic in the center to radially biased in the outskirts (van Albada 1982; Hernquist 1993; Nipoti et al. 2006). A simplified description of the transition can be made with an anisotropy radius parameterization, rani, defining βani as a function of radius r (Osipkov 1979; Merritt 1985)

(51)

To describe the anisotropy distribution on the population level, we explicitly parameterize the profile relative to the measured half-light radius of the galaxy, reff, with the scaled anisotropy parameter

(52)

To account for lens-by-lens differences in the anisotropy configuration, we also introduce a Gaussian scatter in the distribution of aani, parameterized as σ(aani), such that σ(aani)⟨aani⟩ is the standard deviation of aani at sample mean ⟨aani⟩.

3.2.4. Cosmological parameters

All relevant cosmological parameters, π, are part of the hierarchical Bayesian analysis. Wong et al. (2020) and Taubenberger et al. (2019) showed that when adding supernovae of type Ia from the Pantheon (Scolnic et al. 2018) or JLA (Betoule et al. 2014) sample as constraints of an inverse distance ladder, the cosmological-model dependence of strong-lensing H0 measurements is significantly mitigated.

In this work, we assume a flat ΛCDM cosmology with parameters H0 and Ωm. We are using the inference from the Pantheon-only sample of a flat ΛCDM cosmology with Ωm = 0.298 ± 0.022 as our prior on the relative expansion history of the Universe in this work.

3.3. Likelihood calculation

In Sect. 3.1 we presented the generic form of the likelihood ℒ(𝒟i|Dd, s, ds, ξpop) (Eq. (49)) that we need to evaluate for each individual lens for a specific choice of hyper-parameters, and in Sect. 3.2 we provided the specific choices and parameterization of the hyper-parameters used in this work. In this section, we specify the specific likelihood of Eq. (49), ℒ(𝒟i|Dd, s, ds, ξpop, i), that we use, since it is accessible and sufficiently fast to evaluate so that we can sample over a large number of lenses and their population priors.

Specifically, the parameters treated on the population level are ξpop, i = {λint, 0, αλ, σ(λint),⟨aani⟩,σ(aani)}. Our choice of hyper-parameters allows us to reutilize many of the posterior products derived from an independent analysis of single lenses (Eq. (40)). None of the lens model parameters, ξmass, except parameters describing λint and none of the light profile parameters, ξlight, are treated on the population level and thus we can sample those independently for each lens directly from their imaging data

(53)

Furthermore, κext and λint can be merged to a total MST parameter λ according to their definitions (Eq. (34)). All observables and thus the likelihood only respond to this overall MST parameter.

4. Validation on the time-delay lens modeling challenge

Before applying the hierarchical framework to real data, we use the time-delay lens modeling challenge (TDLMC; Ding et al. 2018, 2020) data set to validate the hierarchical analysis and to explore different anisotropy models and priors. The TDLMC was structured with three independent submission rungs. Each of the rungs contained 16 mock lenses with HST-like imaging, time delays and kinematics information. The H0 value used to create the mocks was hidden from the modeling teams. The Rung1 and Rung2 mocks both used PEMD (Eq. (39)) with external shear lens models. The Rung3 lenses were generated by ray-tracing through zoom-in hydrodynamic simulations and reflect a large complexity in their mass profiles and kinematic structure, as expected in the real Universe.

In the blind submissions for Rung1 and Rung2, different teams demonstrated that they could recover the unbiased Hubble constant within their uncertainties under realistic conditions of the data products, uncertainties in the Point Spread Function (PSF) and complex source morphology. In particular, two teams used LENSTRONOMY in their submissions in a completely independent way and achieved precise constraints on H0 while maintaining accuracy. For Rung1 and Rung2, the most precise submissions used the same model parameterization in their inference, thus omitting the problems reviewed in Sect. 2.

It is hard to draw precise conclusions from Rung3 as there are remaining issues in the simulations, such as numerical smoothing scale, sub-grid physics, and a truncation at the virial radius. For more details of the challenge setup we refer to Ding et al. (2018) and on the results and the simulations used in Rung3 to Ding et al. (2020). For a recent study comparing spectroscopic observations with hydrodynamical simulations at z = 0 we refer for instance to van de Sande et al. (2019).

Despite the limitations of the available simulations for accurate cosmology, the application of the hierarchical analysis scheme on TDLMC Rung3 is a stress for the flexibility introduced by the internal MST and the kinematic modeling. Furthermore, the stellar kinematics from the stellar particle orbits provides a self-consistent and highly complex dynamical system. The analysis of TDLMC Rung3 can further help in validating the kinematic modeling aspects in our analysis. However, the removal of substructure in post-processing and truncation effects do not allow, in this regard, conclusions below the 1% level (see Ding et al. 2020). For the effect of substructure on the time delays we refer, for instance, to Mao & Schneider (1998), Keeton & Moustakas (2009) and for a study including the full line-of-sight halo population to Gilman et al. (2020).

We describe the analysis as follow: In Sect. 4.1 we discuss the modeling of the individual lenses. In Sect. 4.2 we describe the hierarchical analysis and priors, and present the inference on H0.

4.1. TDLMC individual lens modeling

For the validation, we make use of the blind submissions of the EPFL team by A. Galan, M. Millon, F. Courbin and V. Bonvin. The modeling of the EPFL team is performed with LENSTRONOMY, including an adaptive PSF reconstruction technique and taking into account astrometric uncertainties explicitly (e.g., Birrer & Treu 2019). Overall, the submissions of the EPFL team follow the standards of the TDCOSMO collaboration. The time that each investigator spent on each lens was substantially reduced due to the homogeneous mock data products, the absence of additional complexity of nearby perturbers and the line of sight, and improvements in the modeling procedure (Shajib et al. 2019). The EPFL team achieved the target precision and accuracy requirement on Rung2, with and without the kinematic constraints, and thus showed reliable inference of lens model parameters within a mass profile parameterization for which the MST does not apply. We refer to the TDLMC paper (Ding et al. 2020) for the details of the performance of all of the participating teams.

We use Rung2 as the reference result for which the MST does not apply, and Rung3 as a test case of the hierarchical analysis. In particular, we make use of the EPFL team’s blind Rung3 submission of the joint time-delay and imaging likelihood (Eq. (C.11)) of their PEMD + external shear models to allow for a direct comparison with the Rung2 results without the kinematics constraints. From the model posteriors of the EPFL team submission, we require the time-delay distance DΔt, Einstein radius θE, power-law slope γpl and half-light radius reff of the deflector. The added external convergence is specified in the challenge setup to be drawn from a normal distribution with mean ⟨κext⟩ = 0 and σ(κext) = 0.025. The EPFL submission of Rung3, which is used in this work, consists of 13 lenses out of the total sample of 16. Three lenses were dropped in their analysis prior to submission due to unsatisfactory results and inconsistency with the submission sample. The uncertainty on the Einstein radius and half-light radius is at subpercent value for all the lenses and the power-law slope reached an absolute precision ranging from below 1% to about 2% for the least constraining lens in their sample from the imaging data alone.

In this work, we perform the kinematic modeling and the likelihood calculation within the hierarchical framework. We use the anisotropy model of Osipkov (1979) and Merritt (1985) (Eq. (51)) with a parameterization of the transition radius relative to the half-light radius (Eq. (52)). We assume a Hernquist light profile with reff in conjunction with the power-law lens model posteriors θE and γpl to model the dimensionless kinematic quantity J (Eqs. (16) and (17)), incorporating the slit mask and seeing conditions (slit 1″ × 1″, seeing FWHM = 06), as specified in the challenge setup.

4.2. TDLMC hierarchical analysis

For the setting of the TDLMC we only sample H0 as a free cosmology-relevant parameter. The matter density Ωm = 0.27 is provided in the challenge setup. We extend the EPFL submission by adding an internal MST distribution with a linear scaling of reff/θE described by λint, 0 and αλ (Eq. (50)) and Gaussian standard deviation σ(λint) of the population at fixed reff/θE. The anisotropy parameter aani is also treated on the population level with mean ⟨aani⟩ and Gaussian standard deviation σ(aani) for the population. In the hierarchical sampling we ignore the covariances between DΔt and the model prediction of the kinematics J. This is justified because of the precise γpl constraints from the imaging data and the inference from the EPFL team.

The summary of the parameters and prior being used in this inference on the TDLMC is presented in Table 1. We chose two different forms of the prior on the anisotropy parameter ⟨aani⟩, one uniform in ⟨aani⟩ and a second one uniform in log(⟨aani⟩), covering the same range in the parameter space, to investigate prior dependences in our inference. To account for the external convergence, we marginalize for each individual lens from the probability distribution p(κext) as specified in the challenge setup12.

Table 1.

Summary of the model parameters sampled in the hierarchical inference on TDLMC Rung3 in Sect. 4.

Figure 5 shows the posteriors of the hierarchical analysis with the priors specified in Table 1.

thumbnail Fig. 5.

Mock data from the TDLMC Rung3 inference with the parameters and prior specified in Table 1. Orange contours indicate the inference with a uniform prior in aani while the purple contours indicate the inference with a uniform priors in log(aani). The thin vertical line indicates the ground truth H0 value in the challenge (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/TDLMC/TDLMC_rung3_inference.ipynb).

Open with DEXTER

We recover the assumed value for the Hubble constant (H0 = 65.413 km s−1 Mpc−1) within the uncertainties of our inference. We find km s−1 Mpc−1 for the 𝒰(log(aani)) prior and km s−1 Mpc−1 for the 𝒰(aani) prior. We note that a uniform prior in log(aani) is a slightly less informative prior than a uniform prior in aani in the same range, as already pointed out by Birrer et al. (2016). In the remaining of this work 𝒰(log(aani)) is the prior of choice in the absence of additional data that constrain the stellar anisotropy of massive elliptical galaxies to provide H0 constraints. The hierarchical analysis and the additional degree of freedom in the mass profile allows us to accurately correct for the insufficient assumptions in the mass profiles on the simulated galaxies. The kinematics modeling indicates that there is more mass in the central part of the galaxies than is modeled with a single power-law profile and infers λint >  1.

We notice a nonzero inferred scatter in the internal MST distribution. One contributing source to this scatter is the fact that the external convergence component was added in post-processing in the TDLMC time delays (Eq. (24)). The rescaling was not applied to the velocity dispersion (Eq. (25)), leading to an artificial scatter in this relation equivalent to the distribution scatter of κext, σ(κext) = 0.025. As the mean in the convergence distribution in the TDLMC is ⟨κext⟩ = 0, we do not expect biases beyond a scatter to occur.

The velocity dispersion measurements allow us to constrain λint and effectively probe a more flexible mass model family. Generally, the velocity dispersion estimates have a 5% relative uncertainty on each individual mock lens. As an ensemble, the 13 lenses of the EPFL submission in the TDLMC Rung3 provide information to infer λint to 2.8% precision (see Eq. (44)) in the limit of a perfect anisotropy model.

The final achieved precision on H0 from the sample of lenses, however, is 8%, dominated by the uncertainty in λint. The fact that, within our chosen priors, the kinematics cannot constrain λint to better than 8% comes from the uncertainty in the anisotropy model. More constraining priors on the anisotropy distribution of the stellar orbits in the lensing galaxies are the key to reducing the uncertainty in the H0 inference (see e.g., Birrer et al. 2016; Shajib et al. 2018; Yıldırım et al. 2020).

5. TDCOSMO mass profile and H0 inference

Having verified the hierarchical approach introduced in Sect. 3 in simultaneously constraining mass profiles and H0 with imaging, kinematics and time-delay observations in the TDLMC (Sect. 4) we employ the inference on the TDCOSMO sample set to measure H0. The inference on the TDCOSMO data is identical to the validation on the TDLMC, apart from some necessary modifications due to the additional complexity in the line-of-sight structure of the real data. In Sect. 5.1 we summarize the data and individual analyses for each single lens of the TDCOSMO sample. In Sect. 5.2 we describe the hierarchical analysis and present the results.

5.1. TDCOSMO sample overview

The analysis presented in this work heavily relies on data and analysis products collected and presented in the literature. We give here a detailed list of the references relevant for our work for the seven lenses of the TDCOSMO sample.

1. B1608+656: The discovery in the Cosmic Lens All-Sky Survey (CLASS) is presented by Myers et al. (1995) with the source redshift by Fassnacht et al. (1996). The imaging modeling is presented by Suyu et al. (2009, 2010). The time-delay measurement is presented by Fassnacht et al. (1999, 2002). The velocity dispersion measurement of 260 km s−1 presented by Suyu et al. (2010) is based on Keck-LRIS spectroscopy. The statistical uncertainty is ±7.7 km s−1 with a systematic spread of ±13 km s−1 depending on wavelength and stellar template solution. The combined uncertainty is 260 ± 15 km s−1. A previous measurement by Koopmans et al. (2003) with 247 ± 35 km s−1 with Echellette Spectrograph and Imager (ESI) on Keck-II is consistent with the more recent one by Suyu et al. (2010). The line-of-sight analysis is presented by Suyu et al. (2010), based on galaxy number counts by Fassnacht et al. (2011).

2. RXJ1131−1231: The discovery is presented by Suyu et al. (2013) and Sluse et al. (2003). The imaging modeling is presented by Suyu et al. (2014) (for HST) and Chen et al. (2019) (for Keck Adaptive Optics data). An independent analysis of the HST data was performed by Birrer et al. (2016). The time-delay measurement is presented by Tewes et al. (2012). The velocity dispersion measurement of 323 ± 20 km s−1 presented by Suyu et al. (2013) is based on Keck-LRIS spectroscopy and includes systematics. The line-of-sight analysis is presented by Suyu et al. (2013).

3. HE0435−1223: The discovery is presented by Wisotzki et al. (2002). The image modeling is presented by Wong et al. (2017) (for HST) and Chen et al. (2019) (for Keck Adaptive Optics data). The time-delay measurement is presented by Bonvin et al. (2016). The velocity dispersion measurement of 222 ± 15 km s−1 presented by Wong et al. (2017) is based on Keck-LRIS spectroscopy and includes systematic uncertainties. An independent measurement of 222 ± 34 km s−1 by Courbin et al. (2011) using VLT is in excellent agreement. The line-of-sight analysis is presented by Rusu et al. (2017).

4. SDSS1206+4332: The discovery is presented by Oguri et al. (2005). The image modeling is presented by Birrer et al. (2019). The time-delay measurement is presented by Eulaers et al. (2013) with an update by Birrer et al. (2019). The velocity dispersion measurement of 290 ± 30 km s−1 presented by Agnello et al. (2016) is based on Keck-DEIMOS spectroscopy and includes systematic uncertainties. The line-of-sight analysis is presented by Birrer et al. (2019).

5. WFI2033−4723: The discovery is presented by Morgan et al. (2004), the image modeling by Rusu et al. (2020) and the time-delay measurement by Bonvin et al. (2019). The velocity dispersion measurement from VLT MUSE is presented by Sluse et al. (2019) with 250 ± 10 km s−1 only accounting for statistical error and 250 ± 19 km s−1 including systematic uncertainties. The line-of-sight analysis is presented by Rusu et al. (2020).

6. DES0408−5354: The discovery is presented by Lin et al. (2017), Diehl et al. (2017). The imaging modeling is presented by Shajib et al. (2020a). A second team within STRIDES and TDCOSMO is performing an independent and blind analysis using a different modeling code (Yildirim et al., in prep.). The time-delay measurement is presented by Courbin et al. (2018). The velocity dispersion measurements are presented by Buckley-Geer et al. (2020). We used the values from Table 3 in Shajib et al. (2020a). The measurements are from Magellan with 230 ± 37 km s−1 (mask A) and 236 ± 42 km s−1 (mask B), from Gemini with 220 ± 21 km s−1 and from VLT MUSE with 227 ± 9 km s−1. The reported values do not include systematic uncertainties and covariances among the different measurements. Following Shajib et al. (2020a) we add a covariant systematic uncertainty of ±17 km s−1 to the reported values. The line-of-sight analysis is presented by Buckley-Geer et al. (2020).

7. PG1115+080: The discovery is presented by Weymann et al. (1980). The image modeling is presented by Chen et al. (2019) using Keck Adaptive Optics. The time-delay measurement is presented by Bonvin et al. (2018), while the line-of-sight analysis by Chen et al. (2019). The velocity dispersion measurement of 281 ± 25 km s−1, presented by Tonry (1998), is based on Keck-LRIS spectroscopy. In this work we add new acquired integral-field spectroscopy obtained with the Multi-Object Survey Explorer (MUSE) on the VLT in March 2019 (0102.A-0600(C), PI Agnello), and we thus go in some detail about the observations. The details and the data will be presented in a forthcoming paper by Agnello et al. (in prep.). At the location of the lens, 3 h of total exposure time were obtained, in clear or photometric conditions and nominal seeing of 0.8″ FWHM. Due to the proximity of the four quasar images to the main galaxy, a dedicated extraction routine was used in order to optimally deblend all components. We followed the same procedure as by Sluse et al. (2019) and Braibant et al. (2014), fitting each spectral channel as a superposition of a Sersic profile (for the main lens) and four point sources as identical Moffat profiles. The separation between the individual components is held fixed to the HST-NICMOS measurements (Sluse et al. 2012).

A nearby star in the MUSE field-of-view was used as a reference PSF. From this direct modeling, the FWHM of the PSF was found to be , with some variation with wavelength that was accounted for in the model-based deblending. This procedure produced an optimal subtraction of the quasar spectra, at least within 1″ from the center of the lens. The lens galaxy 1D spectra were then extracted in two square apertures (, ), and processed with the Penalized PiXel-Fitting (PPXF) code presented in Cappellari & Emsellem (2004) and further upgraded in Cappellari (2017) to obtain velocity dispersions.

The velocity dispersion measurement results from a linear combination of stellar template spectra to which a sum of orthogonal polynomials is added to adjust the continuum shape of the templates to the observed galaxy specttrum. The spectral library used for the fit is the Indo-US spectral library, 1273 stars covering the region from 3460 to 9464 Åat a spectral resolution of 1.35 Å FWHM (Valdes et al. 2004).

We measure for the inner aperture (R <  0.6″) a stellar velocity dispersion value of 277 ± 6.5 km s−1 and for the outer () a value of 241 ± 8.8 km s−1. The uncertainties only include the statistical errors. In order to estimate the systematics, we performed a number of PPXF fits on the smaller aperture, changing each time the wavelength range, the degree of the additive polynomial and the number of stellar templates used to fit the galaxy spectra. We obtained a systematic uncertainty of ±23.6 km s−1 that, as for the case of DES0408, we treat as fully covariant among the two aperture measurements. With the spectral resolution of MUSE, systematic uncertainties are within ≈10% and about three times larger than the nominal, statistical uncertainties thanks to the high signal-to-noise of the spectra.

All the TDCOSMO analyses of lenses used uniform priors on all relevant parameters when performing the inference with a PEMD model13. Six out of the seven lenses were modeled blindly14, that is H0 values were never seen by the modeler at any step of the process.

Detailed line-of-sight analyses for each lens have been performed based on weighted relative number counts of galaxies along the line of sight on deep photometry and spectroscopic campaigns (e.g., Rusu et al. 2017). Furthermore, for a fraction of the lenses, we have used also an external shear constraint inferred by the strong lens modeling to inform the line-of-sight convergence estimate. The weighted galaxy number count and external shear summary statistics have been applied on the Millenium Simulation (Springel et al. 2005) with ray-tracing (Hilbert et al. 2009) to extract a posterior in p(κext) with the prior from the Millenium Simulation and semi-analytic galaxy evolution model with painted synthetic photometry on top (De Lucia & Blaizot 2007)15. The external convergence and shear values from the Millenium simulation are computed from the observer to the source plane, κext ≈ κs. The coupling of the strong lens deflector (e.g., Bar-Kana 1996; McCully et al. 2014; Birrer et al. 2017) is not included in the calculation of κs. Figure 6 shows the κext posteriors for the individual lenses. For the overall sample mean, we get with a scatter of σ(κext) = 0.046 around the mean. Nearby massive galaxies along the line of sight were included explicitly in the modeling where required, and the external convergence term was adapted accordingly in order to not double count mass structure in the analysis. Table 2 presents the redshifts and the relevant lens model posteriors that are used in our analysis.

Table 2.

Overview of the TDCOSMO sample posterior products used in this work.

5.2. TDCOSMO hierarchical inference

We use for each lens the individual time-delay distance likelihood according to Eq. (C.11) that was derived in previous works of this collaboration from a lens model inference on imaging data and the time-delay measurements from the PEMD inference, not including external convergence or internal MST, . We add the same MST transform as a distribution mean λint, 0 and scaling αλ with reff/θE, and with Gaussian scatter across the data set, identical to the TDLMC validation in Sect. 4. The individual p(κext) distributions are added for each lens and in the inference combined with the internal MST parameters.

For the kinematic modeling, we make the same assumptions as for the TDLMC sample (Sect. 4.1) with the anisotropy model of Osipkov (1979), Merritt (1985) (Eq. (51)) with a parameterization of the transition radius relative to the half-light radius (Eq. (52)). The approach is consistent with the previous kinematic analysis and sufficiently verified on the TDLMC to the level of accuracy we can expect from this analysis. We also assume a Hernquist light profile with reff, in conjunction with the power-law lens model posteriors θE and γpl to model the dimensionless kinematic quantity J (Eqs. (16) and (17)), also incorporating the slit mask and seeing conditions of the individual observations.

For each of the lenses in the TDCOSMO sample, we use the distribution p(κext) as derived on the individual blinded analyses and do not invoke an additional population parameter. We leave the hierarchical analysis of the line-of-sight selection to future work. We want to stress that the overall selection bias in this hierarchical approach does not impact the H0 constraints as the kinematics constrains the overall MST (Eq. (34)). An overall shift in the distribution of κs will be compensated by λint in the inference, thus leaving the H0 constraints invariant.

We assume a flat ΛCDM cosmology with a uniform prior on H0 in [0,150] km s−1 Mpc−1. For Ωm we chose the prior based on the Pantheon sample (Scolnic et al. 2018), 𝒩(μ = 0.298, σ = 0.022). We also perform the inference with a flat prior on Ωm in [0.05,0.5] to allow for comparison with the previous work by Wong et al. (2020) and Millon et al. (2020) and to illustrate cosmology dependences in the time-delay cosmography inference. Table 3 summarizes all the hierarchical hyper-parameters sampled in the analysis of this section. The posteriors of the TDCOSMO sample inference are presented in Fig. 7.

Table 3.

Summary of the model parameters sampled in the hierarchical inference on the TDCOSMO sample in Sect. 5 and posteriors presented in Fig. 7.

thumbnail Fig. 7.

Hierarchical analysis of the TDCOSMO-only sample when constraining the MST with kinematic information. Parameter and priors are specified in Table 3. Orange contours correspond to the inference with uniform prior on Ωm, 𝒰([0.05, 0.5]), while the purple contours correspond to the prior based on the Pantheon sample with 𝒩(μ = 0.298, σ = 0.022) (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/TDCOSMO_sample/tdcosmo_sample.ipynb).

Open with DEXTER

For the tight prior on Ωm, we measure km s−1 Mpc−1. For an unconstrained relative expansion history with a prior on Ωm uniform in [0.05,0.5], we measure km s−1 Mpc−1. The 9% precision on H0 is significantly inflated relative to previous studies with the same data set (Wong et al. 2020; Millon et al. 2020). The increase in uncertainty with respect to the H0LiCOW analysis is attributed to two main factors: (1) we relaxed the assumption of NFW+stars or power-law mass density profiles; (2) we considered the impact of covariance between lenses when accounting for uncertainties potentially arising from assumptions about mass profile and stellar anisotropy models. As we show in the next sections, however, this uncertainty can be reduced by adding external information to further constrain the mass profile and anisotropy of the deflectors. The inferred scatter in λint, σ(λint), is consistent with zero. This is a statement on the internally consistent error bars on H0 among the TDCOSMO sample (Wong et al. 2020; Millon et al. 2020).

6. SLACS analysis of galaxy density profiles

Gravitational lenses with imaging and kinematics data can add valuable information about the mass profiles of the lenses. Even though the kinematics data in the current TDCOSMO sample is limited, an additional sufficiently large data set with precise measurements can significantly improve the precision on the mass profiles of the population and thus on the Hubble constant. Resolved kinematics observations may in addition provide constraints on the anisotropy distribution of stellar orbits.

When incorporating external data sets as part of the hierarchical framework, it is important that those external lenses are drawn from the same population as the time-delay lenses – unless explicitly marginalized over population differences. Provided that (i) the lensing sample has a known selection function, (ii) the lens modeling is performed to the same level of precision and with the same model assumptions as the time-delay lenses, (iii) the kinematic modeling assumptions are identical and (iv) the anisotropy uncertainties are mitigated on the population level, we can fold in the extracted likelihood (Eq. (C.12)) into the hierarchical analysis, applying the same population dependence on λint and aani.

Selection biases can arise from different aspects. Ellipticity and shear naturally increase the abundance of quadruple lenses relative to double lenses. Holder & Schechter (2003) use N-body simulations to estimate the level of external shear due to structure near the lens and conclude that the local environment is the dominant contribution that drives the external shear bias in quadruple lenses. Huterer et al. (2005) investigate the external shear bias and conclude that this effect is not sufficient to explain the observed quadruple-to-double ratio. Collett & Cunnington (2016) conclude, based on idealized simulations, that selection based on image brightness and separation leads to significant selection bias in the slope of the mass profiles. In addition, Collett & Cunnington (2016) also find a line-of-sight selection bias in quadruply lensed quasars relative to the overall population on the level of 0.9%. The bias is less prominent for doubly imaged quasars. The specific discovery channel can also lead to selection effects. Dobler et al. (2008) note that a spectroscopically selected search, as performed for the Sloan Lens ACS (SLACS) survey (Bolton et al. 2006), can lead to significant biases on the selected velocity dispersion in the resulting sample. However, Treu et al. (2006) show that, at fixed velocity dispersion, the SLACS sample is indistinguishable from other elliptical galaxies.

In this section we present a hierarchical analysis of the SLACS sample (Bolton et al. 2006, 2008) following the same hierarchical approach as the TDCOSMO sample, based on the imaging modeling by Shajib et al. (2020b). The SLACS sample of strong gravitational lenses is a sample of massive elliptical galaxies selected from the Sloan Digital Sky Survey (SDSS) by the presence in their spectra of emission lines consistent with a higher redshift. Follow-up high-resolution observations with HST revealed the presence of strongly lensed sources. The SLACS data set allows us to further constrain the population distribution in the mass profile parameter λint and the anisotropy distribution aani and, thus, can add significant information to the TDCOSMO sample to be used jointly in Sect. 7 to constrain H0.

In Sect. 6.1 we describe the imaging data and lens model inference. In Sect. 6.2 we describe the spectroscopic data set used and how we model it, including VLT VIMOS IFU data for a subset of the lenses. We analyze the selection effect of the SLACS sample in Sects. 6.3 and 6.4 we constrain the line-of-sight convergence for the individual lenses. In Sect. 6.5 we present the results of the hierarchical analysis of the SLACS sample in regard to mass profile and anisotropy constraints.

6.1. SLACS imaging

To include additional lenses in the hierarchical analysis, we must ensure that the quality and the choices made in the analysis are on equal footing with the TDCOSMO sample. Shajib et al. (2020b) presents a homogeneous lens model analysis of 23 SLACS lenses from HST imaging data. The lens model assumptions are a PEMD model with external shear, identical to the derived products we are using from the TDCOSMO sample. The scaling of the analysis was made possible by advances in the automation of the modeling procedure (e.g., Shajib et al. 2019) with the DOLPHIN pipeline package. The underlying modeling software is LENSTRONOMY (Birrer & Amara 2018; Birrer et al. 2015) for which we also performed the TDLMC validation (Sect. 4).

Shajib et al. (2020b) first select 50 SLACS lenses for uniform modeling from the sample of 85 lenses presented by Auger et al. (2009). The selection criteria for these lenses are: (i) no nearby satellite or large perturber galaxy within approximately twice the Einstein radius, (ii) absence of multiple source galaxies or complex structures in the lensed arcs that require large computational cost for source reconstruction, and (iii) the main deflector galaxy is not disk-like. These criteria are chosen so that the modeling procedure can be carried out automatically and uniformly without tuning the model settings on a lens-by-lens basis. Using the DOLPHIN package on top of LENSTRONOMY, a uniform and automated modeling procedure is performed on the 50 selected lenses with V-band data (Advance Camera for Surveys F555W filter, or Wide Field and Planetary Camera 2 F606W filter).

After the modeling, 23 lenses are selected to have good quality models. The criteria for this final selection are: (i) good fitting to data by visually inspecting the residual between the image and the model-based reconstruction, and (ii) the median of the power-law slope does not diverge to unusual values (i.e., ≲1.5 or ≳2.5)16. For the TDCOSMO sample, iterative PSF corrections have been performed, based on the presence of the bright quasar images, to guarantee a well matched and reliable PSF in the modeling. For the SLACS lenses, such an iterative correction on the image itself cannot be performed due to the absence of quasars in these systems. Nevertheless, extensive tests with variations of the PSF have been performed by Shajib et al. (2020b) and the impact on the resulting power-law slope inference was below ∼0.005 on the population mean of γpl. The half light radius for the deflector galaxies are taken from Auger et al. (2009) in V-band (measured along the intermediate axis).

6.2. SLACS spectroscopy

The constraints on the MST rely on the kinematics observations. In this section we provide details on the data set and reduced products we are using in this work, on top of the already described ones for the TDCOSMO lenses. These include SDSS’s Baryon Oscillation Spectroscopic Survey (BOSS) fiber spectroscopy (Dawson et al. 2013) and VLT VIMOS IFU observations.

6.2.1. SDSS fiber spectroscopy

All the SLACS lenses have BOSS spectra available as part of SDSS-III. The fiber diameter is 3″ and the nominal seeing of the observations are 14 FWHM. The measurements of the velocity dispersion from the SDSS reduction pipeline were originally presented by Bolton et al. (2008). However, in this work, we use improved measurements of the velocity dispersion, determined using an improved set of templates as described in Shu et al. (2015). The SDSS measurements are in excellent agreement with the subsample measured with VLT X-shooter presented by Spiniello et al. (2015).

6.2.2. VLT VIMOS IFU data

The VLT VIMOS IFU data set is described in Czoske et al. (2008) and subsequently used in Barnabè et al. (2009, 2011), Czoske et al. (2012). The VIMOS fibers were in a configuration with spatial sampling of 0.67″, and the seeing was 08 FWHM.

The first moment (velocity) and second moment (velocity dispersion) of the individual VIMOS fibers are fit with a single stellar template for each fiber individually and the uncertainties in the measurements are quantified within Bayesian statistics. Templates were chosen by fitting a random sample of IndoUS spectra to the aperture-integrated VIMOS IFU spectra and selecting one of the best-fitting (in the least-squares sense) template candidates (we refer to details to Czoske et al. 2008). Marginalization over template mismatch adds another 5−10% measurement uncertainties. Within this additional error budget, the integrated velocity dispersion measurements of Czoske et al. (2008) are consistent with the SDSS measured values of Bolton et al. (2008). We bin the fibers in radial bins in steps of 1″ from the center of the deflector. The binning is performed using luminosity weighting and propagation of the independent errors to the uncertainty estimate per bin. Where necessary, we exclude fibers that point on satellite galaxies or line-of-sight contaminants. In this work, we make use of the relative velocity dispersion measurements in radial bins when inferring H0. We do so by introducing a separate internal MST distribution λifu, effectively replacing λint when evaluating the likelihood of the IFU data. λifu is entirely constrained by the IFU data. The MST information that propagates in the joint constraints of TDCOSMO+SLACS analysis, λint, (Sect. 7) is derived from the SDSS velocity dispersion measurements only. In this form, the IFU data informs the anisotropy parameter but not the mass profile directly. We leave the amplitude calibration and usage of this data set to constrain the MST for future work.

From the original sample of 17 SLACS lenses with VIMOS observations, we drop five objects that are fast rotators (when the first moments dominate the averaged dispersion in the outer radius bin) and one slow rotator with velocity dispersion >380 km s−1. This is necessary to match this sample with the TDCOSMO one in velocity dispersion space; the fast rotators are, in fact, all in a lower velocity dispersion range (σP in [185,233] km s−1). Finally, we excluded one more galaxy for which there is no estimate of the Einstein radius, and thus we cannot combine lensing and dynamics. In this way, we end up with a sample of ten lenses, prior to further local environment selection.

6.3. SLACS selection function

The SLACS lenses were preselected from the spectroscopic database of the SDSS based on the presence of absorption-dominated galaxy continuum at one redshift and nebular emission lines (Balmer series, [OII] 3727 A, or [OIII] 5007 A) at another, higher redshift. Details on the method and selection can be found in Bolton et al. (2004, 2006) and Dobler et al. (2008). The lens and source redshifts of the SLACS sample are significantly lower than for the TDCOSMO sample.

Treu et al. (2009) studied the relation between the internal structure of early-type galaxies and their environment with two statistics: the projected number density of galaxies inside the tenth nearest neighbor (Σ10) and within a cone of radius one h−1 Mpc (D1) based on photometric redshifts. It was observed that the local physical environment of the SLACS lenses is enhanced compared to random volumes, as expected for massive early-type galaxies, with 12 out of 70 lenses in their sample known to be in group/cluster environments.

In this study, we are specifically only looking for lenses whose lensing effect can be described as the mass profile of the massive elliptical galaxy and an uncorrelated line-of-sight contribution. Assuming SLACS and TDCOSMO lenses are galaxies within the same homogeneous galaxy population and with the local environment selection of SLACS lenses, the remaining physical mass components in the deflector model are the same physical components of the lensing effect we model in the TDCOSMO sample. The uncorrelated line-of-sight contribution can be characterized based on large scale structure simulations.

6.3.1. Deflector morphology and lensing information selection

Our first selection cut on the SLACS sample is based on Shajib et al. (2020b), which excludes a subset of lenses based on their unusual lens morphology (prominent disks, two main deflectors, or complex source morphology) to derive reliable lensing properties using an automated and uniform modeling procedure.

This first cut reduced the total SLACS sample of 85 lenses, presented by Auger et al. (2009), to 51 lenses17. Out of these 51 lenses, 23 lenses had good quality models from an automated and uniform modeling procedure as described in Sect. 6.1. Producing good quality models for the rest of the SLACS lenses would require careful treatment on a lens-by-lens basis, which was out of the scope of Shajib et al. (2020b).

6.3.2. Mass proxy selection

We want to make sure that the deflector properties are as close as possible to the TDCOSMO sample. To do so without introducing biases regarding uncertainties in the velocity dispersion measurements, we chose a cut based on Singular Isothermal Sphere (SIS) equivalent dispersions, σSIS, derived from the Einstein radius and the lensing efficiency only. The deflectors of the TDCOSMO sample span a range of σSIS in [200,350] km s−1 and we select the same range for the SLACS sample.

6.3.3. Local environment selection

We use the DESI Legacy Imaging Surveys (DLS; Dey et al. 2019) to characterize the environment of the SLACS lenses. We query the DR7 Tractor source photometry catalog (Lang et al. 2016) removing any object that is morphologically consistent with being a point source convolved with the DLS point spread function. We use the R band data to count objects with 18 <  R <  23 within 120″ of the lens galaxy but more than 3″ from the lens.

We quantify the environment with two numbers: N2′, the total number of galaxies within 2 arcmin and an inverse projected distance weighted count N1/r within the same 2 arcmin aperture, defined as (Greene et al. 2013)

(54)

N2′ and N1/r are physically meaningful numbers for our analysis as N2′ should approximately trace the total mass close enough to significantly perturb the lensing (see Collett et al. 2013), and N1/r should be skewed larger by masses close along the line of sight of the lens which are likely to have the most significant perturbative effect. We assess the uncertainty on N2′ and N1/r by taking every object within 120″ of the lens and bootstrap resampling from their R band magnitude errors, before reapplying the 18 <  R <  23 cut. Where the SLACS lens is not in the DLS DR7 footprint we queried the DLS DR8 catalog instead. To put N2′ and N1/r into context, we perform the same cuts centered on 105 random points within the DLS DR7 footprint. Dividing the SLACS N2′ and N1/r by the median ⟨N2′⟩rand and ⟨N1/rrand of the calibration lines of sight allows us to assess the relative over-density of the SLACS lenses as

(55)

and

(56)

We compare this metric on our sample with the overlapping sample of Treu et al. (2009) where local 3-dimensional quantities in the form of D1 are available, and we find good agreement between these two statistics in terms of a rank correlation.

We remove lenses that have ζ1/r >  2.10 within the 2 arcminutes aperture from our sample. This selection cut corresponds to D1 = 1.4 Mpc−3 for the subset by Treu et al. (2009). Independently of the ζ1/r cut, we check and flag all lenses within the Shajib et al. (2020b) sample that have prominent nearby perturbers present in the HST data within 5″. We do not find any additional lenses with prominent nearby perturbers not already removed by the selection cut of ζ1/r >  2.10.

6.3.4. Combined sample selection

With the combined selection on the SLACS sample based on the morphology, mass proxy, local environment, and for the IFU lenses also rotation, we end up with 33 SLACS lenses of which nine lenses have IFU data. 14 lenses out of the sample have quality lens models by Shajib et al. (2020b), including five lenses with IFU data. Figure 8 shows how the individual lenses among the different samples, TDCOSMO, SLACS and the subset with IFU data are distributed in key parameters of the deflector.

thumbnail Fig. 8.

Sample selection of the SLACS lenses being added to the analysis and comparison with the TDCOSMO data set. The comparisons are in lens redshift, zlens, source redshift, zsource, measured velocity dispersion, σP, half light radius of the deflector, reff, Einstein radius of the deflector, θE, the ratio of half light radius to Einstein radius, reff/θE, and the SIS equivalent velocity dispersion estimated from the Einstein radius and a fiducial cosmology, σSIS. Open dots correspond to lenses included in our selection without quality lens models. Red points correspond to SLACS lenses which have VIMOS IFU data (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/sample_selection.ipynb).

Open with DEXTER

We discuss possible differences between the SLACS and TDCOSMO samples and the possibility of trends within the samples impacting our analysis in a systematic way in Sect. 8.3.2 after presenting the results of the hierarchical analysis of the joint sample. We list all the relevant measured values and uncertainties of the 33 SLACS lenses in Appendix E.

6.4. Line of sight convergence estimate

We compute the probability for the external convergence given the relative number counts, P(κext|ζ1, ζ1/r), following Greene et al. (2013) (see e.g., Rusu et al. 2017, 2020; Chen et al. 2019; Buckley-Geer et al. 2020). In brief, we select from the Millennium Simulation (MS; Springel et al. 2005) line of sights which satisfy the relative weighted number density constraints measured above, in terms of both number counts and 1/r weighting (Eq. (54)). While the MS consists only of dark matter halos, we use the catalog of galaxies painted on top of these halos following the semi-analytical models of De Lucia & Blaizot (2007). We implement the same magnitude cut, aperture radius etc. which were employed in measuring the relative weighted number densities for the SLACS lenses, in order to compute ζ1, ζ1/r corresponding to each line of sight in the MS. We then use the κ maps computed by Hilbert et al. (2009) and read off the values corresponding to the location of the selected line of sight, thus constructing the p(κext|ζ1, ζ1/r) probability density function (PDF). The Hilbert et al. (2009) maps were computed for a range of source redshift planes. Over the range spanned by the source redshifts of the SLACS lenses, there are 17 MS redshift planes, with spacing Δz ∼ 0.035−0.095. We used the maps best matching the source redshift of each SLACS lens. For 23 of the SLACS lenses there are available external shear measurements by Shajib et al. (2020b), which we used, optionally, as a third constraint. Compared to previous inferences of p(κ) for the TDCOSMO lenses, we made two computational simplifications to our analysis, in order to be able to scale our technique to the significantly larger number of lenses: (1) We did not resample from the photometry of the MS galaxies, taking into account photometric uncertainties similar to those in the observational data. A toy simulation showed that this step results in negligible differences. (2) We use only 1/8 of the lines of sight in the MS. We then checked that this results in Δκ ≲ 0.001 offsets, negligible for the purpose of our analysis.

Figure 9 shows the p(κext|ζ1, ζ1/r) distributions for the sub-selected sample based on morphology and local environment. As expected from the significantly lower source redshifts of the SLACS sample compared to the TDCOSMO lenses, most of p(κ) PDFs for the individual lenses are very narrow and peak at ∼zero, with dispersion ∼0.01. This is because the volume is smaller and thus there are relatively few structures in the MS at these low redshifts to contribute. In fact, the relative weighted number density constraints have relatively little impact on most of the p(κ) distributions, which resemble the PDFs for all lines of sight. Finally, we note that, while our approach to infer p(κ) for the SLACS lenses is homogeneous, this is not the case for the TDCOSMO lenses. This is by necessity, as the environmental data we used for the TDCOSMO lenses has varied in terms of depth, number of filters and available targeted spectroscopy. Nonetheless, as we have shown through simulations by Rusu et al. (2017, 2020), such differences do not bias the p(κ) inference.

thumbnail Fig. 9.

External convergence posteriors of the 33 SLACS lenses that pass our morphology and local environment selection cut based on the weighted number counts (gray dashed lines) and the population distribution (black solid line) (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/sample_selection.ipynb).

Open with DEXTER

6.5. SLACS inference

Here we present the hierarchical inference on the mass profile and anisotropy parameters from the selected sample of the SLACS lenses. We remind the reader that we use 33 SLACS lenses, of which 14 have imaging modeling constraints on the power-law slope γpl. Nine of the lenses in our final sample have also VLT VIMOS IFU constraints in addition to SDSS spectroscopy (five of which have imaging modeling constraints on the power-law slope). The separate inference presented in this section is meant to provide consistency checks and to gain insights into how the likelihood of the SLACS data set is going to impact the constraints on the mass profiles, and thus H0, when combining with the TDCOSMO data set.

We are making use of the marginalized posteriors in the lens model parameters of Shajib et al. (2020b) in the same way as for the TDLMC and TDCOSMO sample. For SLACS lenses that do not have a model and parameter inference by Shajib et al. (2020b), we use the Einstein radii measured by Auger et al. (2009) derived from a singular isothermal ellipsoid (SIE) lens model. For the power-law slopes of those lenses we apply the inferred Gaussian population distribution prior on γpl from the selected sample which has measured values, with γpl, pop = 2.10 ± 0.16. Figure 10 presents the imaging data inferred γpl for the 14 quality lenses selected in our sample by Shajib et al. (2020b).

thumbnail Fig. 10.

Power-law slope γpl inferences obtained from HST imaging data for the 14 SLACS lenses within our selection cut from Shajib et al. (2020a,b). We derive a population distribution from these lenses to be applied on the subset of lenses without measured γpl from imaging data. The black dashed line indicates the population mean and the blue band the 1-sigma population width based on the 14 individual measurements (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/sample_selection.ipynb).

Open with DEXTER

Table 4 presents the parameters and priors used in the hierarchical inference of this section. In particular, we fix the cosmology to assess constraining power and consistency with the TDCOSMO data set. We separate the inference on λint, 0 of the VIMOS IFU data set from the SDSS measurements to assess systematic differences between the two data products. Further more, we use a uniform prior in aani, 𝒰(aani), rather than a logarithmic prior 𝒰(log(aani)), to assess and illustrate the information on the anisotropy parameter from the IFU data set.

Table 4.

Summary of the model parameters sampled in the hierarchical inference on the SLACS sample of Sect. 6.

For the analysis of the SLACS-only sample in this section, we fix the cosmological model. The cosmological dependence folds in the prediction of the velocity dispersion through the distance ratio Ds/Dds (Eq. (17)). This ratio is not sensitive to H0 and the SLACS-only data set is not constraining H0. When combining the SLACS and TDCOSMO sample in the next section, the cosmology dependence is fully taken into account.

We perform two posterior inferences: one with the SDSS velocity dispersion data only, and one combining SDSS and VIMOS IFU binned dispersions. Figure 11 shows the two different posteriors. The constraints on λint (parameters λint, 0, αλ, σ(λint)) come for all three cases entirely from the kinematics of the SDSS measurements.

thumbnail Fig. 11.

Posterior distribution for the SLACS sample with priors according to Table 4. Orange: inference with the SDSS spectra. Purple: inference with SDSS spectra and VIMOS IFU data set. The posterior of λint, 0 was blinded during the analysis (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/SLACS_sample/SLACS_constraints.ipynb).

Open with DEXTER

All the parameters are statistically consistent with each other and the TDCOSMO analysis of Sect. 5 except the posterior in the scatter of the internal MST, σ(λint). The TDCOSMO constraints of σ(λint) are consistent with zero scatter in the mass profile parameter and 2-sigma bound at 0.1, while the inference of the SLACS sample results in a larger scatter. An underestimation of uncertainties in the velocity dispersion measurements, if not accounted for in the analysis, will directly translate to an increase in σ(λint). We point out the excellent agreement of the anisotropy distribution with the TDLMC Rung3 hydrodynamical simulations (Sect. 4).

7. Hierarchical analysis of TDCOSMO+SLACS

We describe now the final and most stringent analysis of this work, obtained by combining the analysis of the TDCOSMO lenses, presented in Sect. 5, and that of the SLACS sample, presented in Sect. 6. The parameterization and priors have been validated on the TDLMC mock data set in Sect. 4. We remind the reader that the choices of the analyses are identical and thus we can combine the TDCOSMO and SLACS sample on the likelihood level. We define the parameterization and priors of our hierarchical model in Sect. 7.1 and present the result and the H0 measurement in Sect. 7.2.

7.1. Parameterization and priors

For our final H0 measurement, we assume a flat ΛCDM cosmology with uniform prior in H0 in [0, 150] km s−1 Mpc−1 and a narrow prior on Ωm with 𝒩(μ = 0.298, σ = 0.022) from the Pantheon SNIa sample (Scolnic et al. 2018, see Sect. 3.2.4). For λint, we assume an identical distribution for the selected population of the SLACS lenses and the TDCOSMO sample for the scaling in reff/θE (Eq. (50)). We also assume the same stellar anisotropy population distributions for the SLACS and TDCOSMO lenses. To account for potential systematics in the VIMOS IFU measurement (see Sect. 6.2.2), we introduce a separate a separate internal MST distribution λifu, effectively replacing λint when fitting the IFU data. This approach allows us to use the anisotropy constraints from the IFU data while not requiring a perfect absolute calibration of the measurements. For the external convergence we use the individual p(κext) distributions from the two samples.

As discussed in Sect. 6, there is an inconsistency in the inferred spread in the λint distribution between the SLACS and TDCOSMO sample. We attribute this inconsistency to uncertainties that were not accounted for in the velocity dispersion measurements of the SDSS data products. In our joint analysis, we add a parameter that describes an additional relative uncertainty in the velocity dispersion measurements, σσP, sys, such that the total uncertainty in the velocity dispersion measurements is the square of the quoted measurement uncertainty plus this unaccounted term,

(57)

σσP, sys is the same for all the SDSS measured velocity dispersions. Table 5 presents all the parameters being fit for, including their priors, in our joint analysis of the SLACS and TDCOSMO sample18.

Table 5.

Summary of the model parameters sampled in the hierarchical inference on the TDCOSMO+SLACS sample.

7.2. Results

Here we present the posteriors of the joint hierarchical analysis of 33 SLACS lenses (nine of which have IFU data) and the seven quasar time-delay TDCOSMO lenses for the parameterization and priors described in Table 5. To trace back information to specific data sets, we sample different combinations of the TDCOSMO and SLACS data sets under the same priors. The TDCOSMO-only inference was already presented in Sect. 5 and results in km s−1 Mpc−1. Besides the TDCOSMO-only result, we perform the inference for the TDCOSMO+SLACSIFU data set, effectively allowing anisotropy constraints being used on top of the TDCOSMO data set, resulting in km s−1 Mpc−1; the TDCOSMO+SLACSSDSS data set, using the SLACS lenses with their SDSS spectroscopy to inform the analysis, results in km s−1 Mpc−1. For our final inference of this work of the joint data sets of TDCOSMO+SLACSSDSS + IFU, we measure km s−1 Mpc−1.

Figure 12 presents the key parameter posteriors of the TDCOSMO-only, TDCOSMO+SLACSIFU, TDCOSMO+SLACSSDSS, and the TDCOSMO+SLACSSDSS + IFU analyses. Not shown on the plot are the Ωm posteriors (effectively identical to the prior), the σσP, sys posteriors for the SDSS kinematics measurements, the distribution scatter parameters σ(λint and σ(aani), and the IFU calibration nuisance parameter λifu. All the one-dimensional marginalized posteriors, except for the nuisance parameter λifu, of the different combinations of the data sets are provided in Table 6.

thumbnail Fig. 12.

Posterior distributions of the key parameters for the hierarchical inference. Blue: constraints from the TDCOSMO-only sample. Violet: constraints with the addition of IFU data of nine SLACS lenses to inform the anisotropy prior on the TDCOSMO sample, TDCOSMO+SLACSIFU. Orange: constraints with a sample of 33 additional lenses with imaging and kinematics data (HST imaging + SDSS spectra) from the SLACS sample, TDCOSMO+SLACSSDSS. Purple: joint analysis of TDCOSMO and 33 SLACS lenses with SDSS spectra of which nine have VIMOS IFU data, TDCOSMO+SLACSSDSS + IFU. Priors are according to Table 5. The 68th percentiles of the 1D marginalized posteriors are presented in Table 6. The posteriors in H0 and λint, 0 were held blinded during the analysis (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/joint_inference.ipynb).

Open with DEXTER

Table 6.

Marginalized posteriors of our hierarchical Bayesian cosmography inference based on the priors and parameterization specified in Table 5 for a flat ΛCDM cosmology.

We compare the best fit model prediction of the joint TDCOSMO+SLACSSDSS + IFU inference to the time-delay distance and kinematics of the TDCOSMO data set in Fig. 13, to the SDSS velocity dispersion measurements in Fig. 14 and to the IFU data set in Fig. 15. The model prediction uncertainties include the population distributions in λint and aani and the measurement uncertainty in the SDSS and VIMOS velocity dispersion uncertainties include the inferred σσP, sys uncertainty.

thumbnail Fig. 13.

Illustration of the goodness of the fit of the maximum likelihood model of the joint analysis in describing the TDCOSMO data set. Blue points are the measurements with the diagonal elements of the measurement covariance matrix. Orange points are the model predictions with the diagonal elements of the model covariance uncertainties. Left: comparison of measured time-delay distance from imaging data and time delays compared with the predicted value from the cosmological model, the internal and external MST (and their distributions). Right: comparison of the velocity dispersion measurements and the predicted values. In addition to the MST terms, the uncertainty in the model also includes the uncertainty in the anisotropy distribution aani. For lenses with multiple velocity dispersion measurements, the diagonal terms in the error covariance are illustrated (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/joint_inference.ipynb). (a) Fit to the time-delay distance. (b) Fit of velocity dispersion.

Open with DEXTER

thumbnail Fig. 14.

Illustration of the goodness of the fit of the maximum likelihood model of the joint analysis in describing the SDSS velocity dispersion measurements of the 34 SLACS lenses in our sample. Blue points are the measurements with the diagonal elements of the measurement covariance matrix. Orange points are the model predictions with the diagonal elements of the model covariance uncertainties. The measurement uncertainties include the uncertainties in the quoted measurements and the additional uncertainty of σσP, sys. The model uncertainties include the lens model uncertainties and the marginalization over the λint and aani distribution (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/joint_inference.ipynb).

Open with DEXTER

thumbnail Fig. 15.

Illustration of the goodness of the fit of the maximum likelihood model of the joint analysis in describing the VIMOS radially binned IFU velocity dispersion measurements of the nine SLACS lenses with VIMOS data in our sample. Blue points are the measurements with the diagonal elements of the measurement covariance matrix. Orange points are the model predictions with the diagonal elements of the model covariance uncertainties. The measurement uncertainties include the uncertainties in the quoted measurements and the additional uncertainty of σσP, sys. The model uncertainties include the lens model uncertainties and the marginalization over the λint and aani distribution (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/joint_inference.ipynb).

Open with DEXTER

In Fig. 16 we assess trends in the fit of the kinematic data in regards to lensing deflector properties. We see that with the reff/θE scaling by αλ (Eq. (50)) we can remove systematic trends in model predictions. We do not find statistically significant remaining trends in our data set beyond the ones explicitly parameterized and marginalized over.

thumbnail Fig. 16.

Relative difference of the measured vs the predicted velocity dispersion for the SLACS and TDCOSMO sample as a function of different parameters associated with the line-of-sight and the lensing galaxy. In particular, these are the relative inverse distance weighted over density ζ1/r, the ratio of half-light radius to Einstein radius reff/θE, lens redshift zlens, and SIS equivalent velocity dispersion σSIS (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/joint_inference.ipynb).

Open with DEXTER

8. Discussion

In this section19, we discuss the interpretation of our measurement of H0, the robustness of the uncertainties, and present an avenue for further improvements in the precision while maintaining accuracy. We first summarize briefly the key assumptions of this work, and give a physical interpretation of the results (Sect. 8.1). Second, we estimate the contribution of each individual assumption and dataset to the total error budget of the current analysis on H0 (Sect. 8.2). Third, we discuss specific aspects of the analysis that need further investigations to maintain accuracy with increased precision in Sect. 8.3. Fourth, in Sect. 8.4 we present the near future prospects for collecting data sets and revising the analysis to increase further the precision on H0 with strong lensing time-delay cosmography. Finally, in Sect. 8.5, we compare and discuss the H0 measurement of this work with previous work by the TDCOSMO collaboration.

8.1. Physical interpretation of the result

While consistent with the results of Wong et al. (2020), Millon et al. (2020), our inference of H0 has significantly lower precision for the TDCOSMO sample, even with the addition of external datasets from SLACS. The larger uncertainty was expected and is a direct result of relaxing the assumptions on the mass profile. By introducing a mass-sheet degeneracy parameter, we add the maximal degree of freedom in H0 while having minimal constraining power by lensing data on their own. This is the most conservative approach when adding a single degree of freedom in our analysis. While mathematically this result is clearly understood, it is worth discussing the physical interpretation of this choice.

If we had perfect cosmological numerical simulations or perfect knowledge of the internal mass distribution within elliptical galaxies, we would not have to worry about the internal MST. The approach chosen by our collaboration (Wong et al. 2020; Shajib et al. 2019; Millon et al. 2020) was to assume physically motivated mass profiles with degrees of freedom in their parameters. In particular, the collaboration used two different mass profiles, a power-law elliptical mass profile, and a composite mass profile separating the luminous component (with fixed mass-to-light ratio) and a dark component described as a NFW profile. The good fit to the data, the small pixellated corrections on the profiles from the first lens system (Suyu et al. 2010), and the good agreement of H0 inferred with the two mass profiles was a positive sanity check on the result (Millon et al. 2020).

In this paper we have taken a different viewpoint, and asked how much can the mass profiles depart from a power-law and still be consistent with the data. By phrasing the question in terms of the MST we can conveniently carry out the calculations, because the MST leaves the lensing observables unchanged and therefore it corresponds to minimal constraints and assumptions, and thus maximal uncertainties with one additional degree of freedom. However, after the inference, one has to examine the inferred MST transformed profile and evaluate it in comparison with existing and future data to make sure it is realistic. We know that the exact MST cannot be the actual answer because profiles have to go zero density at large radii, but the approximate MST discussed in Sect. 2 provides a convenient interpretation with the addition of a cored mass component.

Figure 17 illustrates a cored mass component approximating the MST inferred from this work, λint = 0.91 ± 0.04, in combination with a power-law model inferred from the population mean of the SLACS analysis by Shajib et al. (2020b). The analysis presented here guarantees that the inferred mass profile is consistent with the properties of TDCOSMO and SLACS lenses. We discuss below how additional data may allow us to constrain the models even further and thus reduce the overall uncertainty while keeping the assumptions at a minimum.

thumbnail Fig. 17.

Illustration of the inferred mass profile of the joint TDCOSMO+SLACSSDSS + IFU analysis. A pure power-law with γpl = 2.10 ± 0.05 is shown in orange. In blue is the result of this work of λint = 0.91 ± 0.045 when interpreted as a cored mass component with Rc uniform in [3″, 10″]. Three dimensional density are illustrated on the left and the lensing convergence on the right. The dashed vertical line on the right panels indicates the Einstein radius. Relative difference in respect to the power-law model are presented in the bottom panels (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER

8.2. Statistical error budget and known systematics

The total error budget of 5% on H0 in our combined TDCOSMO+SLACS analysis can be traced back to specific aspects of the data and the uncertainties in the model components/assumptions. Fixing λint to a single-valued number (i.e., λint = 1) is equivalent to assuming a power-law profile and leads to an uncertainty in H0 of 2% (Millon et al. 2020). By subtracting in quadrature 2% from our total uncertainty, we estimate that the total error contribution of the MST (λint) to the error budget is 4.5%. Once the MST is introduced, the uncertainty in the mass profile is dominated by uncertainties in the measurement and modeling assumptions of the velocity dispersion. The statistical constraints on the combined velocity dispersion measurements of 33 SLACS lenses with SDSS spectroscopy, accounting for the σσP, sys contribution, and the TDCOSMO spectroscopic data set contribute 3% to the total error budget. The remaining 3.5% error contribution (in quadrature) to the total H0 error budget arises in equal parts from the uncertainty in the anisotropy prior distribution (⟨aani⟩, σ(aani)) and the MST dependence with reff/θE (αλ). The uncertainty in the line-of-sight selection effect of the SLACS sample contributes a statistical uncertainty smaller than 0.5%. We note that an overall unaccounted-for shared κext term of the ensemble of lenses in our sample would be mitigated through our MST parameterization and thus not affect our H0 inference.

8.3. Unaccounted-for systematics

Our framework is conservative in the sense that it imposes minimal assumptions of the mass profile in regards to H0. Furthermore, the methods presented here have been internally reviewed and validated on the hydrodynamical simulations used in the TDLMC (Ding et al. 2018, 2020) (Sect. 4). Despite the known limitations of current numerical simulations at the sub-kpc scale, the blind validation on external data corroborates our methodology. In this section, we discuss aspects of our analysis that are not part of our validation scheme. In particular, we discuss uncertainties and potential systematics in the kinematics measurements and selection effects of the different lens samples used in this work. At the current level of precision, these are all subdominant effects, but they may be relevant as we further increase the precision.

8.3.1. Uncertainties in the kinematics measurement and modeling

Under the assumptions of this analysis, aperture stellar kinematic measurements drive the overall precision by providing the information needed to mitigate the MST. Given its crucial role, we highlight here the limitations of our kinematic treatment, in order to point the way to further improvements. First, we used a heterogeneous set of stellar velocity dispersions. The TDCOSMO measurements are based on large telescope high-quality data and were the subject of extensive tests to assess systematic measurements, sometimes through repeated measurements. The nominal uncertainties are thus accurate, resulting in the internal consistency of all the TDCOSMO systems with a scatter on λint consistent with zero20.

The SLACS-only analysis with the reported uncertainties of the stellar velocity dispersions leads to an inferred scatter in λint of about 10%. Assuming the same scatter in λint among the TDCOSMO and SLACS lenses, the discrepancy in the inferred σ(λint) between the two samples indicates that the reported uncertainties of the stellar velocity dispersions of the SLACS lenses do not reflect the total uncertainty. For the present analysis, we have addressed this issue by adding additional terms of uncorrelated errors. However, future work should aim to improve the determination of systematics going back to the original data (or acquiring better data), and contemplate the possibility of correlated calibration errors, as due for example to the choice of stellar library or instrumental setup. Second, our analysis is based on spherical Jeans models, assuming anisotropy of the Osipkov–Merritt form. These approximations are sufficient given the current uncertainties and constraints, but future work should consider at least axis-symmetric Jeans modeling (e.g., Cappellari 2008; Barnabè et al. 2012; Posacki et al. 2015; Yıldırım et al. 2020), and consider alternate parameterizations of anisotropy. Another possibility is the use of axisymmetric modeling of the phase-space distribution function with a two-integral Schwarzschild method by Cretton et al. (1999), Verolme & de Zeeuw (2002) as performed by Barnabè & Koopmans (2007), Barnabè et al. (2009).

The addition of more freedom to the kinematic models will require the addition of more empirical information that can be obtained by spatially resolved data on distant lens galaxies, or from high-quality data (including absorption line shapes) of appropriately selected local elliptical galaxies.

8.3.2. Selection effects of different lens samples

One key pillar in this analysis to improve the precision on the H0 measurement from the TDCOSMO sample is the information on the mass profiles of the SLACS sample. The SLACS sample differs in terms of the redshift distribution and reff/θE relative to the TDCOSMO sample. Beyond our chosen explicit parameterized dependence of the MST parameter λint as a function of reff/θE we do not find trends in the predicted vs measured velocity dispersion within the SLACS sample. However, we do find differences in the external shear contributions between the SLACS and TDCOSMO sample (Shajib et al. 2020b). This is expected because of selection effects. The TDCOSMO sample is composed of quads at higher redshift than SLACS. So it is not surprising that the TDCOSMO lenses tend to be more elongated (to increase the size of the quad cross section) and be more impacted by mass structure along the line of sight than SLACS. Nonetheless, based on previous studies, we have no reason to suspect that the deflectors themselves are intrinsically different between SLACS and TDCOSMO. Complex angular structure of the lenses might also affect the inference in the power-law slope γpl, as the angular degree of freedoms in our model assumptions are, to some degree, limited (Kochanek 2020b). A study with more lenses and particularly sampling the redshift range of the TDCOSMO sample (see Fig. 16) would allow us to better test our current underlying assumption and in case of a significant redshift evolution to correct for it.

8.3.3. Line-of-sight structure

The investigation of the line-of-sight structure of strong gravitational lenses of the TDCOSMO and the SLACS sample follows a specific protocol to provide an individual PDF of the external convergence, p(κext). In our current analysis, the statistical uncertainty of the SLACS line-of-sight structure is subdominant.

In the future – as the other terms of the error budget shrink and this one becomes more relevant – the following steps will be necessary. First, the specific choice of N-body simulation and semi-analytic galaxy evolution model will need to be revisited. Second, it will be necessary to investigate how to improve the comparison with simulation products in order to further mitigate uncertainties. For instance, beyond galaxy number count statistics, weak gravitational lensing observations can also add information on the line-of-sight structure (Tihhonova et al. 2018, 2020).

Ideally, we aim for a validation based on simulations in the full cosmological context. These future simulations should include the presence of the strong lensing deflector, to further quantify nonlinear effects from the line-of-sight structure on the main deflector modeling as well as the main deflector impact on the line-of-sight light path differences (see e.g., Li et al. 2020). Meeting the line-of-sight goal will require large box simulations, and for the main deflector this demands a very high fidelity and resolution at the 10−100 pc scales dominated by baryons in the form of stars and gas.

8.3.4. More flexible lens models and extended hierarchical analysis

Getting the uncertainties right requires careful judgment in the use of theoretical assumptions, validated as much as possible by empirical data. Previous work by TDCOSMO assumed that galaxies were described by power laws or stars plus an NFW profile, leading to a given precision. In this work, we relax this assumption, with the goal to study the impact of the MST. As part of this investigation, we introduce the MST parameter λint in our hierarchical framework and use a PEMD + shear model as baseline. We demonstrate, based on simulations, that these choices are sufficient to the level of precision currently achieved. It is not, however, the end of the story. Additional information will enable better constraints on the mass density profiles. As the precision improves on H0, it will be necessary to keep revisiting our assumptions and validating on a sufficiently large and realistic mock data set.

In the future, additional model flexibility may demand a treatment of more lens model parameters in the full hierarchical context of the inference. Currently, our baseline model is constrained sufficiently by the imaging data of the lensing sample.

However, the development of a hierarchical treatment of additional lensing parameters may also allow us to incorporate lenses with fewer constraints on the lensing nature, such as doubly lensed quasars, or lenses with missing high resolution imaging, or other partially incomplete data products. By pursuing further this development in hierarchical lens modeling, the total number of usable systems can improve, thus, in turn improving the constraints on the Hubble constant.

Substructure adds 0.6%−2% of uncorrelated and un-biased uncertainties on the DΔt inference (Gilman et al. 2020) for individual lenses. Thus, substructure adds a 0.5% uncertainty in quadrature on the combined H0 constraints from the seven TDCOSMO lenses. This effect is highly subdominant to other sources of uncertainties related to the MST in our work and we note that this effect might partially be encapsulated in the scatter in λint, σ(λint), as inferred to be few percent.

8.4. A pathway forward for time-delay cosmography

After having discussed current limitations on the precision and accuracy of our new proposed hierarchical framework applied to time-delay cosmography, we summarize here the key steps to take in the near future, in terms of improvements on the analysis and addition of data, to improve both precision and accuracy in the H0 measurements. Given the new hierarchical context, our largest statistical uncertainty on H0 arises from the stellar anisotropy modeling assumptions and the precision on the velocity dispersion measurements. Multiple and spatially resolved high signal-to-noise velocity dispersion measurements of gravitational lenses are able to further constrain the stellar anisotropy distribution. This can be provided by a large VLT-MUSE and Keck-KCWI campaign of multiple lenses and we expect significant constraining power from JWST (Yıldırım et al. 2020). A complementary approach of studying the mass profile and kinematic structure of the deflector galaxies, is to study the local analogs of those galaxies with high signal-to-noise ratio resolved spectroscopy. Assumptions about potential redshift evolution need to be mitigated and assessed within a lensing sample covering a wide redshift range.

A more straightforward approach in extending our analysis is by incorporating more galaxy-galaxy lenses, in particular lenses that populate a similar distribution to the lensed quasar sample. Such a targeted large sample can reduce potential systematics of our self-similarity assumptions, as well as increase the statistical precision on the mass profiles. Recent searches for strong gravitational lenses in current and ongoing large area imaging surveys, such as the Dark Energy Survey (DES) and the Hyper-Supreme-Cam survey (HSC) have resulted in hundreds of promising galaxy-galaxy scale candidate lenses (see e.g., Jacobs et al. 2019; Sonnenfeld et al. 2020) and dozens of lensed quasars (see e.g., Agnello et al. 2018; Delchambre et al. 2019; Lemon et al. 2020).

With the next generation large ground and space based surveys (Vera Rubin Observatory LSST, Euclid, Nancy Grace Roman Space Telescope), of order 105 galaxy-galaxy lenses and of order 103 quasar-galaxy lenses will be discovered (Oguri & Marshall 2010; Collett 2015). Limited follow-up capabilities with high resolution imaging and spectroscopy will be a key limitation and needs to be mitigated with strategic prioritization of targets to maximize resulting precision and accuracy. We refer to Birrer & Treu (2020) for a forcast based on the precision on H0 we can expect for a current and future lensing sample with spatially resolved kinematics measurement based on the analysis framework presented in this work.

Beyond the addition of external data sets, we emphasize the further demand on the validation of the modeling approach, both in the imaging analysis as well as the stellar anisotropy modeling. Detailed investigation and data challenges based on realistic data with the same complexity level as the real analysis are a useful tool to make progress. To ensure that the requirements are met in the modeling of the deflector galaxy and the local and line-of-sight environment, validation on realistic simulations in the full cosmological context, including selection effects and ray-tracing through the line-of-sight cone of a cosmological box are required. Moreover, we also stress that assessing and tracking systematics at the percent level and the mitigation thereof on the joint inference on H0 would be much facilitated by an automatized and homogenized analysis framework encapsulating all relevant aspects of the analysis of individual lenses.

Finally, a decisive conclusion on the current Hubble tension demands for a rigorous assessment of results by different science collaborations. We stress the importance of conducting the analysis blindly in regard to H0 and related quantities to prevent experimenter bias, a procedure our collaboration has incorporated and followed rigorously. In addition, all measurements of H0 contributing to a decisive conclusion of the tension must guarantee reproducibility. In this work, we provide all software as open-source and release the value-added data products and analysis scripts to the community to facilitate the needed reproducibility.

8.5. Post-blind discussion of the results and comparison with previous time-delay cosmography work

In this section21 we discuss how the measurement presented in this paper related to previous work by members of this collaboration as part of the H0LiCOW, STRIDES, and SHARP projects. We then discuss the relationship between the multiple measurements obtained within the hierarchical framework introduced in this paper. All the relevant measurements are summarized in Fig. 18 for quick visualization.

thumbnail Fig. 18.

Comparison of different blind H0 measurements by the TDCOSMO collaboration, based on different mass profile assumptions and data sets incorporated. All measurements presented on this plot were performed blindly with regard to the inference of H0. The measurement on top is the combined H0LiCOW six lenses constraints presented by Wong et al. (2020), when averaging power-law and composite NFW plus stars (with constant mass-to-light ratio) on a lens-by-lens basis without correlated errors among the lenses. The next two measurements are from Millon et al. (2020) of six TDCOSMO time-delay lenses (five H0LiCOW lenses22 and one STRIDES lens by Shajib et al. 2020a), when performing the inference assuming either a composite NFW plus stars (with constant mass-to-light ratio) or the power-law mass density profile for the galaxy acting as a lens. Lower panel: results from this work. The main difference with respect to previous work is that we have made virtually no assumption on the radial mass density profile of the lens galaxy, and taken into account the covariance between the lenses. The analysis in this work is constrained only by the stellar kinematics and fully accounts for the uncertainty related to the mass sheet transformation (MST). In this framework, we obtain four measurements according to the datasets considered. The TDCOSMO-only inference is based on the same set of seven lenses as those jointly included by Millon et al. (2020) and Wong et al. (2020). The inferred median value is the same, indicating no bias, and the uncertainties, as expected, are larger. The next three measurements rely on external datasets from the SLACS survey, by making the assumption that the lens galaxies in the two surveys are drawn from the same population. The TDCOSMO+SLACSIFU measurements uses, in addition to the TDCOSMO sample, nine lenses from the SLACS sample with IFU observations to inform the anisotropy prior applied on the TDCOSMO lenses. The TDCOSMO+SLACSSDSS measurement comes from the joint analysis of the TDCOSMO sample and 33 SLACS lenses with SDSS spectroscopy. The TDCOSMO+SLACSSDSS + IFU presents the joint analysis of all three data sets, again assuming self-similar distributions of the mass profiles and stellar anisotropy. The TDCOSMO-only and TDCOSMO+SLACSIFU analyses do not rely on self-similar mass profiles of the SLACS and TDCOSMO sample while the TDCOSMO+SLACSSDSS and TDCOSMO+SLACSSDSS + IFU measurements (orange and purple) do. All the measurements shown in this plot are in statistical agreement with each other. See Sect. 8.5 for a discussion and physical interpretation of the results (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/tdcosmo_comparison_plot.ipynb).

Open with DEXTER

The result of our hierarchical TDCOSMO-only analysis is fully consistent with the assumptions on the mass profiles made in previous H0LiCOW/STRIDES/SHARP work (see e.g., Wong et al. 2020; Shajib et al. 2020a; Millon et al. 2020). The consistency is reinforced by Yang et al. (2020) who concluded that the combination of kinematics and time-delay constraints are consistent with General Relativity, an underlying assumptions of time-delay cosmography. The only difference with respect to the H0LiCOW/STRIDES/SHARP analysis is that the uncertainty has significantly increased. This was expected, because we have virtually eliminated the assumptions on the radial mass profile of elliptical galaxies and, due to the MST, the only source of information left to enable a H0 measurement is the stellar kinematics. Without lensing information, due to the well known mass-anisotropy degeneracy, unresolved kinematics has limited power to constrain the mass profiles. Since our parametrization is maximally degenerate with H0 and our assumptions are minimal, this 9% error budget accounts for potential effects of the MST.

Another set of results is obtained within the hierarchical framework with the addition of external information. Under the additional assumption that the galaxies in the external datasets are drawn from the same population as the TDCOSMO deflectors, these results achieve higher precision than TDCOSMO alone. Adding the SLACS dataset shrinks the uncertainty to 5% and shifts the mean inferred H0 to a value about 6 km s−1 Mpc−1 lower than the TDCOSMO-only analysis. This shift is consistent within the uncertainties achieved by the TDCOSMO-only analysis and can be traced back to two factors: (i) the anisotropy constraints prefer a lower aani value and this moves H0 down relative to the chosen prior on aani. The VIMOS+IFU inference is about 2 km s−1 Mpc−1 lower than the equivalent TDCOSMO-only inference. (ii) The SLACS lenses prefer an overall lower – but statistically consistent – λint, 0 value for a given anisotropy model by about 8%. The negative trend of λint with reff/θE (αλ) partially mitigates an even lower λint value preferred by the SLACS sample relative to the TDCOSMO sample.

The shift between the TDCOSMO and TDCOSMO+SLACS results can have two possible explanations (if it is not purely a statistical fluctuation). One option is that elliptical galaxies are more radially anisotropic (and therefore have a flatter mass density profile to reproduce the same velocity dispersion profile) than the prior used to model the TDCOSMO galaxies. The alternative option is that the TDCOSMO and SLACS galaxies are somehow different. Within the observables at disposal, one that may be indicative of a different line of sight anisotropy is the higher ellipticity of the surface brightness and of the projected total mass distribution (Shajib et al. 2020b) of the TDCOSMO deflectors in comparison to the SLACS deflectors. As mentioned in Sect. 6.3, this is understood to be a selection effect because ellipticity increases the cross section for quadruple images and TDCOSMO is a sample of mostly quads (six out of seven), while SLACS is mostly doubles (Treu et al. 2009). Departure from spherical symmetry in elliptical galaxies can arise from rotation or anisotropy. If flattening arises from rotation (which we have neglected in our study) more flattened systems are more likely to be seen edge-on. If it arises from anisotropy, the observed flattening could be due to tangential anisotropy that is not included in our models, or to a smaller degree of radial anisotropy than for other orientations. These two options result in different predictions that can be tested with spatially resolved kinematics of the TDCOSMO lens galaxies. If the shift is just due to an inconsistency between the TDCOSMO prior and the SLACS likelihood, spatially resolved kinematics will bring them in closer alignment. If it is due to intrinsic differences, spatially resolved kinematics will reveal rotation or tangential (less radial) anisotropy. In addition, spatially resolved kinematics of the TDCOSMO sample will reduce the uncertainties of both measurement, and thus resolve whether the shift is a fluctuation or significant.

The other potential way to elucidate the marginal differences between the TDCOSMO and SLACS sample is to obtain precise measurements of mass at scales well beyond the Einstein radius. As seen in Fig. 17, a pure power law and the transformed profile differ by up to 50% in that region (depending on the choice of Rc). Satellite kinematics or weak lensing would help reduce the freedom of the MST, provided they reach sufficient precision.

9. Conclusion

The precision of time-delay cosmography has improved significantly in the past few years, driven by improvement in the quality of the data and methodology. As the precision improves it is critical to revisit assumptions and explore potential systematics, while charting the way forward.

In this work, we relaxed previous assumptions on the mass-profile parameterization and introduced an efficient way to explore potential systematics associated to the mass-sheet degeneracy in a hierarchical Bayesian analysis. In this new approach, the mass density profile of the lens galaxies is only constrained by basic information on stellar kinematics. It thus provides a conservative estimate of how much the mass profile can depart from a power law, and how much the error budget can grow as a result. Based on the consistent results of the power law and stars plus NFW profiles in the inference on H0 (Millon et al. 2020), we expect very similar conclusions had we performed this analysis with a stars plus NFW profile.

We validated our approach on the Time-Delay Lens Modeling Challenge sample of hydrodynamical simulations. We then applied the formalism and assumptions to the TDCOSMO data set in a blind fashion. Based on the TDCOSMO data set alone we infer km s−1 Mpc−1. The uncertainties on H0 are dominated by the precision of the spectroscopic data and the modeling uncertainties therein. To further increase our precision, we added self-consistently to our analysis a set of SLACS lenses with imaging modeling and independent kinematic constraints. We characterized the candidate lenses to be added and explicitly selected only lenses that do not have significantly enhanced local environments. In total, we were able to add 33 additional lenses with no time delay information of which nine have additional 2D kinematics with VIMOS IFU data that allowed us to further constrain uncertainties in the anisotropy profile of the stellar orbits. Our most constrained measurement of the Hubble constant is km s−1 Mpc−1 from the joint TDCOSMO+SLACS analysis, assuming that the two samples are drawn from the same population.

The 5% error budget reported in this work addresses conclusively concerns about the MST (Schneider & Sluse 2013; Sonnenfeld 2018; Kochanek 2020a, 2020b). If the mass density profiles of lens galaxies are not well described by power-laws or stars plus NFW halos, this is the appropriate uncertainty to associate with current time-delay cosmography. Additional effects are very much subdominant for now as compared with the effect of the MST. For example, the small level of pixelated corrections to the elliptical power-law model obtained in our previous work suggests that the departure from ellipticity is not required by the data.

Based on the methodology presented and the results achieved, we lay out a roadmap for further improvements to ultimately enable a 1% precision measurement of the Hubble constant, which is a clear target both for resolving the Hubble tension and to serve as a prior on dark energy studies (Weinberg et al. 2013). The key ingredients required to reduce the statistical uncertainties are (i) spatially resolved high signal-to-noise kinematic measurements; (ii) an increase in the sample size of both lenses with measured time-delays and lenses with high-resolution imaging and precise kinematic measurements. Potential sources of systematic that should be investigated further to maintain accuracy at the target precision are those arising from: (i) measurements of the stellar velocity dispersion; (ii) characterization of the selection function and local environment of all the lenses included in the inference; (iii) mass profile modeling assumptions beyond the MST and stellar anisotropy modeling assumptions.

Upcoming deep, wide-field surveys (such as those enabled by Vera Rubin Observatory, Euclid and the Nancy Grace Roman Observatory) will discover many thousands of lenses of which several hundred will have accurate time delay measurements (see e.g., Oguri & Marshall 2010; Collett 2015; Huber et al. 2019). The analysis framework presented in this work will serve as a baseline for the analysis of these giant samples of lenses; simultaneously enabling precise and accurate constraints on the Hubble constant and the astrophysics of strong lensing galaxies.


2

For the NFW profile parameters, priors on the mass-concentration relation were imposed on the individual analyses.

3

See Birrer et al. (2016) for an analysis explicitly constraining the MST with kinematic data that satisfies the error budget of Kochanek (2020a).

4

Noting however the caveats on the realism of the TDLMC simulations discussed by Ding et al. (2020).

9

The integral between the deflector and the source deviates from the Born approximation as the light paths are significantly perturbed (see e.g., Bar-Kana 1996; Birrer et al. 2017).

10

We note that the mean cosmological background density is already fully encompassed in the background metric and we effectively only require to model the enhancement matter density (see e.g., Wucknitz 2008; Birrer et al. 2017).

12

Alternatively, we could have also transformed the DΔt posteriors accordingly to account for the external convergence for each individual lens.

13

For the composite models, priors on the mass-concentration relation of the dark matter profiles were imposed.

14

The first lens, B1608+656, and the reanalysis of RXJ1131−1231 with AO data were not executed blindly.

15

The Millenium Simulation uses the following flat ΛCDM cosmology: Ωm = 0.25, Ωb = 0.045, H0 = 73 km s−1 Mpc−1, n = 1, and σ8 = 0.9.

16

We note that the prior on the power-law slope γpl is chosen to be uniform in [1,3] during the Bayesian inference with MCMC.

17

To use the IFU data set more optimally, we add the lens SDSSJ0216-0813, which is the remaining lens within the IFU quality sample that was not selected by Shajib et al. (2020b) from the original SLACS sample.

18

The notebooks are publicly available and we facilitate the use of different priors and cosmological models. All choices presented here are made blindly in regard to H0.

19

This section, with the exception of Sect. 8.5, was written before the results of the combined TDCOSMO+SLACS analysis were known to the authors and, thus, reflect the assessment of uncertainties present in our analysis agnostic to its outcome.

20

This statement has been tested with a flat prior on σ(λint).

21

This section was written after the results were known to the authors.

22

Excluding B1608+656 as this lens was only analyzed with a power-law model and not with a composite model and thus not part of the model comparison analysis. Additional lensing potential perturbations on top of the power-law profile lead to only small amounts of corrections Suyu et al. (2010).

Acknowledgments

SB thanks Kfir Blum, Veronica Motta, Timo Anguita, Sampath Mukherjee, Elizabeth Buckley-Geer for useful discussions, Hyungsuk Tak for participating in the TDLMC, the TDLMC team for setting up the challenge and Yiping Shu for providing access to the SDSS velocity dispersion measurements. AJS was supported by the Dissertation Year Fellowship from the UCLA Graduate Division. TT and AJS acknowledge support from NSF through NSF grant NSF-AST-1906976, from NASA through grant HST-GO-15320 and from the Packard Foundation through a Packard Research Fellowship to TT. AA was supported by a grant from VILLUM FONDEN (project number 16599). This project is partially funded by the Danish council for independent research under the project “Fundamentals of Dark Matter Structures”, DFF–6108-00470. SB and MWA acknowledges support from the Kavli Foundation. TC is funded by a Royal Astronomical Society Research Fellowship. C.D.F. and G.C.-F.C. acknowledge support for this work from the National Science Foundation under grant no. AST-1907396 and from NASA through grant HST-GO-15320. This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (COSMICLENS:grant agreement No 787886). CS is supported by a Hintze Fellow at the Oxford Centre for Astrophysical Surveys, which is funded through generous support from the Hintze Family Charitable Foundation. SHS thanks the Max Planck Society for support through the Max Planck Research Group. This work was supported by World Premier International Research Center Initiative (WPI Initiative), MEXT, Japan. This work was supported by JSPS KAKENHI Grant Number JP20K14511. SB acknowledges the hospitality of the Munich Institute for Astro- and Particle Physics (MIAPP) of the Excellence Cluster “Universe”. Based on observations collected at the European Southern Observatory under ESO program 0102.A-0600 (PI Agnello), 075.B-0226 (PI Koopmans), 177.B-0682 (PI Koopmans). This work made use of the following public software packages: HIERARC (this work), LENSTRONOMY (Birrer et al. 2015; Birrer & Amara 2018), DOLPHIN (Shajib et al. 2020b), EMCEE (Foreman-Mackey et al. 2013), CORNER (Foreman-Mackey 2016), PPXF (Cappellari 2012), ASTROPY (Astropy Collaboration 2013, 2018), FASTELL (Barkana 1999) and standard Python libraries.

References

  1. Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2017, Nature, 551, 85 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  2. Abbott, T. M. C., Abdalla, F. B., Annis, J., et al. 2018, MNRAS, 480, 3879 [NASA ADS] [CrossRef] [Google Scholar]
  3. Agnello, A., Evans, N. W., & Romanowsky, A. J. 2014a, MNRAS, 442, 3284 [NASA ADS] [CrossRef] [Google Scholar]
  4. Agnello, A., Evans, N. W., Romanowsky, A. J., & Brodie, J. P. 2014b, MNRAS, 442, 3299 [NASA ADS] [CrossRef] [Google Scholar]
  5. Agnello, A., Sonnenfeld, A., Suyu, S. H., et al. 2016, MNRAS, 458, 3830 [NASA ADS] [CrossRef] [Google Scholar]
  6. Agnello, A., Lin, H., Kuropatkin, N., et al. 2018, MNRAS, 479, 4345 [NASA ADS] [CrossRef] [Google Scholar]
  7. Astropy Collaboration (Robitaille, T. P., et al.) 2013, A&A, 558, A33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  8. Astropy Collaboration (Price-Whelan, A. M., et al.) 2018, AJ, 156, 123 [Google Scholar]
  9. Auger, M. W., Treu, T., Bolton, A. S., et al. 2009, ApJ, 705, 1099 [NASA ADS] [CrossRef] [Google Scholar]
  10. Auger, M. W., Treu, T., Bolton, A. S., et al. 2010, ApJ, 724, 511 [NASA ADS] [CrossRef] [Google Scholar]
  11. Bar-Kana, R. 1996, ApJ, 468, 17 [NASA ADS] [CrossRef] [Google Scholar]
  12. Barkana, R. 1998, ApJ, 502, 531 [NASA ADS] [CrossRef] [Google Scholar]
  13. Barkana, R. 1999, Astrophysics Source Code Library [record ascl:9910.003] [Google Scholar]
  14. Barnabè, M., & Koopmans, L. V. E. 2007, ApJ, 666, 726 [NASA ADS] [CrossRef] [Google Scholar]
  15. Barnabè, M., Czoske, O., Koopmans, L. V. E., et al. 2009, MNRAS, 399, 21 [NASA ADS] [CrossRef] [Google Scholar]
  16. Barnabè, M., Czoske, O., Koopmans, L. V. E., et al. 2011, MNRAS, 415, 2215 [NASA ADS] [CrossRef] [Google Scholar]
  17. Barnabè, M., Dutton, A. A., Marshall, P. J., et al. 2012, MNRAS, 423, 1073 [NASA ADS] [CrossRef] [Google Scholar]
  18. Bernardi, M., Domínguez-Sánchez, H., Margalef-Bentabol, B., Nikakhtar, F., & Sheth, R. K. 2020, MNRAS, 494, 5148 [CrossRef] [Google Scholar]
  19. Betoule, M., Kessler, R., Guy, J., et al. 2014, A&A, 568, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Binney, J., & Mamon, G. A. 1982, MNRAS, 200, 361 [NASA ADS] [CrossRef] [Google Scholar]
  21. Binney, J., & Tremaine, S. 2008, Galactic Dynamics: Second Edition (Princeton: Princeton University Press) [Google Scholar]
  22. Birrer, S., & Amara, A. 2018, Phys. Dark Univ., 22, 189 [NASA ADS] [CrossRef] [Google Scholar]
  23. Birrer, S., & Treu, T. 2019, MNRAS, 489, 2097 [CrossRef] [Google Scholar]
  24. Birrer, S., & Treu, T. 2020, A&A, submitted [arXiv:2008.06157] [Google Scholar]
  25. Birrer, S., Amara, A., & Refregier, A. 2015, ApJ, 813, 102 [NASA ADS] [CrossRef] [Google Scholar]
  26. Birrer, S., Amara, A., & Refregier, A. 2016, JCAP, 2016, 020 [NASA ADS] [CrossRef] [Google Scholar]
  27. Birrer, S., Welschen, C., Amara, A., & Refregier, A. 2017, JCAP, 2017, 049 [CrossRef] [Google Scholar]
  28. Birrer, S., Treu, T., Rusu, C. E., et al. 2019, MNRAS, 484, 4726 [NASA ADS] [CrossRef] [Google Scholar]
  29. Blandford, R., & Narayan, R. 1986, ApJ, 310, 568 [NASA ADS] [CrossRef] [Google Scholar]
  30. Blum, K., Castorina, E., & Simonović, M. 2020, ApJ, 892, L27 [CrossRef] [Google Scholar]
  31. Bolton, A. S., Burles, S., Schlegel, D. J., Eisenstein, D. J., & Brinkmann, J. 2004, AJ, 127, 1860 [NASA ADS] [CrossRef] [Google Scholar]
  32. Bolton, A. S., Burles, S., Koopmans, L. V. E., Treu, T., & Moustakas, L. A. 2006, ApJ, 638, 703 [NASA ADS] [CrossRef] [Google Scholar]
  33. Bolton, A. S., Burles, S., Koopmans, L. V. E., et al. 2008, ApJ, 682, 964 [NASA ADS] [CrossRef] [Google Scholar]
  34. Bonvin, V., Tewes, M., Courbin, F., et al. 2016, A&A, 585, A88 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  35. Bonvin, V., Tihhonova, O., Millon, M., et al. 2018, A&A, 616, A183 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Bonvin, V., Millon, M., Chan, J. H.-H., et al. 2019, A&A, 629, A97 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  37. Braibant, L., Hutsemékers, D., Sluse, D., et al. 2014, A&A, 565, L11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  38. Buckley-Geer, E. J., Lin, H., Rusu, C. E., et al. 2020, MNRAS, 498, 3241 [CrossRef] [Google Scholar]
  39. Cappellari, M. 2008, MNRAS, 390, 71 [NASA ADS] [CrossRef] [Google Scholar]
  40. Cappellari, M. 2012, Astrophysics Source Code Library [record ascl:1210.002] [Google Scholar]
  41. Cappellari, M. 2017, MNRAS, 466, 798 [NASA ADS] [CrossRef] [Google Scholar]
  42. Cappellari, M., & Emsellem, E. 2004, PASP, 116, 138 [NASA ADS] [CrossRef] [Google Scholar]
  43. Cappellari, M., Emsellem, E., Bacon, R., et al. 2007, MNRAS, 379, 418 [NASA ADS] [CrossRef] [Google Scholar]
  44. Chen, G. C. F., Fassnacht, C. D., Suyu, S. H., et al. 2019, MNRAS, 490, 1743 [NASA ADS] [CrossRef] [Google Scholar]
  45. Coles, J. 2008, ApJ, 679, 17 [NASA ADS] [CrossRef] [Google Scholar]
  46. Coles, J. P., Read, J. I., & Saha, P. 2014, MNRAS, 445, 2181 [NASA ADS] [CrossRef] [Google Scholar]
  47. Collett, T. E. 2015, ApJ, 811, 20 [NASA ADS] [CrossRef] [Google Scholar]
  48. Collett, T. E., & Cunnington, S. D. 2016, MNRAS, 462, 3255 [NASA ADS] [CrossRef] [Google Scholar]
  49. Collett, T. E., Marshall, P. J., Auger, M. W., et al. 2013, MNRAS, 432, 679 [NASA ADS] [CrossRef] [Google Scholar]
  50. Courbin, F., Chantry, V., Revaz, Y., et al. 2011, A&A, 536, A53 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Courbin, F., Bonvin, V., Buckley-Geer, E., et al. 2018, A&A, 609, A71 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Cretton, N., de Zeeuw, P. T., van der Marel, R. P., & Rix, H.-W. 1999, ApJS, 124, 383 [NASA ADS] [CrossRef] [Google Scholar]
  53. Czoske, O., Barnabè, M., Koopmans, L. V. E., Treu, T., & Bolton, A. S. 2008, MNRAS, 384, 987 [NASA ADS] [CrossRef] [Google Scholar]
  54. Czoske, O., Barnabè, M., Koopmans, L. V. E., Treu, T., & Bolton, A. S. 2012, MNRAS, 419, 656 [CrossRef] [Google Scholar]
  55. Dawson, K. S., Schlegel, D. J., Ahn, C. P., et al. 2013, AJ, 145, A10 [NASA ADS] [CrossRef] [Google Scholar]
  56. De Lucia, G., & Blaizot, J. 2007, MNRAS, 375, 2 [NASA ADS] [CrossRef] [Google Scholar]
  57. Delchambre, L., Krone-Martins, A., Wertz, O., et al. 2019, A&A, 622, A165 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  58. Dey, A., Schlegel, D. J., Lang, D., et al. 2019, AJ, 157, 168 [NASA ADS] [CrossRef] [Google Scholar]
  59. Diehl, H. T., Buckley-Geer, E. J., Lindgren, K. A., et al. 2017, ApJS, 232, 15 [CrossRef] [Google Scholar]
  60. Ding, X., Treu, T., Shajib, A. J., et al. 2018, MNRAS, submitted [arXiv:1801.01506] [Google Scholar]
  61. Ding, X., Treu, T., Birrer, S., et al. 2020, MNRAS, submitted [arXiv:2006.08619] [Google Scholar]
  62. Dobler, G., Keeton, C. R., Bolton, A. S., & Burles, S. 2008, ApJ, 685, 57 [NASA ADS] [CrossRef] [Google Scholar]
  63. Dobler, G., Fassnacht, C. D., Treu, T., et al. 2015, ApJ, 799, 168 [NASA ADS] [CrossRef] [Google Scholar]
  64. Eulaers, E., Tewes, M., Magain, P., et al. 2013, A&A, 553, A121 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Faber, S. M., & Jackson, R. E. 1976, ApJ, 204, 668 [NASA ADS] [CrossRef] [Google Scholar]
  66. Fadely, R., Keeton, C. R., Nakajima, R., & Bernstein, G. M. 2010, ApJ, 711, 246 [NASA ADS] [CrossRef] [Google Scholar]
  67. Falco, E. E., Gorenstein, M. V., & Shapiro, I. I. 1985, ApJ, 289, L1 [NASA ADS] [CrossRef] [Google Scholar]
  68. Fassnacht, C. D., Womble, D. S., Neugebauer, G., et al. 1996, ApJ, 460, L103 [NASA ADS] [CrossRef] [Google Scholar]
  69. Fassnacht, C. D., Pearson, T. J., Readhead, A. C. S., et al. 1999, ApJ, 527, 498 [NASA ADS] [CrossRef] [Google Scholar]
  70. Fassnacht, C. D., Xanthopoulos, E., Koopmans, L. V. E., & Rusin, D. 2002, ApJ, 581, 823 [NASA ADS] [CrossRef] [Google Scholar]
  71. Fassnacht, C. D., Koopmans, L. V. E., & Wong, K. C. 2011, MNRAS, 410, 2167 [NASA ADS] [CrossRef] [Google Scholar]
  72. Foreman-Mackey, D. 2016, J. Open Sour. Softw., 1, 24 [NASA ADS] [CrossRef] [Google Scholar]
  73. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [NASA ADS] [CrossRef] [Google Scholar]
  74. Freedman, W. L., Madore, B. F., Hatt, D., et al. 2019, ApJ, 882, 34 [NASA ADS] [CrossRef] [Google Scholar]
  75. Freedman, W. L., Madore, B. F., Hoyt, D., et al. 2020, ApJ, 891, 57 [CrossRef] [Google Scholar]
  76. Gerhard, O., Kronawitter, A., Saglia, R. P., & Bender, R. 2001, AJ, 121, 1936 [NASA ADS] [CrossRef] [Google Scholar]
  77. Gilman, D., Birrer, S., & Treu, T. 2020, A&A, 642, A194 [CrossRef] [EDP Sciences] [Google Scholar]
  78. Gorenstein, M. V., Falco, E. E., & Shapiro, I. I. 1988, ApJ, 327, 693 [NASA ADS] [CrossRef] [Google Scholar]
  79. Greene, J. E., Murphy, J. D., Graves, G. J., et al. 2013, ApJ, 776, 64 [NASA ADS] [CrossRef] [Google Scholar]
  80. Grogin, N. A., & Narayan, R. 1996, ApJ, 464, 92 [NASA ADS] [CrossRef] [Google Scholar]
  81. Hernquist, L. 1993, ApJ, 409, 548 [NASA ADS] [CrossRef] [Google Scholar]
  82. Hilbert, S., Hartlap, J., White, S. D. M., & Schneider, P. 2009, A&A, 499, 31 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  83. Holder, G. P., & Schechter, P. L. 2003, ApJ, 589, 688 [NASA ADS] [CrossRef] [Google Scholar]
  84. Huang, C. D., Riess, A. G., Yuan, W., et al. 2020, ApJ, 889, 5 [NASA ADS] [CrossRef] [Google Scholar]
  85. Huber, S., Suyu, S. H., Noebauer, U. M., et al. 2019, A&A, 631, A161 [CrossRef] [EDP Sciences] [Google Scholar]
  86. Huterer, D., Keeton, C. R., & Ma, C.-P. 2005, ApJ, 624, 34 [NASA ADS] [CrossRef] [Google Scholar]
  87. Jacobs, C., Collett, T., Glazebrook, K., et al. 2019, ApJS, 243, 17 [CrossRef] [Google Scholar]
  88. Jee, I., Komatsu, E., & Suyu, S. H. 2015, JCAP, 2015, 033 [NASA ADS] [CrossRef] [Google Scholar]
  89. Keeton, C. R., & Kochanek, C. S. 1997, ApJ, 487, 42 [NASA ADS] [CrossRef] [Google Scholar]
  90. Keeton, C. R., & Moustakas, L. A. 2009, ApJ, 699, 1720 [NASA ADS] [CrossRef] [Google Scholar]
  91. Kochanek, C. S. 2002, ApJ, 578, 25 [NASA ADS] [CrossRef] [Google Scholar]
  92. Kochanek, C. S. 2003, ApJ, 583, 49 [NASA ADS] [CrossRef] [Google Scholar]
  93. Kochanek, C. S. 2006, in Saas-Fee Advanced Course 33: Gravitational Lensing: Strong, Weak and Micro, eds. G. Meylan, P. Jetzer, P. North, et al. (Berlin: Springer), 91 [NASA ADS] [CrossRef] [Google Scholar]
  94. Kochanek, C. S. 2020a, MNRAS, 493, 1725 [CrossRef] [Google Scholar]
  95. Kochanek, C. S. 2020b, MNRAS, submitted [arXiv:2003.08395] [Google Scholar]
  96. Koopmans, L. V. E. 2004, ArXiv e-prints [arXiv:astro-ph/0412596] [Google Scholar]
  97. Koopmans, L. V. E., Treu, T., Fassnacht, C. D., Blandford, R. D., & Surpi, G. 2003, ApJ, 599, 70 [NASA ADS] [CrossRef] [Google Scholar]
  98. Kormann, R., Schneider, P., & Bartelmann, M. 1994, A&A, 284, 285 [NASA ADS] [Google Scholar]
  99. Lang, D., Hogg, D. W., & Mykytyn, D. 2016, Astrophysics Source Code Library [record ascl:1604.008] [Google Scholar]
  100. Lemon, C., Auger, M. W., McMahon, R., et al. 2020, MNRAS, 494, 3491 [CrossRef] [Google Scholar]
  101. Li, N., Becker, C., & Dye, S. 2020, MNRAS, submitted [arXiv:2006.08540] [Google Scholar]
  102. Liao, K., Treu, T., Marshall, P., et al. 2015, ApJ, 800, 11 [NASA ADS] [CrossRef] [Google Scholar]
  103. Lin, H., Buckley-Geer, E., Agnello, A., et al. 2017, ApJ, 838, L15 [NASA ADS] [CrossRef] [Google Scholar]
  104. Mao, S., & Schneider, P. 1998, MNRAS, 295, 587 [NASA ADS] [CrossRef] [Google Scholar]
  105. McCully, C., Keeton, C. R., Wong, K. C., & Zabludoff, A. I. 2014, MNRAS, 443, 3631 [NASA ADS] [CrossRef] [Google Scholar]
  106. Merritt, D. 1985, AJ, 90, 1027 [NASA ADS] [CrossRef] [Google Scholar]
  107. Millon, M., Galan, A., Courbin, F., et al. 2020, A&A, 639, A101 [CrossRef] [EDP Sciences] [Google Scholar]
  108. Morgan, N. D., Caldwell, J. A. R., Schechter, P. L., et al. 2004, AJ, 127, 2617 [NASA ADS] [CrossRef] [Google Scholar]
  109. Myers, S. T., Fassnacht, C. D., Djorgovski, S. G., et al. 1995, ApJ, 447, L5 [NASA ADS] [CrossRef] [Google Scholar]
  110. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1997, ApJ, 490, 493 [NASA ADS] [CrossRef] [Google Scholar]
  111. Nipoti, C., Londrillo, P., & Ciotti, L. 2006, MNRAS, 370, 681 [NASA ADS] [CrossRef] [Google Scholar]
  112. Oguri, M. 2007, ApJ, 660, 1 [NASA ADS] [CrossRef] [Google Scholar]
  113. Oguri, M., & Marshall, P. J. 2010, MNRAS, 405, 2579 [NASA ADS] [Google Scholar]
  114. Oguri, M., Inada, N., Hennawi, J. F., et al. 2005, ApJ, 622, 106 [NASA ADS] [CrossRef] [Google Scholar]
  115. Osipkov, L. P. 1979, Pisma v Astronomicheskii Zhurnal, 5, 77 [Google Scholar]
  116. Paraficz, D., & Hjorth, J. 2009, A&A, 507, L49 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  117. Pesce, D. W., Braatz, J. A., Reid, M. J., et al. 2020, ApJ, 891, L1 [CrossRef] [Google Scholar]
  118. Philcox, O. H. E., Ivanov, M. M., Simonović, M., & Zaldarriaga, M. 2020, JCAP, 2020, 032 [CrossRef] [Google Scholar]
  119. Planck Collaboration VI. 2020, A&A, 641, A6 [CrossRef] [EDP Sciences] [Google Scholar]
  120. Posacki, S., Cappellari, M., Treu, T., Pellegrini, S., & Ciotti, L. 2015, MNRAS, 446, 493 [NASA ADS] [CrossRef] [Google Scholar]
  121. Rathna Kumar, S., Stalin, C. S., & Prabhu, T. P. 2015, A&A, 580, A38 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  122. Read, J. I., Saha, P., & Macciò, A. V. 2007, ApJ, 667, 645 [NASA ADS] [CrossRef] [Google Scholar]
  123. Refsdal, S. 1964, MNRAS, 128, 307 [NASA ADS] [CrossRef] [Google Scholar]
  124. Riess, A. G., Casertano, S., Yuan, W., Macri, L. M., & Scolnic, D. 2019, ApJ, 876, 85 [NASA ADS] [CrossRef] [Google Scholar]
  125. Romanowsky, A. J., & Kochanek, C. S. 1999, ApJ, 516, 18 [NASA ADS] [CrossRef] [Google Scholar]
  126. Rusu, C. E., Fassnacht, C. D., Sluse, D., et al. 2017, MNRAS, 467, 4220 [NASA ADS] [CrossRef] [Google Scholar]
  127. Rusu, C. E., Wong, K. C., Bonvin, V., et al. 2020, MNRAS, 498, 1440 [CrossRef] [Google Scholar]
  128. Saha, P., & Williams, L. L. R. 2006, ApJ, 653, 936 [NASA ADS] [CrossRef] [Google Scholar]
  129. Saha, P., Coles, J., Macciò, A. V., & Williams, L. L. R. 2006, ApJ, 650, L17 [NASA ADS] [CrossRef] [Google Scholar]
  130. Schechter, P. L., Bailyn, C. D., Barr, R., et al. 1997, ApJ, 475, L85 [NASA ADS] [CrossRef] [Google Scholar]
  131. Schneider, P. 1985, A&A, 143, 413 [Google Scholar]
  132. Schneider, P., & Sluse, D. 2013, A&A, 559, A37 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  133. Schneider, P., & Sluse, D. 2014, A&A, 564, A103 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  134. Schneider, P., Ehlers, J., & Falco, E. E. 1992, Gravitational Lenses (Berlin, Heidelberg: Springer-Verlag), 560 [Google Scholar]
  135. Schwarzschild, M. 1979, ApJ, 232, 236 [NASA ADS] [CrossRef] [Google Scholar]
  136. Scolnic, D. M., Jones, D. O., Rest, A., et al. 2018, ApJ, 859, 101 [NASA ADS] [CrossRef] [Google Scholar]
  137. Sereno, M., & Paraficz, D. 2014, MNRAS, 437, 600 [NASA ADS] [CrossRef] [Google Scholar]
  138. Shajib, A. J., Treu, T., & Agnello, A. 2018, MNRAS, 473, 210 [NASA ADS] [CrossRef] [Google Scholar]
  139. Shajib, A. J., Birrer, S., Treu, T., et al. 2019, MNRAS, 483, 5649 [NASA ADS] [CrossRef] [Google Scholar]
  140. Shajib, A. J., Birrer, S., Treu, T., et al. 2020a, MNRAS, 494, 6072 [CrossRef] [Google Scholar]
  141. Shajib, A. J., Treu, T., Birrer, S., & Sonnenfeld, A. 2020b, MNRAS, submitted [arXiv:2008.11724] [Google Scholar]
  142. Shu, Y., Bolton, A. S., Brownstein, J. R., et al. 2015, ApJ, 803, 71 [NASA ADS] [CrossRef] [Google Scholar]
  143. Sluse, D., Surdej, J., Claeskens, J.-F., et al. 2003, A&A, 406, L43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  144. Sluse, D., Chantry, V., Magain, P., Courbin, F., & Meylan, G. 2012, A&A, 538, A99 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  145. Sluse, D., Rusu, C. E., Fassnacht, C. D., et al. 2019, MNRAS, 490, 613 [CrossRef] [Google Scholar]
  146. Sonnenfeld, A. 2018, MNRAS, 474, 4648 [NASA ADS] [CrossRef] [Google Scholar]
  147. Sonnenfeld, A., Verma, A., More, A., et al. 2020, A&A, 642, A148 [CrossRef] [EDP Sciences] [Google Scholar]
  148. Spiniello, C., Koopmans, L. V. E., Trager, S. C., et al. 2015, MNRAS, 452, 2434 [NASA ADS] [CrossRef] [Google Scholar]
  149. Springel, V., White, S. D. M., Jenkins, A., et al. 2005, Nature, 435, 629 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  150. Suyu, S. H., Marshall, P. J., Hobson, M. P., & Blandford, R. D. 2006, MNRAS, 371, 983 [NASA ADS] [CrossRef] [Google Scholar]
  151. Suyu, S. H., Marshall, P. J., Blandford, R. D., et al. 2009, APJ, 691, 277 [NASA ADS] [CrossRef] [Google Scholar]
  152. Suyu, S. H., Marshall, P. J., Auger, M. W., et al. 2010, ApJ, 711, 201 [NASA ADS] [CrossRef] [Google Scholar]
  153. Suyu, S. H., Auger, M. W., Hilbert, S., et al. 2013, ApJ, 766, 70 [NASA ADS] [CrossRef] [Google Scholar]
  154. Suyu, S. H., Treu, T., Hilbert, S., et al. 2014, ApJ, 788, L35 [NASA ADS] [CrossRef] [Google Scholar]
  155. Suyu, S. H., Bonvin, V., Courbin, F., et al. 2017, MNRAS, 468, 2590 [NASA ADS] [CrossRef] [Google Scholar]
  156. Taubenberger, S., Suyu, S. H., Komatsu, E., et al. 2019, A&A, 628, L7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  157. Tewes, M., Courbin, F., Meylan, G., et al. 2012, The Messenger, 150, 49 [NASA ADS] [Google Scholar]
  158. Tewes, M., Courbin, F., & Meylan, G. 2013, A&A, 553, A120 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  159. Tihhonova, O., Courbin, F., Harvey, D., et al. 2018, MNRAS, 477, 5657 [NASA ADS] [CrossRef] [Google Scholar]
  160. Tihhonova, O., Courbin, F., Harvey, D., et al. 2020, MNRAS, 498, 1406 [CrossRef] [Google Scholar]
  161. Tonry, J. L. 1998, AJ, 115, 1 [NASA ADS] [CrossRef] [Google Scholar]
  162. Treu, T., & Koopmans, L. V. E. 2002, MNRAS, 337, L6 [NASA ADS] [CrossRef] [Google Scholar]
  163. Treu, T., & Koopmans, L. V. E. 2004, ApJ, 611, 739 [NASA ADS] [CrossRef] [Google Scholar]
  164. Treu, T., Koopmans, L. V., Bolton, A. S., Burles, S., & Moustakas, L. A. 2006, ApJ, 640, 662 [NASA ADS] [CrossRef] [Google Scholar]
  165. Treu, T., Gavazzi, R., Gorecki, A., et al. 2009, ApJ, 690, 670 [NASA ADS] [CrossRef] [Google Scholar]
  166. Unruh, S., Schneider, P., & Sluse, D. 2017, A&A, 601, A77 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  167. Valdes, F., Gupta, R., Rose, J. A., Singh, H. P., & Bell, D. J. 2004, ApJS, 152, 251 [NASA ADS] [CrossRef] [Google Scholar]
  168. Vanderriest, C., Schneider, J., Herpe, G., et al. 1989, A&A, 215, 1 [NASA ADS] [Google Scholar]
  169. Verolme, E. K., & de Zeeuw, P. T. 2002, MNRAS, 331, 959 [NASA ADS] [CrossRef] [Google Scholar]
  170. Vuissoz, C., Courbin, F., Sluse, D., et al. 2008, A&A, 488, 481 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  171. Weinberg, D. H., Mortonson, M. J., Eisenstein, D. J., et al. 2013, Phys. Rep., 530, 87 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  172. Wertz, O., Orthen, B., & Schneider, P. 2018, A&A, 617, A140 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  173. Weymann, R. J., Latham, D., Angel, J. R. P., et al. 1980, Nature, 285, 641 [NASA ADS] [CrossRef] [Google Scholar]
  174. Wisotzki, L., Schechter, P. L., Bradt, H. V., Heinmüller, J., & Reimers, D. 2002, A&A, 395, 17 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  175. Wong, K. C., Suyu, S. H., Auger, M. W., et al. 2017, MNRAS, 465, 4895 [NASA ADS] [CrossRef] [Google Scholar]
  176. Wong, K. C., Suyu, S. H., Chen, G. C.-F., et al. 2020, MNRAS, 498, 1420 [CrossRef] [Google Scholar]
  177. Wucknitz, O. 2008, MNRAS, 386, 230 [CrossRef] [Google Scholar]
  178. Xu, D., Sluse, D., Schneider, P., et al. 2016, MNRAS, 456, 739 [NASA ADS] [CrossRef] [Google Scholar]
  179. Yang, T., Birrer, S., & Hu, B. 2020, MNRAS, 497, L56 [CrossRef] [Google Scholar]
  180. Yıldırım, A., Suyu, S. H., & Halkola, A. 2020, MNRAS, 493, 4783 [CrossRef] [Google Scholar]
  181. van Albada, T. S. 1982, MNRAS, 201, 939 [NASA ADS] [CrossRef] [Google Scholar]
  182. van de Sande, J., Lagos, C. D. P., Welker, C., et al. 2019, MNRAS, 484, 869 [NASA ADS] [CrossRef] [Google Scholar]
  183. van der Marel, R. P. 1994, MNRAS, 270, 271 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Internal MST + PEMD

Figure A.1 shows different approximate MST’s with a core radius of 10 arcsec on top of a power-law profile (see also Blum et al. 2020). Figure A.2 shows the mock lens used in Sect. 2.6.1 to perform the imaging modeling inference on the lens model parameters, including the cored component resembling the MST.

thumbnail Fig. A.1.

Illustration of the power-law profile (Eq. (39)) in three dimensions (left panel) and in projection (right panel) under an approximate MST with a cored mass component (Eq. (38)). The transforms presented here were indistinguishable by the mock imaging data of Fig. A.2 (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER

thumbnail Fig. A.2.

Mock HST image with a power-law mass profile for which we perform the inference on the detectability of an approximate MST (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER

Appendix B: Mass-anisotropy degeneracy

Figure B.1 shows the predicted projected velocity dispersions (Eq. (16)) in radial bins form the center for PEMD profiles with different logarithmic mass-profile slopes and half-light radii. We chose a fiducial seeing of FWHM = 10. Alternatively, we display the results assuming a constant anisotropy βani(r) = const in Fig. B.2. In Fig. B.3 we plot, without seeing and under fixed anisotropy model, the predicted radial change in the velocity dispersion for different core masses, λc, and core radii, Rc.

thumbnail Fig. B.1.

Radial dependence on the projected velocity dispersion measurement for an Osipkov–Merritt anisotropy profile (Eq. (51)). Top to bottom: increase in the half light radius of the deflector. Left to right: change in the mass profile slope (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/anisotropy_ifu.ipynb).

Open with DEXTER

thumbnail Fig. B.2.

Radial dependence on the projected velocity dispersion measurement for a constant anisotropy βani. Top to bottom: increase in the half light radius of the deflector. Left to right: change in the mass profile slope (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/anisotropy_ifu.ipynb).

Open with DEXTER

thumbnail Fig. B.3.

Radial dependence on the projected velocity dispersion measurement for different cored components (38) on top of a PEMD profile approximating a pure MST, with normalization λc and core radii, Rc. The projected radius from the center of the galaxy is extended to 5 arcsec to visibly see the impact on the kinematic of larger cored components (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/anisotropy_ifu.ipynb).

Open with DEXTER

Appendix C: Likelihood calculation

In this section, we provide the specifics of the likelihood calculation for individual lenses and how we efficiently evaluate the likelihood in the hierarchical context. This includes the imaging likelihood (Appendix C.1), time-delay likelihood (Appendix C.2) and velocity dispersion likelihood (Appendix C.3). Appendix C.4 describes our formalism to track covariances and the marginalization as implemented in HIERARC.

C.1. Imaging likelihood

The likelihood and the lens model inference is not prominently featured in this work, as we are making use of products being derived by our collaboration presented in other work. Nevertheless, the high resolution imaging data and lens model inferences on the likelihood level are essential parts of the analysis.

Given a lens model with parameters ξmass and surface brightness model with parameters ξligth, a model of the imaging data can be constructed, dmodel. The likelihood is computed at the individual pixel level accounting for the noise properties from background and other noise properties, such as read-out, as well as the Poisson contribution from the sources. The imaging likelihood is given by

(C.1)

where k is the number of pixels used in the likelihood and Σpixel is the error covariance matrix. Current analyses assume uncorrelated noise properties in the individual pixels and the covariance matrix becomes diagonal. The model of the surface brightness of the lensed galaxy requires high model flexibility. The surface brightness components can be captured with linear components and solved for and marginalized over analytically. TDCOSMO uses pixelized grids as well as smooth basis sets (see e.g., Suyu et al. 2006; Birrer et al. 2015, for the current methods in use).

C.2. Time-delay likelihood

The likelihood of the time delay data 𝒟td given a model prediction is

(C.2)

with Δtdata is the data vector of relative time delays, is the measurement covariance between the relative delays and

(C.3)

is the model predicted time-delay vector (Eq. (5)) with ΔϕFermat is the relative Fermat potential vector (Eq. (6)). Effectively, the time-delay distance posterior transform according to Eq. (26) under an MST.

C.3. Velocity dispersion likelihood

The model prediction of the velocity dispersion transforms under MST according to Eq. (25) and cosmological distance ratio relevant for the kinematics is Ds/Dds and scales according to Eq. (17). We can write the likelihood of the spectroscopic data, 𝒟spec, given a model as

(C.4)

where is a vector of velocity dispersion measurements, is the measurement error covariance between the measurements (including, for example, stellar template fitting, calibration systematics etc.) and

(C.5)

is the model prediction. The impact of the anisotropy distribution depends on the specific lens and light configuration. We can compute numerically the change in the model predicted dimensionless velocity dispersion component for each individual aperture 𝒜j, J𝒜j(ξmass, ξlight, βani)

(C.6)

C.4. Marginalization and covariances

The marginalization over ξmass and ξlight (Eq. (53)) affects the relative Fermat potential ΔϕFermat in the time-delay likelihood (Eq. (C.3)) and the dimensionless factors (Eqs. (C.5) and (C.6)). We can compute the marginalized likelihood over ξmass and ξlight under the assumption that the posteriors in ξmass and ξlight transform to covariant Gaussian distributions in ΔϕFermat and as a model addition to the error covariances, such that

(C.7)

The model covariance matrix for the time delays can be expressed as

(C.8)

the covariance matrix on the kinematics as

(C.9)

and the cross-covariance between the kinematics and the time delays as

(C.10)

In this form, the model covariances are explicitly dependent on the anisotropy model, the MST and the cosmology.

The covariance between the kinematics and the time delays, , above in Eq. (C.10) is primarily impacted by the average density slope parameter γ of the mass model. γ affects both the kinematics and the Fermat potential and uncertainty in γ can lead to covariances. However, if the density slope parameter is well constrained by imaging data (modulo explicit MST), the covariance in Eq. (C.10) becomes subdominant relative to the uncertainty in the measurement of the kinematics.

When setting , we can separate the inference of DΔt/λ from the kinematics likelihood and can work directly on the DΔt/λ posteriors from the inference from the image data, Dimage, and the time-delay measurement, Dtd,

(C.11)

This allows us to use individually sampled angular diameter distance posteriors (expression (40)) without sampling an additional MST and then transform them in post-processing. This is applicable for both, external convergence and internal MST and we effectively evaluate the likelihood on the one-dimensional posterior density in DΔt/λ.

In the same way as for the time-delay likelihood, we can perform the marginalization of the kinematics likelihood over the imaging data constraints

(C.12)

Appendix D: TDLMC inference with more general anisotropy models

thumbnail Fig. D.1.

TDLMC Rung3 inference with fixed Ωm to the correct value and a generalized Osipkov–Merritt anisotropy profile (Eq. (D.1)). Blue contours indicate the inference with a uniform prior in aani while the red contours indicate the inference with uniform priors in log(aani). The thin vertical line indicates the ground truth H0 value in the challenge (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/TDLMC/TDLMC_rung3_inference.ipynb).

Open with DEXTER

thumbnail Fig. D.3.

TDLMC Rung3 inference on the profile and anisotropy parameter when assuming the correct cosmology for a generalized Osipkov–Merritt anisotropy profile (Eq. (D.1)) (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/TDLMC/TDLMC_rung3_inference.ipynb).

Open with DEXTER

In this work, we presented inferences based on the anisotropy parameterization by Osipkov (1979), Merritt (1985) (Eq. (51)). In this appendix we perform the inference on the TDLMC with a more general anisotropy parameterization. Agnello et al. (2014a) introduced a generalization of the Osipkov–Merritt profile with an asymptotic anisotropy value, β, different than radial

(D.1)

We perform the identical analysis as presented in Sect. 4 except for the addition of one free parameter, β. Table D.1 presents the parameters and priors used in the hierarchical analysis on the TDLMC data set. Figure D.1 shows the results of this inference for the two different priors in aani. The additional degree of freedom in the anisotropy is not constrained by the mock data and leads to a prior-volume effect. The constraining power on the mass profile relies on the mean anisotropy in the orbits within the aperture of the measurement, and not particularly on the parameterization of the radial dependence (see also e.g., Agnello et al. 2014b). It is more challenging to find uninformative priors in higher dimension. As we found an uninformative prior in a simpler parameterization that leads to a consistent result on the TDLMC data set, we do not explore more degrees of freedom in the anisotropy parameterization in this work.

Table D.1.

Summary of the model parameters sampled in the hierarchical inference on TDLMC Rung3 with the anisotropy model of Eq. (D.1).

On the mock data with known input cosmology, we can also reverse the problem and ask which anisotropy parameter configurations result in statistically consistent cosmologies. To do so, we fix the cosmology to the input values and only perform the inference on the anisotropy parameters. Figure D.2 presents the results for the Osipkov–Merritt model of Sect. 4 and Fig. D.3 presents the results for the generalized Osipkov–Merritt profile of this appendix. The posterior on the anisotropy parameter can be interpreted as an informative prior on the anisotropy model parameters from the hydrodynamical simulations of the TDLMC. We do not make use of such a prior in this work but note the consistent inference of the anisotropy parameters for the TDCOSMO+SLACS analysis with this exercise performed on the TDLMC.

Appendix E: SLACS sample details

In this appendix we provide the detailed numerical numbers used in this analysis for the SLACS lenses. Table E.1 lists the data derived from external works that are used in our analysis for the 33 lenses of the SLACS sample. Redshifts are from SDSS presented by Auger et al. (2009), Einstein radii from Auger et al. (2009) and Shajib et al. (2020b) (where available), half-light radii, reff, from Auger et al. (2009), power-law slopes from Shajib et al. (2020b) (where available) and velocity dispersions are based on Bolton et al. (2008) and Shu et al. (2015). Local environment statistics ζ1/r and external shear κext are derived in this work (see Sects. 6.3 and 6.4).

Table E.1.

Summary of the parameters being used of the individual 33 SLACS lenses selected in Sect. 6 to infer mass profile constraints in combination of imaging and kinematics.

All Tables

Table 1.

Summary of the model parameters sampled in the hierarchical inference on TDLMC Rung3 in Sect. 4.

Table 2.

Overview of the TDCOSMO sample posterior products used in this work.

Table 3.

Summary of the model parameters sampled in the hierarchical inference on the TDCOSMO sample in Sect. 5 and posteriors presented in Fig. 7.

Table 4.

Summary of the model parameters sampled in the hierarchical inference on the SLACS sample of Sect. 6.

Table 5.

Summary of the model parameters sampled in the hierarchical inference on the TDCOSMO+SLACS sample.

Table 6.

Marginalized posteriors of our hierarchical Bayesian cosmography inference based on the priors and parameterization specified in Table 5 for a flat ΛCDM cosmology.

Table D.1.

Summary of the model parameters sampled in the hierarchical inference on TDLMC Rung3 with the anisotropy model of Eq. (D.1).

Table E.1.

Summary of the parameters being used of the individual 33 SLACS lenses selected in Sect. 6 to infer mass profile constraints in combination of imaging and kinematics.

All Figures

thumbnail Fig. 1.

Illustration of a composite profile consisting of a stellar component (Hernquist profile, dotted lines) and a dark matter component (NFW + cored component (Eq. (38)), dashed lines) which transform according to an approximate MST (joint as solid lines). The stellar component gets rescaled by the MST while the cored component transforms the dark matter component. Left: profile components in three dimensions. Right: profile components in projection. The transforms presented here cannot be distinguished by imaging data alone and require i.e., stellar kinematics constraints (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_composite_cored.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 2.

Illustration of the constraining power of imaging data on a cored mass component (Eq. (35)). Shown are the parameter inference of the power-law profile mock quadruply lensed quasar of Fig. A.2 when including a marginalization of an additional cored power law profile (Eq. (38)). Orange lines indicate the input truth of the model without a cored component. λc is the scaled core model parameter (Eq. (35)) resembling the pure MST for large core radii (λc ≈ λint) (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 3.

Constraints on an approximate internal MST transform with a cored component, λc, of an NFW profile as a function of core radius. In gray are the 1-σ exclusion limits that imaging data can provide. In orange is the region where the total mass of the core within a three-dimensional radius exceeds the mass of the NFW profile in the same sphere. In blue is the region where the transformed profile results in negative convergence at the core radius. The white region is effectively allowed by the imaging data and simple plausibility considerations and where we can use the mathematical MST as an approximation (λc ≈ λint). The halo mass, concentration and the redshift configuration is displayed in the lower left box (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 4.

Comparison of the actual predicted kinematics from the modeling of the physical three-dimensional mass distribution κλint (Eq. (35)) for varying core sizes (solid) and the analytic relation of a perfect MST (Eq. (25), dashed) for the mock lens presented in Fig. A.2. Lower panel: fractional differences between the exact prediction and a perfect MST calculation. The MST prediction matches to < 1% in the considered range. Minor numerical noise is present at the subpercent level (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 5.

Mock data from the TDLMC Rung3 inference with the parameters and prior specified in Table 1. Orange contours indicate the inference with a uniform prior in aani while the purple contours indicate the inference with a uniform priors in log(aani). The thin vertical line indicates the ground truth H0 value in the challenge (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/TDLMC/TDLMC_rung3_inference.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 7.

Hierarchical analysis of the TDCOSMO-only sample when constraining the MST with kinematic information. Parameter and priors are specified in Table 3. Orange contours correspond to the inference with uniform prior on Ωm, 𝒰([0.05, 0.5]), while the purple contours correspond to the prior based on the Pantheon sample with 𝒩(μ = 0.298, σ = 0.022) (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/TDCOSMO_sample/tdcosmo_sample.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 8.

Sample selection of the SLACS lenses being added to the analysis and comparison with the TDCOSMO data set. The comparisons are in lens redshift, zlens, source redshift, zsource, measured velocity dispersion, σP, half light radius of the deflector, reff, Einstein radius of the deflector, θE, the ratio of half light radius to Einstein radius, reff/θE, and the SIS equivalent velocity dispersion estimated from the Einstein radius and a fiducial cosmology, σSIS. Open dots correspond to lenses included in our selection without quality lens models. Red points correspond to SLACS lenses which have VIMOS IFU data (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/sample_selection.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 9.

External convergence posteriors of the 33 SLACS lenses that pass our morphology and local environment selection cut based on the weighted number counts (gray dashed lines) and the population distribution (black solid line) (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/sample_selection.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 10.

Power-law slope γpl inferences obtained from HST imaging data for the 14 SLACS lenses within our selection cut from Shajib et al. (2020a,b). We derive a population distribution from these lenses to be applied on the subset of lenses without measured γpl from imaging data. The black dashed line indicates the population mean and the blue band the 1-sigma population width based on the 14 individual measurements (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/sample_selection.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 11.

Posterior distribution for the SLACS sample with priors according to Table 4. Orange: inference with the SDSS spectra. Purple: inference with SDSS spectra and VIMOS IFU data set. The posterior of λint, 0 was blinded during the analysis (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/SLACS_sample/SLACS_constraints.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 12.

Posterior distributions of the key parameters for the hierarchical inference. Blue: constraints from the TDCOSMO-only sample. Violet: constraints with the addition of IFU data of nine SLACS lenses to inform the anisotropy prior on the TDCOSMO sample, TDCOSMO+SLACSIFU. Orange: constraints with a sample of 33 additional lenses with imaging and kinematics data (HST imaging + SDSS spectra) from the SLACS sample, TDCOSMO+SLACSSDSS. Purple: joint analysis of TDCOSMO and 33 SLACS lenses with SDSS spectra of which nine have VIMOS IFU data, TDCOSMO+SLACSSDSS + IFU. Priors are according to Table 5. The 68th percentiles of the 1D marginalized posteriors are presented in Table 6. The posteriors in H0 and λint, 0 were held blinded during the analysis (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/joint_inference.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 13.

Illustration of the goodness of the fit of the maximum likelihood model of the joint analysis in describing the TDCOSMO data set. Blue points are the measurements with the diagonal elements of the measurement covariance matrix. Orange points are the model predictions with the diagonal elements of the model covariance uncertainties. Left: comparison of measured time-delay distance from imaging data and time delays compared with the predicted value from the cosmological model, the internal and external MST (and their distributions). Right: comparison of the velocity dispersion measurements and the predicted values. In addition to the MST terms, the uncertainty in the model also includes the uncertainty in the anisotropy distribution aani. For lenses with multiple velocity dispersion measurements, the diagonal terms in the error covariance are illustrated (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/joint_inference.ipynb). (a) Fit to the time-delay distance. (b) Fit of velocity dispersion.

Open with DEXTER
In the text
thumbnail Fig. 14.

Illustration of the goodness of the fit of the maximum likelihood model of the joint analysis in describing the SDSS velocity dispersion measurements of the 34 SLACS lenses in our sample. Blue points are the measurements with the diagonal elements of the measurement covariance matrix. Orange points are the model predictions with the diagonal elements of the model covariance uncertainties. The measurement uncertainties include the uncertainties in the quoted measurements and the additional uncertainty of σσP, sys. The model uncertainties include the lens model uncertainties and the marginalization over the λint and aani distribution (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/joint_inference.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 15.

Illustration of the goodness of the fit of the maximum likelihood model of the joint analysis in describing the VIMOS radially binned IFU velocity dispersion measurements of the nine SLACS lenses with VIMOS data in our sample. Blue points are the measurements with the diagonal elements of the measurement covariance matrix. Orange points are the model predictions with the diagonal elements of the model covariance uncertainties. The measurement uncertainties include the uncertainties in the quoted measurements and the additional uncertainty of σσP, sys. The model uncertainties include the lens model uncertainties and the marginalization over the λint and aani distribution (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/joint_inference.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 16.

Relative difference of the measured vs the predicted velocity dispersion for the SLACS and TDCOSMO sample as a function of different parameters associated with the line-of-sight and the lensing galaxy. In particular, these are the relative inverse distance weighted over density ζ1/r, the ratio of half-light radius to Einstein radius reff/θE, lens redshift zlens, and SIS equivalent velocity dispersion σSIS (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/joint_inference.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 17.

Illustration of the inferred mass profile of the joint TDCOSMO+SLACSSDSS + IFU analysis. A pure power-law with γpl = 2.10 ± 0.05 is shown in orange. In blue is the result of this work of λint = 0.91 ± 0.045 when interpreted as a cored mass component with Rc uniform in [3″, 10″]. Three dimensional density are illustrated on the left and the lensing convergence on the right. The dashed vertical line on the right panels indicates the Einstein radius. Relative difference in respect to the power-law model are presented in the bottom panels (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER
In the text
thumbnail Fig. 18.

Comparison of different blind H0 measurements by the TDCOSMO collaboration, based on different mass profile assumptions and data sets incorporated. All measurements presented on this plot were performed blindly with regard to the inference of H0. The measurement on top is the combined H0LiCOW six lenses constraints presented by Wong et al. (2020), when averaging power-law and composite NFW plus stars (with constant mass-to-light ratio) on a lens-by-lens basis without correlated errors among the lenses. The next two measurements are from Millon et al. (2020) of six TDCOSMO time-delay lenses (five H0LiCOW lenses22 and one STRIDES lens by Shajib et al. 2020a), when performing the inference assuming either a composite NFW plus stars (with constant mass-to-light ratio) or the power-law mass density profile for the galaxy acting as a lens. Lower panel: results from this work. The main difference with respect to previous work is that we have made virtually no assumption on the radial mass density profile of the lens galaxy, and taken into account the covariance between the lenses. The analysis in this work is constrained only by the stellar kinematics and fully accounts for the uncertainty related to the mass sheet transformation (MST). In this framework, we obtain four measurements according to the datasets considered. The TDCOSMO-only inference is based on the same set of seven lenses as those jointly included by Millon et al. (2020) and Wong et al. (2020). The inferred median value is the same, indicating no bias, and the uncertainties, as expected, are larger. The next three measurements rely on external datasets from the SLACS survey, by making the assumption that the lens galaxies in the two surveys are drawn from the same population. The TDCOSMO+SLACSIFU measurements uses, in addition to the TDCOSMO sample, nine lenses from the SLACS sample with IFU observations to inform the anisotropy prior applied on the TDCOSMO lenses. The TDCOSMO+SLACSSDSS measurement comes from the joint analysis of the TDCOSMO sample and 33 SLACS lenses with SDSS spectroscopy. The TDCOSMO+SLACSSDSS + IFU presents the joint analysis of all three data sets, again assuming self-similar distributions of the mass profiles and stellar anisotropy. The TDCOSMO-only and TDCOSMO+SLACSIFU analyses do not rely on self-similar mass profiles of the SLACS and TDCOSMO sample while the TDCOSMO+SLACSSDSS and TDCOSMO+SLACSSDSS + IFU measurements (orange and purple) do. All the measurements shown in this plot are in statistical agreement with each other. See Sect. 8.5 for a discussion and physical interpretation of the results (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/JointAnalysis/tdcosmo_comparison_plot.ipynb).

Open with DEXTER
In the text
thumbnail Fig. A.1.

Illustration of the power-law profile (Eq. (39)) in three dimensions (left panel) and in projection (right panel) under an approximate MST with a cored mass component (Eq. (38)). The transforms presented here were indistinguishable by the mock imaging data of Fig. A.2 (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER
In the text
thumbnail Fig. A.2.

Mock HST image with a power-law mass profile for which we perform the inference on the detectability of an approximate MST (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/MST_pl_cored.ipynb).

Open with DEXTER
In the text
thumbnail Fig. B.1.

Radial dependence on the projected velocity dispersion measurement for an Osipkov–Merritt anisotropy profile (Eq. (51)). Top to bottom: increase in the half light radius of the deflector. Left to right: change in the mass profile slope (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/anisotropy_ifu.ipynb).

Open with DEXTER
In the text
thumbnail Fig. B.2.

Radial dependence on the projected velocity dispersion measurement for a constant anisotropy βani. Top to bottom: increase in the half light radius of the deflector. Left to right: change in the mass profile slope (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/anisotropy_ifu.ipynb).

Open with DEXTER
In the text
thumbnail Fig. B.3.

Radial dependence on the projected velocity dispersion measurement for different cored components (38) on top of a PEMD profile approximating a pure MST, with normalization λc and core radii, Rc. The projected radius from the center of the galaxy is extended to 5 arcsec to visibly see the impact on the kinematic of larger cored components (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/MST_impact/anisotropy_ifu.ipynb).

Open with DEXTER
In the text
thumbnail Fig. D.1.

TDLMC Rung3 inference with fixed Ωm to the correct value and a generalized Osipkov–Merritt anisotropy profile (Eq. (D.1)). Blue contours indicate the inference with a uniform prior in aani while the red contours indicate the inference with uniform priors in log(aani). The thin vertical line indicates the ground truth H0 value in the challenge (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/TDLMC/TDLMC_rung3_inference.ipynb).

Open with DEXTER
In the text
thumbnail Fig. D.3.

TDLMC Rung3 inference on the profile and anisotropy parameter when assuming the correct cosmology for a generalized Osipkov–Merritt anisotropy profile (Eq. (D.1)) (https://github.com/TDCOSMO/hierarchy_analysis_2020_public/blob/6c293af582c398a5c9de60a51cb0c44432a3c598/TDLMC/TDLMC_rung3_inference.ipynb).

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.