Inference of the cold dark matter substructure mass function at z = 0.2 using strong gravitational lenses

Abstract

1 INTRODUCTION

One of the most robust predictions of the cold dark matter (CDM) theory is that dark matter haloes of all scales should be populated by thousands of low-mass substructures (Kauffmann, White & Guiderdoni 1993; Klypin et al. 1999; Diemand et al. 2008; Springel et al. 2008). While an equivalent number of satellites has not yet been detected in observations of the Local Group galaxies (e.g. Kravtsov 2010, and references therein), reconciling the observations with the theoretical expectations could be possible by invoking the existence of a large population of low surface brightness or even dark satellites. Alternatively, solving the problem would require modifying the nature of dark matter (e.g. Nierenberg et al. 2013a).

Strong gravitational lensing, being sensitive to any form of matter, provides a unique opportunity to address this issue. It was initially suggested by Mao & Schneider (1998) that flux-ratio anomalies of multiply imaged lensed quasars can be interpreted as the signature of mass substructure in the lens galaxies at scales smaller than the lensed image separation. Indeed, Metcalf & Madau (2001), Chiba (2002), Dalal & Kochanek (2002), Metcalf & Zhao (2002), Keeton, Gaudi & Petters (2003) and Nierenberg et al. (2014) showed that mass substructure could be responsible for a violation of the magnification relations.

While early comparisons with N-body numerical simulations have indicated an agreement between the predicted and the gravitational lensing inferred amount of substructure (e.g. Bradač et al. 2002), more recently Xu et al. (2009) have shown that the probability of reproducing the observed rate of anomalies with the subhalo population of the Aquarius simulation is of the order of 10⁻³. The inclusion of line-of-sight structures from the Millennium II simulations only increases the chance of flux-ratio anomalies to about 10–30 per cent (Xu et al. 2012). At the same time however, gravitational lens galaxies are more commonly massive early-type galaxies and are hosted therefore by dark matter haloes that are more massive and potentially different (e.g. with a higher ellipticity) than those of the Aquarius simulation, which are primarily Milky-Way-type haloes. In particular, Metcalf & Amara (2012) using more realistic simulations of lens galaxies have found consistent results between the expectations and the observations. A similar conclusion was more recently obtained by Xu et al. (2013).

Schechter & Moore (1993), McKean et al. (2007) and More et al. (2009) have found that in some cases, the flux-ratio anomaly is related to the presence of a single massive luminous satellite rather than a large population of small mass CDM substructures. It has also been suggested that evidence for mass substructure in lens galaxies can be obtained from astrometric measurements alone (Chen et al. 2007) or in combination with flux-ratio anomalies (Fadely & Keeton 2012) and time-delay measurements (Keeton & Moustakas 2009).

Koopmans (2005) and Vegetti & Koopmans (2009a) have developed an alternative technique of direct gravitational imaging of mass substructure based on perturbations to the surface brightness distribution of lensed arcs and Einstein rings. While the gravitational imaging technique is complementary to the above methods, it has the major advantage of being able to directly measure the substructure mass and position, and therefore is less degenerate in the mass models. This feature of the gravitational imaging technique makes it particularly powerful in constraining the substructure mass function (Vegetti & Koopmans 2009b).

The gravitational imaging technique has previously been applied to high resolution Hubble Space Telescope (HST; Vegetti et al. 2010) and Keck adaptive optics imaging (Vegetti et al. 2012). This method has already led to the detection of two substructures at redshifts of 0.222 and 0.881, with respective masses of (3.51 ± 0.15) × 10⁹ and (1.9 ± 0.1) × 10⁸ M_⊙. The substructure mass here quoted is the total 3D mass measured under the assumption of a spherical pseudo-Jaffe mass profile (Section 4.3). As shown in Vegetti et al. (2012) the mass so measured is consistent with the projected mass derived in a model-independent fashion from the pixelized density corrections (Section 4.2). This indicates that the measured lensing mass is not significantly affected by our choice of mass density profile (see also Vegetti & Vogelsberger 2014). Under the assumption that these detections are associated with the host lens galaxy, they imply a substructure projected (within ∼5 kpc) mass fraction of |$\langle f \rangle =0.033^{+0.036}_{-0.018}$| and a slope of the substructure mass function of |$\alpha =1.1^{+0.6}_{-0.4}$| (where dN/dm ∝ m^−α). Within the large errors, these results are consistent with the theoretical expectations at the 95 per cent confidence level (CL; Diemand et al. 2008; Springel et al. 2008; Xu et al. 2013). The mass fraction is inferred by assuming a mass function normalized between 4 × 10⁶ and 4 × 10⁹ M_⊙.

In this paper, we present a statistical analysis of a sample of 11 gravitational lens galaxies with the aim of deriving tighter constraints on the substructure mass function. The layout of the paper is as follows. In Section 2 we provide a short description of the analysed data. In Section 3 we provide an overview of the full analysis from the lens modelling to the derivation of the substructure mass function. In Section 4 we describe the gravitational imaging technique, while in Section 5 we present the main lens modelling results. In Section 6 we derive statistical constraints on the substructure mass function and we finally conclude in Section 7.

Throughout the paper we assume the following cosmology, H₀ = 73 km s⁻¹ Mpc⁻¹, Ω_m = 0.25 and |$\Omega _\Lambda =0.75$|⁠.

2 DATA

The gravitational lenses analysed in this paper are all taken from the Sloan Lens ACS Survey (SLACS; e.g. Bolton et al. 2008a). These were selected from the main (Strauss et al. 2002) and luminous red galaxy (LRG; Eisenstein et al. 2001) samples of the Sloan Digital Sky Survey (SDSS; York et al. 2000). The SLACS gravitational lenses are generally massive early-type galaxies, typical of the parent samples from which they are selected (e.g. Bolton et al. 2006).

The data used in the current analysis were obtained using the wide-field channel (WFC) of the Advanced Camera for Surveys (ACS) aboard the HST. A single orbit was used for each unique filter, split into four dithered exposures. The data were rectified and co-added on to a uniform image grid using custom software written in the idl language. The surface brightness profiles of the lensing foreground galaxies were modelled and subtracted using the radial B-spline formalism described by Bolton et al. (2006, 2008a). Initial models for the lens galaxy surface brightness distribution were optimized non-linearly to determine the centroid, axis ratio and position angle of the lens galaxies. These parameters were then used to define an elliptical coordinate system for the final linear model fitting, including sufficient multipole orders to provide model-subtracted images free of appreciable systematic residuals. Masks for the lensed image features were generated manually and applied during the B-spline model fitting. In some cases, these masks were updated based on initial B-spline subtracted residual images, and the fitting procedures were repeated. Vegetti et al. (2010) have explicitly tested the effect of the galaxy subtraction procedure on the significance of substructure detections using a Sérsic profile rather than a B-spline surface brightness profile and found no significant influence. Moreover, the B-spline subtraction was done independently for each filter and we find in this paper that the sensitivity function to substructures (see Section 6.1) is consistent across filters. A similar test with an equivalent output was done by Vegetti et al. (2012). We are confident, therefore, that the galaxy subtraction procedure is not a major source of systematic error on our inference of the substructure mass function. This can be understood because the effect of a substructure on the brightness distribution of a lensed image occurs on spatial scales over which the galaxy surface brightness model varies very little (i.e. there is very little mixing of spatial scales and early-type galaxies show very little structure on the scales where substructure causes anomalies).

In order to increase our sensitivity to lower mass substructures, the subsample of lenses investigated in this paper was selected purely on the basis of the lensed images signal-to-noise ratio (S/N) in the I and V band, whenever ACS images are available. Specifically, all the images in the full 11 lens sample have at least 200 pixels with a S/N ≥ 2 and at least 50 pixels with S/N ≥ 3. This ensures that there is a region on the lensed images large enough for substructures to be detected when present with a mass above the detection threshold. This selection does not bias in favour or against substructure in the lens, since substructure in the probed mass range has a negligible effect on the source magnification. Since no substructure detection was obtained for the sample here considered, the lack of V-band data for some of the systems does not represent a major issue for the interpretation of the lens modelling. When constraining the substructure mass function, we also consider our previous detection in the gravitational lens system SDSS J0946+1006 that was similarly analysed by Vegetti et al. (2010). In Table 1 we list all of the lens systems in our sample with the lens and source redshifts, the filter used for the imaging and other relevant properties. The gravitational lens system B1938+666, previously analysed by Vegetti et al. (2012), is not included in the current sample as this system is not part of the SLACS survey but of the Cosmic Lens All-sky Survey (CLASS), and follows therefore a different selection function (source based rather than lensed based). Moreover, since the amount of substructure in lens galaxies can significantly depend on the properties of the host galaxy (e.g. Xu et al. 2013), by limiting our sample to the SLACS lenses, we can avoid possible issues related to trends of the substructure mass fraction with the host redshift and mass.

3 ANALYSIS OVERVIEW

We summarize here all of the steps involved in our analysis and we refer to the following sections and to Koopmans (2005) and Vegetti & Koopmans (2009a,b) for details on the lens modelling technique and the statistical inference on the substructure mass function:

• all of the lens galaxies in the sample are initially modelled under the assumption of a smooth mass distribution with an elliptical power-law density profile and with the inclusion of external shear (Section 4.1);

• for each lens, the smooth model is improved with linear and localized potential corrections, which are examined for the presence of mass substructures (Section 4.2);

• all of the lens systems that have shown evidence for non-negligible (at the 10σ level) potential corrections are remodelled with the same smooth mass model family of step (i), plus one or more mass substructures with a given analytical mass density profile (Section 4.3);

• for each lens we calculate the 10σ mass detection threshold at each pixel in the lens plane and derive the mass sensitivity as a function of image-plane position (Section 6.1, Fig. 1);

• we determine the joint likelihood for the set of detections and non-detections and then transform it, in a Bayesian way, into posterior probability density functions (PDFs) for the CDM substructure mass fraction and mass function slope (Section 6.3, Fig. 2).

4 LENS MODELLING TECHNIQUE

The gravitational lens modelling is carried out in a three step procedure using the Bayesian grid-based technique of Vegetti & Koopmans (2009a). Below, we describe each step of this technique in more detail.

4.1 Smooth models

The free parameters of this mass model are the mass density normalization κ₀, its position angle θ, the axis ratio q, the coordinate centre x and y, the logarithmic slope γ, the external shear strength Γ and its position angle Γ_θ. The mass distribution is assumed to have a negligible core radius (i.e. r_c ≡ 0 arcsec).

4.1.1 Model optimization and inference

4.2 Lens-potential corrections

4.3 Clumpy models

For a pseudo-Jaffe substructure in a spherically symmetric isothermal galaxy, with mass normalization κ₀, the truncation radius in arcseconds is therefore |$r_{\rm t} = r\sqrt{\pi \kappa _{0,{{\rm sub}}}/(2\beta \kappa _{0}})$|⁠, which reduces to |$r_{\rm t}=\sqrt{\pi \kappa _{0,{{\rm sub}}}\kappa _{0}/(2\beta )}$| if the substructure is located at the Einstein radius of the lens. In addition to the free parameters defined in Section 4.1, the free parameters of the clumpy model now also include the substructure position and mass. In order to allow for a proper comparison, we use the same source regularization form and the same number of pixels n_src characterizing the source grid that we used for the smooth model analysis. The posterior probability density distribution |$P(\lambda _{{\rm s}},\boldsymbol \eta \,|\,\boldsymbol d,{\bf M}_{\rm c},{\bf R}_{{\rm s}})$|⁠, is now maximized for |$\boldsymbol \eta =\lbrace \kappa _0,\theta ,q,x,y,\gamma ,\Gamma ,\Gamma _{\theta },M_{{\rm sub}},x_{{\rm sub}},y_{{\rm sub}}\rbrace$| and λ_s, following the same three-step procedure outlined above. This analytic clumpy model can also be used to assess the significance of the substructure detection via the total marginalized Bayesian evidence, as described in the following section.

4.4 Model comparison

5 LENS MODELLING RESULTS

In this section, we provide more detailed comments on the individual systems and give comparisons with previous results where appropriate.

5.1 Smooth models

The MAP models under the assumption of a smooth mass distribution are listed in Table 2 for each system and plotted in Figs 3–12. Models for SDSS J0946+1006 can be found in Vegetti et al. (2010). Smooth models for all of the SLACS lenses, under the assumption of a singular isothermal elliptical (SIE) mass distribution, were derived by Bolton et al. (2008a) (see also Auger et al. 2009), while for a subset of the lenses considered here, a joint lensing and 2D kinematics analysis can be found in the paper by Barnabè et al. (2011). A measure of the slope of the mass density profile γ was also derived by Koopmans et al. (2009) and Auger et al. (2010), using a combination of lensing and the stellar velocity dispersion for each galaxy. We note that the comparison between our results and those by Koopmans et al. (2009), Bolton et al. (2008a) and Barnabè et al. (2011) is not a straightforward one and that agreement should not be necessarily expected. For example, the mass density slope derived in this paper is a measure of the local slope near the Einstein radius and could be therefore different from the average mass density slope measured by Koopmans et al. (2009) and Barnabè et al. (2011), which are mostly set to fit the dynamical data (see discussion by Sonnenfeld et al. 2013). We do find that the Einstein radius derived in this paper is for every lens consistent with those determined by Bolton et al. (2008a), although our lenses are often not (but close to) isothermal.

The lens system SDSS J1430+4105 was modelled by Eichner, Seitz & Bauer (2012) by making use of different assumptions on the mass model for the main lens and its environment, and providing consistent results within the error among several models. Here, we can only directly compare with the power-law plus external shear model from Eichner et al. (2012) derived by modelling the surface brightness distribution of the lensed images with the lensing code lensview (Wayth & Webster 2006). We find that the Einstein radius derived in this paper is  ∼1 per cent smaller than their Einstein radius, while our mass distribution has an axis ratio which is  ∼15 per cent smaller and a mass density slope which is  ∼10 per cent steeper. Finally, the external shear strength is  ∼8 per cent larger than the value recovered by Eichner et al. (2012) and its position angle is not consistent with their model, but it is consistent with the location of the brightest cluster galaxy of the lens environment (Koester et al. 2007). These differences are consistent with the scatter on the lens parameters derived by Eichner et al. (2012) under different model assumptions. We believe that they are related to a different choice of the surface mass density normalization, which in our case explicitly encodes the mass distribution axis ratio and logarithmic slope (equation 1). Also, the power-law plus external shear model by Eichner et al. (2012) does not seem to recover the external shear correctly.

In the case of SDSS J0956+5100 a luminous companion was identified in the HST images. This luminous satellite was included in the lens model as an additional mass component with a truncated pseudo-Jaffe (see equation 8) mass density profile.

5.2 Potential corrections and clumpy models

We claim the presence of a substructure detected at a significant statistical level when the following conditions are simultaneously satisfied:

• a positive convergence correction that improves the image residuals is found independently from the potential regularization, number of source pixels, PSF rotations and galaxy subtraction procedure;

• a clumpy model is preferred over a smooth model with a Bayes factor |$\Delta \log {{\cal E}}=\log {{\cal E}_{{\rm smooth}}}-\log {{\cal E}_{{\rm clumpy}}}\ge -50$| (to first order equivalent to a 10σ detection, under the assumption of Gaussian noise);

• the mass and the position of the substructure obtained via the nested sampling analysis are consistent with those independently obtained by the potential corrections and the MAP parametric clumpy model;

• the results are consistent among the different HST filters, where available.

We find that none of the systems considered in this paper satisfies the above conditions, implying that no significant detection of a mass substructure was made above the mass threshold set by the S/N and resolution of the HST images. Importantly, we notice that some lenses satisfy the condition (ii) without satisfying conditions (iii) and (iv). We believe that this is due partly to degeneracies in the lens model parameters, a relatively poor fitting of the smooth model and the fact that only a clumpy model was considered as an alternative to the smooth one. Generalizing, we stress that is important to verify the uniqueness of the substructure detection by taking fully into account the degeneracies with the macromodel and the source model and by testing other models than clumpy ones as an alternative to smooth models, both with gravitational imaging and with flux ratio anomalies.

5.2.1 SDSS J0956+5100

For the case of SDSS J0956+5100, since the luminous satellite was explicitly included as part of the smooth model, this implies that no substructure other than the luminous satellite was found. By repeating the analysis of this lens without the inclusion of the luminous satellite, we find that this satellite is also detected via the gravitational imaging technique in both I and V bands at the  ∼37σ and  ∼20σ level, respectively. Including our previous detection in SDSS J0946+1006 (Vegetti et al. 2010), we detected in total, two satellites in 11 gravitational lens systems. In both cases the satellites have a projected distance of about ≤5 kpc from the centre of the host galaxy. Although we do not know its redshift, the colours of this luminous clump are consistent with those of the lens galaxy, indicating that it is most likely a satellite. By fitting the surface brightness of this satellite with a Gaussian profile, we find that the total magnitude is 22.70 and 22.93 in I- and V band, respectively, which is about 6 mag fainter than the host galaxy. From Auger et al. (2009) we derive a stellar mass of M_* = 1.4 × 10⁹ M_⊙ for a Chabrier initial mass function and M_* = 2.5 × 10⁹ M_⊙ for a Salpeter one. Under the assumption of a Pseudo-Jaffe total mass density profile, we infer a total mass of M_sub = (1.69 ± 0.01) × 10¹⁰ and (1.89 ± 0.05) × 10¹⁰ M_⊙ in I- and V band, respectively. This yields a stellar- to total-mass ratio of  ∼10, which is not unexpected considering that the satellite has probably been significantly tidally stripped. Finally, we note that the MAP position is located in both bands offset from the centre of the light distribution by 3 pixels (0.15 arcsec) in the x direction and 8 pixels (0.4 arcsec) in the y direction.

5.2.2 SDSS J0946+1006

Vegetti et al. (2010) reported the detection of a dark-matter-dominated substructure in the gravitational lens system SDSS J0946+1006. When modelled as a pseudo-Jaffe profile, the total mass of this substructure was inferred to be M_sub = (3.51 ± 0.15) × 10⁹ M_⊙. The substructure was located very close to the Einstein radius but in a position where the surface brightness is not dominated by the lensed images. This allowed a tight constraint to be set on the upper limit of the substructure luminosity, hence its mass-to-light ratio lower limit is (M/L)_{V, ⊙} ≥ 120 M_⊙/L_{V, ⊙} (3σ) inside a sphere of radius 0.3 kpc.

5.2.3 Systematic error on the substructure masses

The major source of systematic error on the substructure total mass is related to the de-projection of the substructure position. In particular, we have assumed that the detected mass substructure and satellites are at the same redshift as the lens galaxy, and when using a pseudo-Jaffe profile, that the observed projected position is the true physical position. As discussed previously, in the case of SDSS J0956+5100, the luminous satellite has colours consistent with the lens galaxy redshift. The mass clump in SDSS J0946+1006 is a high mass-to-light ratio object for which direct light was not detected and colours could not be measured. In principle, we can estimate the probability that this clump is at the same redshift of the lens galaxy by making use of numerical simulations, in practice however, current N-body simulations of massive galaxies do not reach the necessary resolution, especially at distances so close to the centre of the host halo. Moreover, one would need to take into account for the different lensing effects of clumps at different redshifts. This can significantly depend on the density profile that is assumed and, in the case of extended lensed images, it cannot be computed by linearly rescaling the clump mass. These are all important issues that we will extensively address in a future paper.

6 SUBSTRUCTURE MASS FUNCTION

In this section, we use the results of this paper, in combination with the mass substructure detected by Vegetti et al. (2010), to constrain the substructure mass function slope and normalization. For this purpose, we assume that the detected mass density corrections are due to mass clumps associated with the lens galaxy, and not with line-of-sight contaminations. In the latter case, the inferred substructure fractions in the galaxy would decrease. Although only one out of 11 lens galaxies shows evidence for mass substructure, the remaining 10 non-detections also provide valuable information on the mass function, once the sensitivity to substructure for each lens galaxy has been quantified. For each lens in the considered sample, we quantify the substructure sensitivity in the following subsection.

6.1 Substructure sensitivity function

Vegetti & Koopmans (2009b) and Vegetti et al. (2010, 2012) assumed that the sensitivity to substructure was constant (i.e. it only depended on the substructure mass but not its position) over an annulus of a few arcseconds width around the Einstein radius and that substructures were randomly distributed within this annulus with a uniform probability density distribution. In practice, this means assuming a constant value for the lowest detectable substructure mass within the considered region. In this paper, we quantify the sensitivity to substructure for each lens system in the considered sample as a function of the substructure position and mass. This is particularly important to properly take into account for the large number of non-detections.

As shown by Koopmans (2005) and Vegetti & Koopmans (2009a), the effect of a potential correction on the surface brightness distribution of the lensed images is given by |$\delta I(\boldsymbol {x})=\nabla S(\boldsymbol {y})\cdot \nabla \delta \psi (\boldsymbol {x})$|⁠, where |$\delta I(\boldsymbol {x})$| is the difference in the surface brightness distribution of the lensed images between a smooth lens and a model that contains substructure, |$\nabla S(\boldsymbol {y})$| is the gradient of the source surface brightness distribution and |$\nabla \delta \psi (\boldsymbol {x})$| is the gradient of the potential correction |$\delta \psi (\boldsymbol {x})$|⁠. For a given substructure, therefore, its effect is different on different parts of the lensed images and larger on those parts which are more structured (see Rau, Vegetti & White 2013). This, in combination with the fact that the effect of a substructure can be partly re-absorbed with a change of the source surface brightness structure and/or with a change in the lens macromodel, implies that the most rigorous way to quantify the sensitivity to a given substructure is via the probability density of a specific substructure model given the data, marginalized over all possible macromodel parameters and source regularization constant (i.e. via equation 10). This would require one to compute the Bayes factor and marginalize over the lens parameters and source parameters for a large number of clumpy models, each defined by a different combination of substructure mass and position (ideally one for every image pixel). However, this approach is computationally prohibitive, and the MultiNest technique does not sample the posterior probability density distribution densely enough to allow us to derive the sensitivity function directly from the analysis of Section 4.4. At the same time, we find that most of the re-absorption is done by changes in the source surface brightness distribution rather than changes in the macromodel parameters. We can therefore reduce the multidimensional integral of the Bayesian evidence to a one-dimensional integral over the source regularization constant, while keeping the lens (and substructure) parameters fixed. Because this integral implicitly involves also an integral over the source surface brightness distribution via equation (4), it properly takes into account for the degeneracy between the substructure properties and the structure of the source. We, therefore, properly quantify the sensitivity function by creating 100 N_pix mocked lens realizations, each of which is based on the most probable smooth lens model and most probable source model as derived in Section 4.1, and has a substructure of given mass between 4 × 10⁶ and 4 × 10¹⁰ M_⊙ (uniformly sampled in log-scale) and located at one of the image pixels. This assumes that we are not sensitive to changes in the substructure position within a single pixel. For each of these realizations, we then compute the marginalized Bayesian evidence under the assumption of a smooth and a clumpy lens model and define the mass detection threshold at each position on the lens plane as the lowest mass for which the Bayes factor is |$\Delta \log {{\cal E}}\ge -50$| (see Fig. 1), as this would correspond to a 10σ detection, under the assumption of Gaussian noise. We find that given the data quality of the sample in this paper, the lowest detectable mass typically varies on average between 0.04 × 10¹⁰ and 0.14 × 10¹⁰ M_⊙ depending on the lens. In some cases, for specific substructure positions the detection threshold, however, can be as low as 0.01 × 10¹⁰ M_⊙.

Figure 1.

Lowest detectable mass (in log [M_low/10¹⁰ M_⊙]) as function of the substructure position for the lens system SDSS J0252+0039.

Open in new tabDownload slide

Finally, we note that re-absorption effects were properly taken into account in the detection of the mass substructure in the lens system SDSS J0946+1006, as its significance was properly calculated via the integral in equation (10). We refer to Vegetti et al. (2010) for a more detailed discussion on this particular lens systems.

6.2 Likelihood of substructure detections

6.2.1 Considered region

Given a substructure mass, its effect on the surface brightness of the lensed images is a function of the substructure distance from the lensed images, the data angular resolution and S/N, and the level of structure of the source surface brightness distribution. Thus, for a given data quality (S/N and angular resolution) the gravitational effect of a mass substructure is maximized when located right on top of the lensed images of a highly structured source. We therefore restrict our measurement of the mass function to a small region where the S/N of the data surface brightness distribution is larger than 3. This approximately corresponds to a small area around the critical curve that can vary between  ∼4 and  ∼200 kpc² (i.e. between 2 and 13 arcsec²) around the lensed images. In this way, we are only considering substructure positions where the effect of the substructure is large and unlikely to be re-absorbed by changes in the lens macromodel. We note that this region does not include the luminous satellite of SDSS J0956+5100, after having subtracted the light from the satellite itself, hence this luminous satellite is not included in the Bayesian inference of the substructure mass function.

Our mass function constraints are therefore solely derived by making use of the substructure detected in SDSS J0946+1006 and the 10 other non-detections. We note that changes in the window in which we allow a detection does not bias the inference since we properly account for the mass detection threshold and dark matter mass inside the window. Shrinking the window, however, does lower the detection probability, increasing the error bars. But, since the expectation value of the number of substructures (per unit area) strongly peaks around the lensed images, the precise window choice (if not too small) has little impact and choosing it relatively narrow significantly reduces the computation effort in deriving the mass thresholds as a function of position in the image plane.

6.2.2 Lens galaxy dark matter mass

Vegetti & Koopmans (2009b) and Vegetti et al. (2010, 2012) calculated the host lens galaxy dark matter mass assuming an isothermal total mass density distribution, a Jaffe profile for the stellar mass distribution and the scaling relations derived for the SLACS lenses by Bolton et al. (2008b). Here, we make use instead of the dark matter fraction derived by Auger et al. (2010). In particular, we evaluate the de Vaucoleurs surface brightness profiles from Bolton et al. (2008b) and Auger et al. (2009) at the Einstein radius and then multiply by the mass-to-light ratio implied by Auger et al. (2009) under the assumption of a Salpeter initial mass function (consistent with results from lensing and dynamics, stellar population models and stellar kinematics; e.g. Auger et al. 2010; Treu et al. 2010; Spiniello et al. 2011; Cappellari et al. 2012; Barnabè et al. 2013). The dark matter mass within the considered region is then calculated by subtracting the stellar mass so derived from the total mass measured with the gravitational lens modelling.

6.3 Posterior probability of |$\boldsymbol {\alpha }$| and |$\boldsymbol {f}$|

6.4 Results

In this section, we use the non-detections obtained in this paper in combination with the previously reported detection in the lens system SDSS J0946+1006 to derive statistical constraints on the substructure mass function within a region of about 32 kpc² on average around the Einstein radius, for substructure masses between 4 × 10⁶ and 4 × 10⁹ M_⊙. The results of this analysis are shown in Fig. 2.

Figure 2.

Joint probability P(α, f|n_s, m, p) contours and marginalized probabilities P(f|n_s, m, p) and P(α|n_s, m, p) assuming a uniform prior (blue lines) and assuming a Gaussian prior in α (black lines). Contours (inside out) are set at levels Δln (P) = −1, −4, −9 from the peak of the posterior probability density. The cross represents the prediction from the CDM model as derived by Xu et al. (2013), the solid grey line is the range of mass fractions allowed by flux-ratio anomalies in the multiply imaged quasar sample analysed by Dalal & Kochanek (2002), while the results derived by Vegetti et al. (2010) from the single detection in the lens system SDSS J0946+1006 with a uniform prior and with a Gaussian prior on the slope α are shown by the orange and magenta points, respectively.

Open in new tabDownload slide

Figure 3.

Best smooth models for the lens gravitational system J0252+0039 observed with HST in the F814W filer.

Open in new tabDownload slide

Figure 4.

Same as Fig. 3 for J0737+3216 as observed in the HST F814W filter (top) and in the HST F555W filter (bottom).

Open in new tabDownload slide

Figure 5.

Same as Fig. 3 for J0956+5100 as observed in the HST F814W filter (top) and in the HST F555W filter (bottom).

Open in new tabDownload slide

Figure 6.

Same as Fig. 3 for J0959+4416.

Open in new tabDownload slide

Figure 7.

Same as Fig. 3 for J1023+4230.

Open in new tabDownload slide

Figure 8.

Same as Fig. 3 for J1205+4910 as observed in the HST F814W filter (top) and in the HST F555W filter (bottom).

Open in new tabDownload slide

Figure 9.

Same as Fig. 3 for J1430+4105.

Open in new tabDownload slide

Figure 10.

Same as Fig. 3 for J1627−0053 as observed in the HST F814W filter (top) and in the HST F555W filter (bottom).

Open in new tabDownload slide

Figure 11.

Same as Fig. 3 for J2238−0754 as observed in the HST F814W filter (top) and in the HST F555W filter (bottom).

Open in new tabDownload slide

Figure 12.

Same as Fig. 3 for J2300+0022 as observed in the HST F814W filter (top) and in the HST F555W filter (bottom).

Open in new tabDownload slide

For a uniform prior probability density on α, we find a mean substructure projected mass fraction of f = |$0.0076_{-0.0052}^{+0.0208}$| at the 68 per cent CL and a substructure mass function slope α <2.93 at the 95 per cent CL. For a Gaussian prior of mean 1.9 and standard deviation 0.1, we infer a fraction of f = |$0.0064^{+0.0080}_{-0.0042}$| and a slope of α = |$1.90^{+0.098}_{-0.098}$| at the 68 per cent CL (see Table 3).

With a Bayes factor of |$2\times \Delta \ln {\cal E} = 0.36$|⁠, we find that there is no positive preference for the model with a uniform prior on the slope (Kass & Raftery 1995). This reflects the fact that with only one detection we have poor constraints on the mass function slope, hence only an upper limit is derived. The mass fraction is instead constrained at the lower limit by the single detection and at the upper limit by the non-detections. These results are robust against assumptions on the prior probability density P( f). Thanks to the larger sample of lenses considered in this paper, we can derive tighter constraints on the substructure mass fraction than previously obtained. However, owing to the large number of non-detections, potentially related to the rather low sensitivity of the data analysed here, we still suffer from a large uncertainty on the substructure mass function slope.

Since detailed studies of the luminous satellites of SLACS lens galaxies have shown that the substructure found in lenses is representative of that measured around non-lens galaxies (Jackson et al. 2010; Nierenberg, Oldenburg & Treu 2013b) we are confident that our results can be generalized to the massive elliptical galaxy population.

7 DISCUSSION AND CONCLUSIONS

Vegetti et al. (2010) reported from the analysis of the single detection in the lens system SDSS J0946+1006 a substructure mass fraction of |$f = 0.0256^{+0.0326}_{-0.0150}$| and |$0.0215^{+0.0205}_{-0.0125}$| for a uniform prior and a Gaussian prior on α, respectively. This is higher but consistent within the error with the fraction reported in this paper. The mean value of our current analysis is lower due to the large number of non-detections which tightly constrains the fraction upper limit. From the joint analysis of the lens system SDSS J0946+1006 and B1938+666, Vegetti et al. (2012) derived instead a fraction of |$f=0.039^{+0.036}_{-0.024}$| and |$0.015^{+0.015}_{-0.009, }$| for a uniform prior and a Gaussian prior on α, respectively. This results is also larger but consistent with our constraints and this difference is also due to the different number of detections and non-detections considered. Thanks to a larger number of detections, Vegetti et al. (2012) were also able to set a constraint on the mass function slope, |$\alpha = 1.0^{+0.8}_{-0.1}$|⁠. In this paper, we have not included the detection in the lens system B1938+666 in order to keep the lens sample uniform in both lens galaxy mass and redshift. In the future, larger samples of strong gravitational lens galaxies will allow us to constrain the substructure mass function as a function of the host galaxy mass and redshift and to tightly constraint both the mass function slope and normalization.

Our current results imply a fraction which is consistent within the error with the mass fraction derived by Dalal & Kochanek (2002) (0.006 < f < 0.07) using measurements of flux-ratio anomalies in multiply imaged quasars. It should be kept in mind, however, that the lens sample analysed by Dalal & Kochanek (2002) has a mean redshift of 0.6, hence a larger mass fraction should be expected relatively to our z = 0.2 sample (Xu et al. 2013).

Using rescaled versions of the Aquarius (Springel et al. 2008) and Phoenix simulations (Gao et al. 2012), Xu et al. (2013) have explored the effect of mass substructure on the flux ratio of multiply imaged quasars in galaxies with masses, ellipticity and redshift comparable to the observed lens galaxy population. While the effect of substructures on the flux of point sources is different from their effect on the surface brightness distribution of extended sources (the former is a function of the second derivative of the potential, while the latter is a function of its first derivative) we can still use some of the Xu et al. (2013) results to compare our findings with the expectations from the CDM model. In particular, Xu et al. (2013) have found the number density of subhaloes in the inner regions of the host haloes to be constant and the subhaloes surface mass density to be roughly 4.72 × 10⁶ h⁻¹ M_⊙ (h⁻¹ kpc)⁻², for subhalo masses between 5.48 × 10⁶ and 5.48 × 10⁹ M_⊙ h⁻¹, a host halo mass of about 10¹³ h⁻¹ M_⊙ and a host halo redshift of z = 0.2. This implies an average substructure mass fraction of f_CDM ≈ 0.0046, for substructure masses between 4 × 10⁶ and 4 × 10⁹ M_⊙ and for substructure positions within the same regions considered in this paper. This is consistent with our results at the 68 per cent CL. We caution, however, that higher quality data are required in order to extend our measurements to lower substructure masses (Lagattuta et al. 2012; Vegetti et al. 2012), where larger discrepancies between the data and the model are expected. Also, high-resolution hydrodynamical numerical simulations are required for a proper comparison that takes into account for the effect of baryons on the survival of mass substructure within massive host galaxies.

Vegetti et al. (2012) have shown that the substructure mass derived under the assumption of a pseudo-Jaffe profile is consisted with the mass derived from the pixelized density corrections. This clearly indicates that the measured substructure mass is a robust quantity that is not significantly affected by the choice of mass profile. Using mock observations of the lens system B1938+666, Vegetti & Vogelsberger (2014) have shown that the assumed profile does not lead to a measurement error on the lensing mass, but that profiles with too low concentrations would not be detectable with the gravitational imaging technique. In practice, this implies that while the detections are not affected by the choice of mass profile, the interpretation of the non-detections could be profile dependent. This is an important issue and in a future paper we will thoroughly investigate the potential biases related to assumptions on the substructure mass profile in different observational scenarios. Vegetti et al. (2010, 2012) have shown instead that the galaxy subtraction procedure does not affect the substructure detections and their interpretation.

In conclusion, we have applied the gravitational imaging technique by Vegetti & Koopmans (2009a) to a sample of 11 gravitational lens galaxies from the SLACS survey with the aim of constraining the substructure mass function for host lens galaxies at a mean redshift of 0.2 and with a mean velocity dispersion of 270 km s⁻¹. Our main results can be summarized as follows:

• given our criteria for the significant detection of a mass substructure, no new dark mass substructure has been identified;

• for each lens we have calculated the mass sensitivity function as a function of the substructure position on the lensed image plane within relatively tight regions around the critical curves. We have found that, given the quality of the considered data, the lowest detectable mass typically varies on average between 0.04 × 10¹⁰ and 0.14 × 10¹⁰ M_⊙ depending on the lens. In some cases, for specific substructure positions the detection threshold, it can be as low as 0.01 × 10¹⁰ M_⊙;

• by statistically combining these non-detections with the detection of a dark-matter-dominated substructure in the lens system SDSS J0946+1006 (Vegetti et al. 2010), we have inferred a projected mass fraction in substructure which is consistent with observations based on flux-ratio anomalies of multiply imaged quasars (Dalal & Kochanek 2002) and with the predictions from N-body CDM simulations (Xu et al. 2013) at the 68 per cent CL.

In the near future, thanks to high-resolution adaptive optics (e.g. from the Strong lensing at High Angular Resolution Programme, SHARP; Lagattuta et al. 2012) and to highly structured Einstein rings observed with the HST in the ultraviolet (e.g. HST observing program 12898), we expect to extend our sensitivity to substructure masses as low as 10⁶ M_⊙, and therefore test the CDM paradigm in a mass regime where most of the discrepancy between the CDM predicted subhalo mass function and the observed luminosity function is found.

During part of this work SV acknowledged the support of a Pappalardo Fellowship at the Massachusetts Institute of Technology. TT acknowledges support from the Packard Foundation in the form of a Packard Research Fellowship. Support for programs 10494, 10798 and 11202 was provided by NASA through grants from the Space Telescope Science Institute. SV is grateful to Paul Schechter and John McKean for insightful discussions, Anna Nierenberg, Sherry Suyu for discussions on the substructure spatial distribution and truncation radius and Dandan Xu for discussions on the comparison between observations and expectations. LVEK, SV and TT thank the Kavli Institute for Theoretical Physics for hosting the 2012 First Galaxies and Faint Dwarfs workshop, during which important parts of this work were made. This work is based on observations made with the NASA/ESA Hubble Space Telescope, obtained from the data archive at the Space Telescope Science Institute. STScI is operated by the Association of Universities for Research in Astronomy, Inc. under NASA contract NAS 5-26555.

REFERENCES