Open Access
Issue
A&A
Volume 703, November 2025
Article Number A285
Number of page(s) 14
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202555967
Published online 21 November 2025

© The Authors 2025

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

A concordance of theoretical and observational evidence suggests that most matter in the Universe is non-baryonic dark matter. The nature of dark matter is unknown, but it is understood to be fundamental to the formation and evolution of structures in the Universe (see Wechsler & Tinker 2018 for review).

Cold dark matter (CDM) models of structure formation, in which dark matter consists of a very weakly interacting particle beyond the Standard Model, succeed in predicting large-scale structure and the global properties of galaxies at z = 0 (see Somerville & Davé 2015 for review). Hypothetical particle candidates include axions and neutralinos (Bergström 2000), which are yet undetected in ground-based experiments (e.g. Bartram et al. 2021; Aalbers et al. 2023). CDM models of structure formation are largely untested on mass scales ≲1010 M. In this regime, many low-mass dark matter haloes are predicted to form and follow a well-defined mass function (Moore et al. 1999; Springel et al. 2008; Frenk & White 2012) and density profile (Wang et al. 2020).

Possible alternatives to CDM include warm dark matter (WDM), in which the hypothetical dark matter particle has a non-negligible thermal velocity, and self-interacting dark matter (SIDM), in which the particles have a non-zero self-interaction cross-section. In contrast to CDM, WDM produces low-mass haloes that are less concentrated and less abundant (Lovell et al. 2012; Lovell 2020), and their haloes have different density profiles (Despali et al. 2020). On the other hand, SIDM produces haloes with more diverse properties than CDM (Despali et al. 2019; Nadler et al. 2023). Therefore, the abundance and properties of low-mass haloes can constrain possible dark matter particle candidates.

In the absence of a detectable electromagnetic signature (e.g. see Sawala et al. 2016), a promising way to quantify low-mass haloes at cosmological distances is through their influence on gravitationally lensed images. This is now a well-established method to constrain the shape of the halo mass function and the properties of hypothetical dark matter particles (Birrer et al. 2017; Vegetti et al. 2012, 2014, 2018; Vegetti et al. 2024; Ritondale et al. 2019; Gilman et al. 2020a; Powell et al. 2023). Recent work suggests it may also be possible to constrain the concentration of low-mass perturbers (Vegetti & Vogelsberger 2014; Gilman et al. 2020b; Minor et al. 2021; Zhang et al. 2022; Amorisco et al. 2022; Despali et al. 2025; Minor 2025; Enzi et al. 2025; Tajalli et al. 2025).

Here, we report on our investigation of SDP 81, a dusty star-forming galaxy at z = 3.042 (Eales et al. 2010) that is lensed by a massive elliptical galaxy at z = 0.2999 (Negrello et al. 2014). Two previous studies have reported finding evidence for two different low-mass haloes in (or along the line of sight to) the lensing galaxy halo using Atacama Large Millimetre/sub-millimetre Array (ALMA) observations at ≈30 mas angular resolution (Hezaveh et al. 2016; Inoue et al. 2016). Currently, Hezaveh et al. (2016) represents the only widely cited claim of detection of a low-mass halo in interferometric data, and the third such claim overall (Vegetti et al. 2010, 2012). While the two earlier reported detections in optical/infrared data have been independently verified with different methodologies (Minor et al. 2021; Sengül et al. 2022; Nightingale et al. 2024; Ballard et al. 2024; Minor 2025; Despali et al. 2025; Enzi et al. 2025; Tajalli et al. 2025), the Hezaveh et al. (2016) sub-halo inference has not yet been reproduced.

In this paper, we re-analyse the ALMA observations of SDP 81 to test the robustness of previously reported low-mass haloes. We begin by giving an overview of the ALMA observations and the calibration procedure in Section 2. In Section 3, we explain the lens model parameterisation and methodology. In Section 4, we compare our lens models and perform several tests of our methodology before discussing the implications of our findings in Sections 6 and 7. Throughout, we assume a Planck Collaboration VI (2020) flat ΛCDM cosmology1.

2. Data

We obtained band 6 and 7 calibrated visibilities of SDP 81 from the ALMA archive taken as part of the 2014 ALMA Long Baseline Campaign (ALMA Partnership 2015). The pipeline calibration and subsequent reduction of the band 7 data were performed using Common Astronomy Software Applications (CASA; McMullin et al. 2007) and described by Stacey et al. (2024). The band 6 data were reduced with the same procedure.

We performed no further time or spectral averaging, resulting in data sets of 16 spectral channels of 0.5 GHz2 and 2 s integration per visibility for the band 7 data and 8 channels of 1 GHz and 6 s integration for the band 6 data. This channel bandwidth results in ≤1% and ≤3% amplitude loss, respectively, at 3 arcsec from the phase centre, when adopting the formula from Bridle & Schwab (1989) for a square bandpass. The time integrations of 2 and 6 s result in ≪1% amplitude loss at 3 arcsec (Bridle & Schwab 1989). We set the absolute data weights and noise standard deviations by subtracting subsequent visibilities in a time series per baseline and computing the average in 1200 s chunks for each spectral channel.

For reference, ‘Clean’ images of the data are shown in Fig. 1, which were made by convolving the maximum a posteriori sky models (detailed in the subsequent section) with a Gaussian fit to the central lobe of the dirty beam (i.e. the synthesised beam) and adding it to an image of the residual visibilities. This forward modelling approach results in fewer artefacts in the imaging than deconvolution with the commonly adopted CLEAN algorithm (see also Stacey et al. 2025).

thumbnail Fig. 1.

ALMA images of SDP 81 with natural weighting of the visibilities created via forward modelling. The synthesised beam FWHM with natural weighting is 26 × 33 mas and 33 × 42 mas for bands 7 (left) and 6 (right), respectively.

All the data are used in the subsequent analysis, although the data from the two bands are treated separately. For comparison, Hezaveh et al. (2016) binned the data in 12 m bins and only used visibilities with lower measured instrumental noise levels, amounting to 40% of the total data. No additional time averaging was done, and the minimum and maximum uv distance of the data are the same (Y. Hezaveh, private communication). Inoue et al. (2016) performed their analysis on images produced with the CLEAN algorithm after tapering the data to an effective synthesised beam size of 0.17 arcsec.

3. Lens modelling

We employed PRONTO, a Bayesian pixellated lens modelling technique appropriate for interferometric data introduced by Powell et al. (2021, 2022; see also Vegetti & Koopmans 2009; Rybak et al. 2015a,b; Rizzo et al. 2018; Ritondale et al. 2019; Ndiritu et al. 2025). An analysis of the smooth model for these data is detailed by Stacey et al. (2024) and extended here to include additional lens model components. Table 1 gives a summary explanation of these model components. The Bayesian inference methodology is detailed in Stacey et al. (2024).

Table 1.

Summary explanation of the various models tested on the band 7 data.

3.1. Smooth model

As described by Stacey et al. (2024), the global mass distribution of the lens is parameterised by a singular power law ellipsoidal density profile with external shear. The convergence of the power law density profile is described by

κ l ( R ) = 3 γ 2 ( κ 0 , l R ) γ 1 , $$ \begin{aligned} \kappa _l(R) = \frac{3-\gamma }{2} \left(\frac{\kappa _{0,l}}{R}\right)^{\gamma -1}, \end{aligned} $$(1)

where R is the elliptical radius (R2 = q2x2 + y2), γ is the 3D logarithmic density slope, q is the axis ratio, and κ 0 , l = q θ E $ \kappa_{0,l} = \sqrt{q}\ \theta_{E} $, where θE is the Einstein radius. The external shear is described by a strength, Γ, and position angle, θΓ. The elliptical power law plus external shear component is denoted as ‘PL’.

Following Stacey et al. (2024), we allowed for multipole perturbations that add angular variations to the PL model, which are described by

κ m ( R 0 , θ ) = κ 0 , m ( R 0 1 arcsec ) ( γ 1 ) [ A m sin ( m θ ) + B m cos ( m θ ) ] $$ \begin{aligned} \kappa _m (R_0,\theta ) = \kappa _{0,m} \left(\frac{R_{0}}{1\,\mathrm{arcsec}}\right)^{-(\gamma -1)} \left[A_m \sin (m\theta ) + B_m\cos (m\theta )\right] \end{aligned} $$(2)

in polar coordinates for multipole order m, where γ is the slope of the density profile in Eq. (1), R0 is the (spherical) radius (R02 = x2 + y2) in terms of a scale radius of 1 arcsec, and κ0, m is the normalisation. Am and Bm are dimensionless free parameters. We consider m = 3, 4, and denote this multipole model component as ‘MP’.

3.2. Sub-halo model

We modelled a low-mass halo as a sub-halo of the primary lensing galaxy (i.e. at the same redshift as the lens), considering a spherical pseudo-Jaffe (denoted as ‘PJ’; Dalal & Kochanek 2002) profile. We used this model for two reasons. Firstly, we wished to compare our results to the previous analysis of these data by Hezaveh et al. (2016), who used a PJ model; secondly, because a combination of tidal stripping and feedback is predicted to result in a mass density profile close to isothermal (Heinze et al. 2024).

The PJ component is a truncated isothermal density profile. We followed the definition of Despali et al. (2018), whereby the 2D truncation radius (Rt) is a function of the sub-halo mass (Ms) and its distance from the main halo (i.e. the primary lens). The convergence is defined as

κ s ( R 0 ) = κ 0 , s R t [ R 0 1 ( R 0 2 + R t 2 ) 1 / 2 ] . $$ \begin{aligned} \kappa _{s}(R_0) = \kappa _{0,s}\ R_t [R_0^{-1} - (R_0^2 + R_t^2)^{-1/2}]. \end{aligned} $$(3)

The mass is defined as

M s = 2 π Σ c R t 2 κ 0 , s , $$ \begin{aligned} M_{s} = 2\pi \Sigma _{c} R_t^2 \kappa _{0,s}, \end{aligned} $$(4)

where Σc is the critical surface mass density at the redshift of the lens and κ0, s is the sub-halo normalisation. The truncation radius is described by

R t = R 0 π κ 0 , s 2 ξ κ 0 , l , $$ \begin{aligned} R_t = R_0 \sqrt{\frac{\pi \kappa _{0,s}}{2 \xi \kappa _{0,l}}}, \end{aligned} $$(5)

where κ0, l is the normalisation of the PL component and ξ is the impact parameter (here assumed to be 3).

4. Analysis of real data

The goal of this work is to compare our lens modelling analysis with previous studies by Hezaveh et al. (2016) and Inoue et al. (2016) (H16 and I16, hereafter). First, we report the results of our smooth modelling, including flexibility for angular and radial variations in the mass distribution (Section 4.1). Second, we test the reported detections of low-mass haloes in these data by including sub-halo model components (Section 4.2). We then test the effects of regular and adaptive source gridding (Section 4.3).

The summary statistics for the posterior distributions of the lens model parameters from MULTINEST are given in Table 2, and the maximum a posteriori models are shown in Fig. 2. Fig. A.1 is a corner plot showing the two-dimensional posterior probability distributions. Heat maps comparing the log Bayes factors of the various models tested are shown in Figs. 3a–5b. Extended tables with summary statistics for the various models tested are in the appendix.

thumbnail Fig. 2.

Maximum a posteriori lens models for the band 7 data. First column: Lens model. Crosses show lens and sub-halo positions; critical curves are in black. Second column: Difference in convergence relative to the PL model. Contours of the sky model are in grey. Third column: Difference in magnification relative to the PL model; contours of the sky model in black. Fourth column: Image of the residuals (data–model) normalised by the RMS noise (σ), with the mask shown by the dashed outline. Last column: Reconstructed source with Gouraud interpolation and caustics in black. The locations of I16 and H16 purported sub-haloes are shown with a triangle and a circle, respectively.

Table 2.

Lens model parameters for the band 7 data, their 1σ uncertainties, and priors.

As detailed in Stacey et al. (2024), to prevent overfitting, it is necessary to impose regularisation when solving for the source surface brightness. We considered three forms of regularisation: curvature, gradient and area-weighted gradient methods (Suyu et al. 2006; Stacey et al. 2024). Additionally, we considered different numbers of pixels cast into the source plane to impose another form of regularisation (Bayer et al. 2023). We performed Multinest sampling of the PL model for the band 7 data, with 1, 2, 3 and 4 cast pixels, for each of the different regularisation types. We found that 3 pixels and gradient regularisation resulted in the highest Bayesian evidence, and this is what we used for all the subsequent analysis.

4.1. Smooth model

4.1.1. Power law model

Our inferred PL model parameters are similar to the PL models inferred from previous studies of this lens using the visibility data from baselines < 2 km (Rybak et al. 2015a) or ‘Clean’ images of the data at ≈30 mas resolution (Dye et al. 2015; Inoue et al. 2016; Tamura et al. 2015; Wong et al. 2015). We find a power law slope, ellipticity, ellipticity position angle, shear, and shear position angles consistent within 1σ with Dye et al. (2015) who used a pixellated modelling of the Briggs-weighted Clean band 6 and band 7 images. We also find ellipticity, ellipticity position angle, shear, and shear position angle consistent within 1σ of those inferred by Inoue et al. (2016) for a singular isothermal ellipsoid and external shear lens model, and using point-sources extracted from a Clean image as constraints.

H16 reported the summary statistics only for PL+MP+PJ, not for PL and PL+MP individually, so a fair comparison with our PL model cannot be made. We compare our PL+MP component parameters in the following section.

4.1.2. Multipoles

PL+MP was preferred with a log Bayes factor (ln𝒦) of 28 relative to PL (i.e. ln𝒦 = lnℰPL + MP − lnℰPL, where ℰ is the Bayesian model evidence), suggesting angular complexity in the lens model beyond the PL model. The parameter posteriors and summary statistics differ between PL and the PL component of PL+MP (see Fig. A.1 and Table 2). Several PL model parameters shift by ≥1σ with the inclusion of MP, most significantly the lens centre and position angles of the external shear (by > 3σ). A comparison of PL and PL+MP for this lens was previously reported by Stacey et al. (2024) for these band 7 data with the same methodology used here. We refer to Stacey et al. (2024) for a more detailed analysis of the multipole structure.

Besides this work and that of Stacey et al. (2024), the only other study in the literature to include an MP component in their lens model was H16. Despite shifts in PL component parameters with the inclusion of multipoles, we find significant disagreement with the parameters found by H16. We compute the mass enclosed within 10 kpc for our inferred lens model posterior and find log(M10 kpc/M) = 11.64 (see Table 2), which is consistent within 3σ of the enclosed mass reported by H16, accounting for differences in cosmology. However, after accounting for differences in parameter definitions, we find that the PL model we infer has significantly different (> 10σ) external shear strength and ellipticity, as well as a significantly different ellipticity position angle. The external shear of the H16 PL is an order of magnitude smaller, the ellipticity is significantly higher (> 10σ), and the ellipticity position angle is 70 deg different. Furthermore, H16 inferred significantly higher m = 3 multipole coefficient amplitudes than we find, albeit with much larger uncertainties. After accounting for different multipole parameterisations (Y. Hezaveh, private communication), the H16 mean multipole parameters are A3 = 0.003 ± 0.003; B3 = 0.013 ± 0.003; A4 = 0.006 ± 0.005; B4 = 0.003 ± 0.006 for our parameterisation (Eq. 2). We find m = 3 multipole amplitudes ( η m = A m 2 + B m 2 $ \eta_m=\sqrt{A_m^2+B_m^2} $), a factor of 3.3 ± 0.8 times smaller than H16. The m = 4 amplitudes and position angles (ϕm = (1/m)arctanBm/Am) are consistent within their 1σ uncertainties.

To confirm that we did not converge to a local minimum in the posterior parameter space, we re-ran our initial PL optimisation beginning from the lens parameters inferred by H16. The optimisation ultimately converged to the same maximum a posteriori parameters that we initially inferred. Additionally, we re-ran the MULTINEST sampling with all uniform priors large enough to encompass the mean and uncertainties of the MP and PL parameters inferred by H16, but it did not change the inferred evidence or summary statistics.

H16 reported that a PL model similar to ours and others reported in the literature (Dye et al. 2015; Inoue et al. 2016; Tamura et al. 2015; Wong et al. 2015) gave a comparably good fit to the data and did not significantly alter their sub-halo inference. Therefore, as an additional test, we fixed the lens model to our maximum a posteriori PL and the mean MP parameters reported by H16 (denoted PL+MP) and optimised only the regularisation hyperparameter to allow the source to focus. This allowed us to test these multipoles directly and investigate how they might impact the sub-halo inference. This model resulted in larger image-plane residuals (> 5σ) compared to our best-fit PL+MP (< 3σ, see Fig. 2). We then considered whether a different PL component could produce a better model in combination with MP by allowing the PL component to optimise. Still, this resulted in > 5σ features in the image-plane residuals. To quantify this goodness of fit, we ran MULTINEST for PL+MP and our best fit PL+MP model with all model parameters fixed except for the source regularisation hyperparameter so that we can compare the model evidence. We inferred a log Bayes factor of −697 for PL+MP relative to our best PL+MP, confirming that the H16 MP parameters are a poorer fit to the data.

4.2. Sub-halo search

4.2.1. Sub-halo search with our smooth model

To determine whether the data favour a sub-halo, we included a parametric sub-halo component in our models (PJ, as detailed in Section 3.2). We explored the posterior parameter space with MULTINEST, allowing the PJ to be anywhere on the lens plane with a log-prior on the mass of log(Msub/M)∈[6, 10.4]. The smooth model parameters were free with the same priors as for the smooth model inference. We found that the PL+PJ is preferred by a log Bayes factor of 22 relative to PL. While the PJ component is constrained in mass and location, it is in a different location than reported in the previous studies and coincides with a positive convergence feature in the PL+MP model (Fig. 2). The parameter summary statistics for PL+PJ show that the PL component parameters move closer to the PL component parameters of PL+MP, suggesting that the PL and PJ components conspire to induce a similar, non-localised lensing effect. This is also indicated by the apparent convergence and magnification differences in Fig. 2.

Previous work has found that incorrect conclusions about the presence of sub-haloes may be inferred if angular structure is not included in the smooth model (O’Riordan & Vegetti 2024). To avoid this, we include a multipole component in subsequent sub-halo searches to account for angular as well as localised mass structure. For posterior sampling with MULTINEST, the smooth model and the PJ component parameters were free with the same priors as before. Comparing the two models, PL+MP and PL+MP+PJ, we found no difference in Bayesian evidence (see Table 2), meaning the data do not favour one model over the other. The PL and MP component parameters were consistent across both models, with differences falling within 1σ uncertainties (see Fig. A.1 and Table 2). The PJ component position was not constrained (see Fig. A.1), and we infer an upper mass limit (at the 86th percentile) of 108.1 M. Therefore, we infer no evidence in favour of a sub-halo when multipoles are included in the model.

This is in contradiction to the inference of the ≈109 M sub-halo by H16. Therefore, to more directly test the sub-halo reported by H16, and to ensure we sufficiently explore the relevant parameter space, we fixed the sub-halo’s position to the coordinates given by H16 (this model is denoted PL+MP+PJ(p)). The sub-halo mass was left as a free parameter, with the same prior bounds as before. The PL+MP parameters and regularisation remained free with the same priors as before. In this case, we found that the PJ mass tended to the lower edge of the prior and resulted in no improvement in the Bayesian evidence over PL+MP. We infer an 86th percentile upper limit of 107.8 M, indicating that the data do not favour a sub-halo in this location.

To better quantify the level of exclusion of the H16 sub-halo, we directly fixed both the sub-halo mass and position to that found by H16 while leaving the PL+MP parameters free to allow them to compensate for the sub-halo component. This model is denoted PL+MP+PJ(p,m). We note that the PJ parameterisation implemented by H16 is slightly different from our analysis: the truncation radius used by H16 is defined relative to the ratio of velocity dispersions of isothermal spheres of the masses of the primary lens and PJ. This results in a deflection equivalent to a less massive sub-halo in the Despali et al. (2018) definition. To ensure a fair comparison, we implemented the exact parameterisation used by H16 (Y. Hezaveh, private communication). With MULTINEST posterior sampling, we found that the PL+MP+PJ(p,m) model is disfavoured relative to PL+MP, with a log Bayes factor of 13. Still, we do not find agreement with the PL or MP component parameters of H16. A heat map showing the log Bayes factors of the models tested on the band 7 data is shown in Fig. 3a. Summary statistics for PL+MP+PJ(p) and PL+MP+PJ(p,m) lens models are reported in Table A.1.

thumbnail Fig. 3.

Heat map showing the log Bayes factors of the models tested on the real band 7 data (y-axis model relative to x-axis model). a: Models tested on the real band 7 data. b: Models tested on the real data with H16 multipole coefficients (MP). The smooth models are fixed, and the regularisation is the only free parameter.

While the band 7 data have a better signal-to-noise ratio and angular resolution, we note that H16 found improved evidence for a sub-halo by combining the band 6 and 7 data sets. We explicitly choose not to combine these data sets due to the different flux density, source structure (see Rybak et al. 2015a) and possible astrometric offsets. Still, comparing our findings for the two data sets helps test the sensitivity of the sub-halo detection to phase and amplitude errors in the data. Therefore, we also searched for evidence of a sub-halo in the band 6 data, but we found no improvement in the Bayesian evidence of PL+MP+PJ or PL+MP+PJ(p) relative to PL+MP (i.e. ln𝒦 = 0). We infer 86th percentile upper limits of 108.0 M for the mass of a sub-halo in the location inferred by H16, indicating that neither band 6 nor 7 data sets favour the H16 sub-halo.

4.2.2. Sub-halo search with the H16 smooth model

As a final test, we investigated whether the different multipole coefficients used by H16 could impact the sub-halo inference. Specifically, we considered whether the combination of the H16 smooth model parameters and sub-halo might conspire to improve the overall model fit. Using the band 7 data, we fixed the smooth model components to our maximum a posteriori PL and multipoles reported by H16 (PL+MP), leaving only the regularisation free to allow the source to focus (as described in Section 4.1.2). We added a PJ component with mass prior bounds of log(Msub/M)∈[6, 10.4] and position prior bounds of the image plane. This allowed us to test whether a sub-halo with similar properties to that found by H16 could be recovered under these assumptions.

With this setup, PL+MP+PJ is strongly favoured relative to PL+MP with a log Bayes factor of 365. This is not unsurprising, because we established that PL+MP produces > 5σ residuals, so it is expected that the additional degrees of freedom will act to reduce these residuals. However, the inferred mass and position of the PJ component are log Msub = 9.93 ± 0.03, xs = −1.25 ± 0.04, ys = −1.64 ± 0.04. These differ from the H16 sub-halo by more than 10σ, indicating that it is not the same model inferred by H16. Therefore, we also tested the H16 sub-halo directly by fixing both its mass and position. This model, PL+MP+PJ(p,m), was preferred by a log Bayes factor of 75 relative to PL+MP. This is equivalent to 12σ, higher than the 6σ inferred by H16.

Nevertheless, as shown in Fig. 3b, all these models based on the H16 smooth lens model (PL+MP) are very strongly disfavoured relative to our maximum a posteriori PL+MP model, with log Bayes factors in the range −332 to −697. Therefore, these results still indicate no statistical preference for a sub-halo. The relative increase in evidence for the H16 sub-halo in combination with the H16 MP component indicates that the MP component may have led to a spurious sub-halo inference.

4.3. Regular vs adaptive source gridding

A key methodological difference between this work and H16 lies in the choice of the source reconstruction grid. In our analysis, we use a magnification-adaptive source grid, as described in Vegetti & Koopmans (2009), whereas H16 employed a regular grid with an arbitrary 10 mas pixel scale. Each approach has its advantages and limitations. Adaptive grids provide greater flexibility in highly magnified regions, potentially leading to overfitting if not properly regularised. Regular grids, on the other hand, may either overfit or underfit depending on the selected pixel scale, which is not dynamically adjusted to the lensing geometry. As a result, the sensitivity to sub-halo perturbations may differ between the two methods: an adaptive grid might smooth over small perturbations that a regular grid could capture, or vice versa, depending on the local structure of the lensing potential and the chosen source grid resolution.

To assess whether the choice of source grid influences sub-halo detection, we repeated our analysis using a regular source grid with a fixed 10 mas pixel scale, matching the setup used by H16. This allowed us to test whether similar evidence for the H16 sub-halo (or any sub-halo) could be recovered in this scenario.

Using the regular grid, we found that the PL+MP+PJ model was disfavoured relative to PL+MP, with a log Bayes factor of 5. In this case, the PJ component was unconstrained in both mass and position. We also tested a model with the PJ fixed at the position of the H16 sub-halo. This configuration was similarly disfavoured, with a log Bayes factor of 6, and yielded an 86th percentile upper limit on the PJ mass of 107.6 M. The model where both the mass and position of the PJ were fixed to the H16 values (PL+MP+PJ(p,m)) was more strongly disfavoured than for the adaptive grid, with a log Bayes factor of 70 relative to PL+MP.

The larger overall Bayes factors found using the regular grid (compared to the adaptive grid) likely reflect the stronger effective regularisation imposed by the fixed 10 mas resolution, in contrast to the more flexible Delaunay-based adaptive grid. For visual comparison, Fig. 5a shows a heat map of the log Bayes factors across all tested models for a regular grid.

As a further investigation, we performed these same tests with the smooth model inferred by H16 (PL+MP), to exclude the possibility that this model succeeds for a regular grid but not an adaptive grid. We fixed the smooth model parameters, with the only free parameter being the regularisation to allow the source to focus (i.e. the same test telegraphed by Fig. 3b). As shown in Fig. 5b, these models are also all very strongly disfavoured relative to our best smooth model for a regular source grid (ln𝒦 = −484 to −879).

As for the adaptive grid, the evidence increases with the inclusion of the H16 sub-halo relative to the H16 smooth model (by a factor of 92). While the Bayes factors differ, the overall result is the same as for the adaptive grid. Therefore, we do not find evidence that the H16 sub-halo was an artefact produced by the regular source grid, nor do we find evidence that the difference in source grid type could change our interpretation.

5. Analysis of mock data

We tested a number of our findings with realistic mock data. One of our findings was angular structure, and an apparent degeneracy between that and a single massive sub-halo, which we explore in Sections 5.1 and 5.2. Additionally, we do not find agreement with the previous study by H16 with a similar approach, so we tested whether we could detect such a sub-halo (Section 5.2) and test the effects of phase errors (Section 5.3).

We made realistic mock data by overwriting the real band 7 visibility data with a Fourier transform of the maximum a posteriori lens model and adding Gaussian noise per channel per correlation as measured in the real visibility data in 1200 s chunks. We modelled the mock data with the same procedure as for the real band 7 data, including the same priors for Multinest sampling (Table 2).

5.1. Smooth model

With mock data generated from our maximum a posteriori PL+MP lens model, we verify that the inferred multipole amplitudes can be reliably recovered. The posterior distributions obtained from the mock data are consistent with those from the real data, within 1σ (see Fig. 6 and Table A.3). We also find a log Bayes factor of 27 in favour of the PL+MP model over the PL model, closely matching the value of 28 obtained from the real observations. This agreement indicates that the detected multipole features are not the result of noise artefacts but are instead genuine structures present in the data.

5.2. Sub-halo search

Firstly, we assessed whether the similar evidence for PL+MP and PL+PJ for the real data could be reproduced. We suspected that this apparent degeneracy could be because a single massive sub-halo can induce a non-localised lensing effect and mimic a multipole effect. Therefore, we tested this using the mock data with ground truth PL+MP, introduced in Section 5.1. From sampling with MULTINEST, we found PL+PJ is also favoured by 24 relative to PL; a similar, but slightly lower significance than for PL+MP (Fig. 4a; summary statistics are given in Table A.3). This is consistent with what we found for the real data. As there is no PJ in the ground truth, this confirms that a single massive sub-halo can mimic the lensing effect of multipoles. In this case, the PJ is constrained at the north-east corner of the image plane (edge of the prior) and is very massive (≈1010 M), where it has a non-localised effect. Although these are not the same mass and location inferred for the real data, most of the PL component parameters move closer to those in the PL+MP model, as we observed for the real data.

thumbnail Fig. 4.

Heat map showing the log Bayes factors of the models tested on mock band 7 data (y-axis relative to x-axis). a: Models tested on mock band 7 data where the ground truth is PL+MP. b: Models tested on mock band 7 data where the ground truth is PL+MP+PJ (p,m).

We then considered whether our lack of evidence for the H16 sub-halo was due to methodological differences (e.g. source regularisation or posterior inference) by performing the sub-halo search, as in Section 4.2, on mock data with ground truth PL+MP+PJ(p,m). The mock lens model was generated from the maximum a posteriori PL+MP source, lensed by the maximum a posteriori PL+MP plus a PJ component. The PJ sub-halo had the mass, location and parameterisation as reported by H16. We performed sampling with MULTINEST of the mock data with and without the H16 sub-halo (including multipoles) to compute the models’ evidence, marginalising over the source and smooth model parameters. A Heat map comparing the log Bayes factors is shown in Fig. 4b.

For the mock data where the ground truth is PL+MP+PJ(p,m), we found that the PL+MP+PJ model has the highest evidence relative to PL and PL+MP. PL+MP is favoured by a log Bayes factor of 13 relative to PL, while PL+MP+PJ(p,m) is favoured by a log Bayes factor of 29 relative to PL+MP. This suggests that we would find evidence in favour of the H16 PJ with a log Bayes factor of 29 relative to PL+MP, should it exist in the real data. The posterior mass and position of the sub-halo are consistent within 3σ uncertainties of the ground truth (see Table A.4 for summary statistics). Therefore, we conclude that our inability to find evidence in favour of the H16 sub-halo in the real data is not a consequence of different methodologies.

5.3. Phase errors

In the analysis thus far, we have only considered the effect of instrumental noise on the model evidence. However, antennas may have residual systematic errors as a result of transferring gain solutions from the phase calibrator to the target. The magnitude of these phase errors3 depends on the specifics of the observations (e.g. distance from target to calibrator, weather conditions, flux density of the calibrator and target, baselining precision). Residual phase errors can be corrected through self-calibration, provided there is a sufficient signal-to-noise ratio per antenna (e.g. Spingola et al. 2018). However, in the case of SDP 81, the low signal-to-noise ratio per antenna means that self-calibration solutions are only possible for long time intervals (see Stacey et al. 2024). Therefore, any phase errors on shorter timescales (≲1 h) are unknown and uncorrected.

With this in mind, it is plausible that phase errors could induce artefacts in the data, reducing the sensitivity of the data to a sub-halo, as this method relies on the surface brightness accuracy of lensed images. Therefore, we tested whether phase errors could reduce the evidence for a sub-halo by generating mock data sets. These mock data sets had a ground truth of PL+MP+PJ(p,m) (as detailed in Section 5.2) with systematic phase offsets added to each antenna, drawn from a Gaussian distribution. This assumption may not be representative of true antenna-based phase errors for ALMA that, in addition to being observation-dependent, are also time-dependent and baseline-length-dependent (see Maud et al. 2022, 2023). However, this is a simple test of the effect of phase errors on our mock data and is the same assumption made by H16.

As the actual phase errors in the real data are unknown, we tested standard deviations of 5 and 10 deg. We modelled the mock data with PL+MP and PL+MP+PJ. In both cases, we recover strong evidence in favour of the sub-halo. The log Bayes factor in favour of PL+MP+PJ is 29 relative to PL+MP with 5 deg phase errors, consistent with the mock without phase errors (see Fig. 4b), while the evidence decreases to 25 for 10 deg. We also created mock data for a standard deviation of 20 deg, but the resulting decoherence made it visibly inconsistent with the real data.

As the evidence in favour of PL+MP+PJ decreases only mildly when random phase errors are included, we conclude that these do not have a significant impact on any sub-halo inference and differences in approach to self-calibration cannot account for the lack of evidence we find for the H16 sub-halo.

On the other hand, one could ask whether random phase errors could lead to spurious evidence in favour of a sub-halo. We tested this scenario by generating mock band 7 data with a ground truth of PL+MP (i.e. no sub-halo) and adding random phase errors, as before. We modelled these mock data with PL+MP and PL+MP+PJ and performed posterior sampling with MULTINEST. However, the log Bayes factors for PL+MP+PJ were < 1 relative to PL+MP for these mock data sets, indicating no preference for a sub-halo. This suggests that random phase errors would not lead to spurious inferences of sub-haloes for the real data, also indicating that the H16 inference cannot be attributed to this origin.

6. Discussion

6.1. Sub-halo or multipoles

The models of the real band 7 data with the highest evidence include multipoles. PL+PJ is favoured relative to the PL model, but not relative to PL+MP+PJ. The uncertainties in the log Bayes factors reported by MULTINEST are typically < 1, but, as explored by Nelson et al. (2020), these uncertainties may be underestimated, and log Bayes factors of ≈3 may not be conducive to reliable model selection. The log Bayes factor of the PL+MP and PL+PJ models is 6, more than what is commonly ascribed to stochasticity.

However, we find that the multipole coefficient posteriors are consistent within 1σ for models tested on the real data and mock data with ground truths PL+MP+PJ and PL+MP (Fig. 6), suggesting that the multipoles are largely stable to the inclusion of a sub-halo in both the model and ground truth. Additionally, with mock data with ground truth PL+MP, we found PL+PJ produced similar (but slightly lower) log Bayes factors than PL+MP, consistent with the real data, suggesting that a massive PJ at the edge of the field can induce a non-localised effect to compensate for angular structure. We note that the convergence and magnification differences of PL+MP and PL+PJ shown in Fig. 2 have some similar features, which further hints that the PJ produces similar deflections in combination with the smooth model parameters. These findings favour the multipole scenario and suggest that a large sub-halo (in combination with the smooth model parameters) can induce angular perturbations in the lens potential similar to multipoles. Our results are in line with recent theoretical work that has found that the presence of angular structure induces spurious ‘detections’ of low-mass haloes if not accounted for in the lens model (O’Riordan & Vegetti 2024) and can fully account for the flux ratio anomalies in lensed quasar images commonly attributed to a population of low-mass haloes (Cohen et al. 2024).

thumbnail Fig. 5.

Heat map showing the log Bayes factors of the models tested on the real band 7 data (y-axis model relative to x-axis model), using a regular source grid of 10 mas. a: Same as Fig. 3a, using a regular source grid. b: Same as Fig. 3b, using a regular source grid.

thumbnail Fig. 6.

Multipole coefficient mean and 1σ uncertainties for the real and mock data. Yellow and red circles show the posteriors for the mock data with ground truth PL+MP and PL+MP+PJ, respectively; cyan and black triangles show the posteriors for PL+MP and PL+MP+PJ models tested on the real data. The multipole coefficients are consistent within their 1σ uncertainties.

6.2. Differences with previous work

We do not find evidence in favour of the sub-haloes reported by I16 or H16 in our re-analysis of these data. However, our tests on mock data indicate that we would have found evidence in favour of the H16 sub-halo with a log Bayes factor of 29 (equivalent to 7.6σ significance) if it were present, a higher significance than the 6σ found by H16.

I16 used a PL smooth model and a sub-halo parameterised as a homogeneous spherical clump of uniform density with a negative density shell localised to the southern point of the eastern arc. They inferred this model by fitting a PL and external shear to lensed source positions measured from a deconvolved image (made by tapering the data to an effective synthesised beam of 0.1 arcsec) and adopting a trial-and-error approach. While the constrained image position is consistent with a positive multipole convergence feature (Fig. 2), it is very plausible that the differences in our findings are because of the different model assumptions and the number of degrees of freedom available using the full surface brightness of the lensed emission. For example, Powell et al. (2022) noted significant differences in the lens model inferred for the system MG J0751+2716 using a pixellated source reconstruction compared to a previous analysis using point-like image positions (Spingola et al. 2018). Furthermore, recent work by Cohen et al. (2024) show that multipole perturbations with amplitudes of 1% (consistent with our best-fit multipole amplitudes; see also Stacey et al. 2024) can perfectly fit quasar image positions and flux ratios, and O’Riordan & Vegetti (2024) showed that multipoles, if not included in the model, can readily lead to artificial sub-halo constraints.

On the other hand, H16 performed a full forward modelling of the source surface brightness using 40% of the visibility data. We ruled out that the treatment of phase errors led to different results; evidence in favour of a sub-halo did not significantly change for mock data with random phase errors. Another key difference with our analysis is how the source is constructed. We used a Delaunay source grid adapted to the lensing magnification, but H16 used a regular source grid with a pixel size of 10 mas and a different surface brightness prior. Nightingale & Dye (2015) found that a regular source grid can induce systematic biases in the inferred lens model parameters. We did not test whether a 10 mas regular source grid sufficiently samples the source structure of SDP 81, but Nightingale & Dye (2015) found that a finer grid would not resolve these systematic biases. While we found that a 10 mas regular source grid led to inconsistent inferred lens parameters for these data, it did not produce evidence in favour of a sub-halo in any location, so we cannot ascribe the H16 sub-halo inference to regular grid bias.

We inferred a PL model consistent with those inferred with previous work (e.g. Dye et al. 2015), but different multipole coefficients from those inferred by H16. We find that our multipole coefficients are more consistent with the data, and the H16 coefficients are strongly disfavoured. Furthermore, the multipole structure was previously modelled and analysed by Stacey et al. (2024), who found a strikingly similar structure in the stellar isophotes and mass isodensity shapes despite the stellar light not having been used as a prior for the lens model.

We find that the most plausible origin of the H16 sub-halo is these inferred smooth model parameters. Their different multipole coefficients suggest that the statistical analysis performed by H16 did not identify the global minimum in the smooth model posterior. This is the only scenario we tested where the H16 sub-halo was favoured relative to its corresponding smooth model. In this scenario, the positive inference of a sub-halo is produced as a result of a poor smooth model fit: the additional degrees of freedom afforded by the sub-halo component act to reduce residuals.

Notably, this scenario we ascribe to the H16 sub-halo inference, is not the same problem discussed in Section 6.1 concerning the similar lensing signatures of multipoles and massive sub-haloes. Rather, it exemplifies how the inference of low-mass haloes requires a good smooth model fit. A poor smooth model fit can be characterised by the level of residuals in multiple lensed images For example, Nightingale et al. (2024) set a criterion for the suitability of a lens for substructure analysis based on the level of smooth model residuals: they assessed that > 3σ residuals over multiple lensed images indicated a poor smooth model fit, such that inclusion of a sub-halo component in the lens model would not produce a reliable inference. In general, a poor smooth model fit may indicate non-localised complexity in the mass distribution beyond what can be described by PL+MP (e.g. see Powell et al. 2022), but we find no need for complexity beyond PL+MP to fit these data.

7. Conclusions

The ‘detection’, or lack thereof, of low-mass dark matter haloes via gravitational lensing is a practical approach to testing dark matter models. Several recent works have explored systematic biases in inferring the presence of low-mass haloes, such as the complexity of the primary lens mass distribution (Nightingale et al. 2024; O’Riordan & Vegetti 2024) and source resolution (Vernardos & Koopmans 2022; Nightingale et al. 2024; Ballard et al. 2024; Minor 2025; Ephremidze et al. 2025). It is nonetheless essential that reported detections are independently verified.

Here, we have been unable to reproduce the inference of a sub-halo by Hezaveh et al. (2016) in ALMA data of SDP 81. Using mock data, we have confirmed that we would have found evidence in favour of that sub-halo should it exist with the prescribed properties. After testing several systematics (smooth model, source gridding, phase errors), we find that the most likely cause of this spurious inference was a poorly fitting multipole component, suggesting that their model inference did not find the global minimum in the parameter space. Another sub-halo reported by Inoue et al. (2016) can be explained by several limitations in the earlier methodology and the presence of angular structure in the lens. Therefore, we conclude that these previous results are not robust.

We find that the presence of angular structure in the lens described by multipoles of third and fourth order (see also Stacey et al. 2024) is sufficient to explain the data, and we find no evidence in favour of any individual sub-halo. However, we find agreement with recent works indicating degeneracy between massive sub-haloes and multipoles (O’Riordan & Vegetti 2024). We suggest that future searches for low-mass perturbers via lensed arcs must allow for complex angular structure to avoid false inferences, which may be achieved by exploiting the full surface brightness of the lensed emission. The goodness of fit of this smooth model is an essential prerequisite to sub-halo inference, which may be assessed via the significance and distribution of image-plane residuals (e.g. Nightingale et al. 2024).

Nevertheless, here, we have shown the ability to differentiate between models including a ≈109 M sub-halo, demonstrating that ALMA data of this quality are useful to test the prevalence of sub-haloes in the < 1010 M mass range (see also Despali et al. 2022). In this case, the inferred multipole coefficients are largely stable in the presence of sub-haloes. Forthcoming analysis from 10 strongly lensed dusty star-forming galaxies observed with ALMA at similar angular resolution will give insight into the generalisability of this result.


1

While previous work by Hezaveh et al. (2016) and Inoue et al. (2016) used different cosmologies, the different parameters do not change the angular diameter distances by more than 1%.

2

Four spectral windows of four channels were manipulated into 16 spectral windows of 1 channel so the data could be converted to UVFITS format. We refer to these as 16 spectral channels for simplicity.

3

We consider only phase errors here, which we expect to be more relevant than amplitude errors since phase errors redistribute the surface brightness on the sky.

Acknowledgments

We thank Yashar Hezaveh and Luke Maud for useful discussions. H. R. S., D. M. P. and S. V. have received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement No 758853). S. V. thanks the Max Planck Society for support through a Max Planck Lise Meitner Group. This research was carried out on the High-Performance Computing resources of the Freya cluster at the Max Planck Computing and Data Facility (MPCDF) in Garching, operated by the Max Planck Society (MPG). JPM and DW acknowledge support from the Netherlands Organization for Scientific Research (NWO) (Project No. 629.001.023) and the Chinese Academy of Sciences (CAS) (Project No. 114A11KYSB20170054). This work is based on the research supported in part by the National Research Foundation of South Africa (Grant Number: 128943). The modified bibliography style is courtesy of Chentao Yang (https://github.com/yangcht/AA-bibstyle-with-hyperlink). This research used Astropy, SciPy, NumPy, and Matplotlib packages for Python (Astropy Collaboration 2013, 2018; Virtanen et al. 2020; Harris et al. 2020; Hunter 2007). We made use of ALMA data with project code 2011.0.00016.SV. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ.

References

  1. Aalbers, J., Akerib, D. S., Akerlof, C. W., et al. 2023, Phys. Rev. Lett., 131, 041002 [Google Scholar]
  2. ALMA Partnership (Vlahakis, C., et al.) 2015, ApJ, 808, L4 [Google Scholar]
  3. Amorisco, N. C., Nightingale, J., He, Q., et al. 2022, MNRAS, 510, 2464 [NASA ADS] [CrossRef] [Google Scholar]
  4. Astropy Collaboration (Robitaille, T. P., et al.) 2013, A&A, 558, A33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Astropy Collaboration (Price-Whelan, A. M., et al.) 2018, AJ, 156, 123 [Google Scholar]
  6. Ballard, D. J., Enzi, W. J. R., Collett, T. E., Turner, H. C., & Smith, R. J. 2024, MNRAS, 528, 7564 [Google Scholar]
  7. Bartram, C., Braine, T., Burns, E., et al. 2021, Phys. Rev. Lett., 127, 261803 [CrossRef] [PubMed] [Google Scholar]
  8. Bayer, D., Koopmans, L. V. E., McKean, J. P., et al. 2023, MNRAS, 523, 1326 [NASA ADS] [CrossRef] [Google Scholar]
  9. Bergström, L. 2000, Rep. Progr. Phys., 63, 793 [Google Scholar]
  10. Birrer, S., Amara, A., & Refregier, A. 2017, JCAP, 2017, 037 [CrossRef] [Google Scholar]
  11. Bridle, A. H., & Schwab, F. R. 1989, Wide Field Imaging I: Bandwidth and Time-Average Smearing (San Francisco: Astronomical Society of the Pacific) [Google Scholar]
  12. Cohen, J. S., Fassnacht, C. D., O’Riordan, C. M., & Vegetti, S. 2024, MNRAS, 531, 3431 [NASA ADS] [CrossRef] [Google Scholar]
  13. Dalal, N., & Kochanek, C. S. 2002, ApJ, 572, 25 [Google Scholar]
  14. Despali, G., Vegetti, S., White, S. D. M., Giocoli, C., & van den Bosch, F. C. 2018, MNRAS, 475, 5424 [NASA ADS] [CrossRef] [Google Scholar]
  15. Despali, G., Sparre, M., Vegetti, S., et al. 2019, MNRAS, 484, 4563 [NASA ADS] [CrossRef] [Google Scholar]
  16. Despali, G., Lovell, M., Vegetti, S., Crain, R. A., & Oppenheimer, B. D. 2020, MNRAS, 491, 1295 [NASA ADS] [CrossRef] [Google Scholar]
  17. Despali, G., Vegetti, S., White, S. D. M., et al. 2022, MNRAS, 510, 2480 [NASA ADS] [CrossRef] [Google Scholar]
  18. Despali, G., Heinze, F. M., Fassnacht, C. D., et al. 2025, A&A, 699, A222 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Dye, S., Furlanetto, C., Swinbank, A. M., et al. 2015, MNRAS, 452, 2258 [NASA ADS] [CrossRef] [Google Scholar]
  20. Eales, S. A., Raymond, G., Roseboom, I. G., et al. 2010, A&A, 518, L23 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Enzi, W. J. R., Krawczyk, C. M., Ballard, D. J., & Collett, T. E. 2025, MNRAS, 540, 247 [Google Scholar]
  22. Ephremidze, N., Chandrashekar, C., Şengül, A. Ç., & Dvorkin, C. 2025, MNRAS, 542, 2610 [Google Scholar]
  23. Frenk, C. S., & White, S. D. M. 2012, Ann. Phys., 524, 507 [NASA ADS] [CrossRef] [Google Scholar]
  24. Galan, A., de Vyvere, L. V., Gomer, M. R., Vernardos, G., & Sluse, D. 2023, J. Open Source Softw., 8, 5567 [NASA ADS] [CrossRef] [Google Scholar]
  25. Gilman, D., Birrer, S., Nierenberg, A., et al. 2020a, MNRAS, 491, 6077 [NASA ADS] [CrossRef] [Google Scholar]
  26. Gilman, D., Du, X., Benson, A., et al. 2020b, MNRAS, 492, L12 [NASA ADS] [CrossRef] [Google Scholar]
  27. Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357 [NASA ADS] [CrossRef] [Google Scholar]
  28. Heinze, F. M., Despali, G., & Klessen, R. S. 2024, MNRAS, 527, 11996 [Google Scholar]
  29. Hezaveh, Y. D., Dalal, N., Marrone, D. P., et al. 2016, ApJ, 823, 37 [Google Scholar]
  30. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [NASA ADS] [CrossRef] [Google Scholar]
  31. Inoue, K. T., Minezaki, T., Matsushita, S., & Chiba, M. 2016, MNRAS, 457, 2936 [NASA ADS] [CrossRef] [Google Scholar]
  32. Lovell, M. R. 2020, MNRAS, 493, L11 [NASA ADS] [CrossRef] [Google Scholar]
  33. Lovell, M. R., Eke, V., Frenk, C. S., et al. 2012, MNRAS, 420, 2318 [NASA ADS] [CrossRef] [Google Scholar]
  34. Maud, L. T., Asaki, Y., Dent, W. R. F., et al. 2022, ApJS, 259, 10 [Google Scholar]
  35. Maud, L. T., Asaki, Y., Nagai, H., et al. 2023, ApJS, 267, 24 [Google Scholar]
  36. McMullin, J. P., Waters, B., Schiebel, D., Young, W., & Golap, K. 2007, ASP Conf. Ser., 376, 127 [Google Scholar]
  37. Minor, Q. E. 2025, ApJ, 981, 2 [Google Scholar]
  38. Minor, Q., Gad-Nasr, S., Kaplinghat, M., & Vegetti, S. 2021, MNRAS, 507, 1662 [NASA ADS] [CrossRef] [Google Scholar]
  39. Moore, B., Ghigna, S., Governato, F., et al. 1999, ApJ, 524, L19 [Google Scholar]
  40. Nadler, E. O., Yang, D., & Yu, H.-B. 2023, ApJ, 958, L39 [NASA ADS] [CrossRef] [Google Scholar]
  41. Ndiritu, S. W., Vegetti, S., Powell, D. M., & McKean, J. P. 2025, MNRAS, 538, 671 [Google Scholar]
  42. Negrello, M., Hopwood, R., Dye, S., et al. 2014, MNRAS, 440, 1999 [NASA ADS] [CrossRef] [Google Scholar]
  43. Nelson, B. E., Ford, E. B., Buchner, J., et al. 2020, AJ, 159, 73 [Google Scholar]
  44. Nightingale, J. W., & Dye, S. 2015, MNRAS, 452, 2940 [Google Scholar]
  45. Nightingale, J. W., He, Q., Cao, X., et al. 2024, MNRAS, 527, 10480 [Google Scholar]
  46. O’Riordan, C. M., & Vegetti, S. 2024, MNRAS, 528, 1757 [CrossRef] [Google Scholar]
  47. Planck Collaboration VI. 2020, A&A, 641, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  48. Powell, D., Vegetti, S., McKean, J. P., et al. 2021, MNRAS, 501, 515 [Google Scholar]
  49. Powell, D. M., Vegetti, S., McKean, J. P., et al. 2022, MNRAS, 516, 1808 [NASA ADS] [CrossRef] [Google Scholar]
  50. Powell, D. M., Vegetti, S., McKean, J. P., et al. 2023, MNRAS, 524, L84 [NASA ADS] [CrossRef] [Google Scholar]
  51. Ritondale, E., Vegetti, S., Despali, G., et al. 2019, MNRAS, 485, 2179 [Google Scholar]
  52. Rizzo, F., Vegetti, S., Fraternali, F., & Di Teodoro, E. 2018, MNRAS, 481, 5606 [Google Scholar]
  53. Rybak, M., McKean, J. P., Vegetti, S., Andreani, P., & White, S. D. M. 2015a, MNRAS, 451, L40 [Google Scholar]
  54. Rybak, M., Vegetti, S., McKean, J. P., Andreani, P., & White, S. D. M. 2015b, MNRAS, 453, L26 [NASA ADS] [CrossRef] [Google Scholar]
  55. Sawala, T., Frenk, C. S., Fattahi, A., et al. 2016, MNRAS, 456, 85 [Google Scholar]
  56. Sengül, A. Ç., Dvorkin, C., Ostdiek, B., & Tsang, A. 2022, MNRAS, 515, 4391 [CrossRef] [Google Scholar]
  57. Somerville, R. S., & Davé, R. 2015, ARA&A, 53, 51 [Google Scholar]
  58. Spingola, C., McKean, J. P., Auger, M. W., et al. 2018, MNRAS, 478, 4816 [NASA ADS] [CrossRef] [Google Scholar]
  59. Springel, V., Wang, J., Vogelsberger, M., et al. 2008, MNRAS, 391, 1685 [Google Scholar]
  60. Stacey, H. R., Powell, D. M., Vegetti, S., et al. 2024, A&A, 688, A110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  61. Stacey, H. R., Kaasinen, M., O’Riordan, C. M., et al. 2025, A&A, 693, L17 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  62. Suyu, S. H., Marshall, P. J., Hobson, M. P., & Blandford, R. D. 2006, MNRAS, 371, 983 [Google Scholar]
  63. Tajalli, M., Vegetti, S., O’Riordan, C. M., et al. 2025, MNRAS, 543, 540 [Google Scholar]
  64. Tamura, Y., Oguri, M., Iono, D., et al. 2015, PASJ, 67, 72 [Google Scholar]
  65. Vegetti, S., & Koopmans, L. V. E. 2009, MNRAS, 392, 945 [Google Scholar]
  66. Vegetti, S., & Vogelsberger, M. 2014, MNRAS, 442, 3598 [NASA ADS] [CrossRef] [Google Scholar]
  67. Vegetti, S., Koopmans, L. V. E., Bolton, A., Treu, T., & Gavazzi, R. 2010, MNRAS, 408, 1969 [Google Scholar]
  68. Vegetti, S., Lagattuta, D. J., McKean, J. P., et al. 2012, Nature, 481, 341 [NASA ADS] [CrossRef] [Google Scholar]
  69. Vegetti, S., Koopmans, L. V. E., Auger, M. W., Treu, T., & Bolton, A. S. 2014, MNRAS, 442, 2017 [NASA ADS] [CrossRef] [Google Scholar]
  70. Vegetti, S., Despali, G., Lovell, M. R., & Enzi, W. 2018, MNRAS, 481, 3661 [NASA ADS] [CrossRef] [Google Scholar]
  71. Vegetti, S., Birrer, S., Despali, G., et al. 2024, Space Sci. Rev., 220, 58 [CrossRef] [Google Scholar]
  72. Vernardos, G., & Koopmans, L. V. E. 2022, MNRAS, 516, 1347 [CrossRef] [Google Scholar]
  73. Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nat. Meth., 17, 261 [Google Scholar]
  74. Wang, J., Bose, S., Frenk, C. S., et al. 2020, Nature, 585, 39 [Google Scholar]
  75. Wechsler, R. H., & Tinker, J. L. 2018, ARA&A, 56, 435 [NASA ADS] [CrossRef] [Google Scholar]
  76. Wong, K. C., Suyu, S. H., & Matsushita, S. 2015, ApJ, 811, 115 [NASA ADS] [CrossRef] [Google Scholar]
  77. Zhang, G., Mishra-Sharma, S., & Dvorkin, C. 2022, MNRAS, 517, 4317 [Google Scholar]

Appendix A: Additional figure and tables

thumbnail Fig. A.1.

Two-dimensional posterior probability distributions of parameters of four models fit to the band 7 data of SDP 81. Contours show 1σ and 2σ volumes of the posterior.

Table A.1.

Lens model parameters for the band 7 data and their 1σ uncertainties (extended version of Table 2).

Table A.2.

Lens model parameters for the band 7 data and their 1σ uncertainties, using a regular source grid.

Table A.3.

Lens model parameters and their 1σ uncertainties for mock band 7 data with ground truth that is our maximum a posteriori PL+MP model inferred for the real band 7 data.

Table A.4.

Lens model parameters and their 1σ uncertainties for mock band 7 data with ground truth PL+MP+PJ (p,m), where PL+MP is our maximum a posteriori model inferred for the real band 7 data.

All Tables

Table 1.

Summary explanation of the various models tested on the band 7 data.

Table 2.

Lens model parameters for the band 7 data, their 1σ uncertainties, and priors.

Table A.1.

Lens model parameters for the band 7 data and their 1σ uncertainties (extended version of Table 2).

Table A.2.

Lens model parameters for the band 7 data and their 1σ uncertainties, using a regular source grid.

Table A.3.

Lens model parameters and their 1σ uncertainties for mock band 7 data with ground truth that is our maximum a posteriori PL+MP model inferred for the real band 7 data.

Table A.4.

Lens model parameters and their 1σ uncertainties for mock band 7 data with ground truth PL+MP+PJ (p,m), where PL+MP is our maximum a posteriori model inferred for the real band 7 data.

All Figures

thumbnail Fig. 1.

ALMA images of SDP 81 with natural weighting of the visibilities created via forward modelling. The synthesised beam FWHM with natural weighting is 26 × 33 mas and 33 × 42 mas for bands 7 (left) and 6 (right), respectively.

In the text
thumbnail Fig. 2.

Maximum a posteriori lens models for the band 7 data. First column: Lens model. Crosses show lens and sub-halo positions; critical curves are in black. Second column: Difference in convergence relative to the PL model. Contours of the sky model are in grey. Third column: Difference in magnification relative to the PL model; contours of the sky model in black. Fourth column: Image of the residuals (data–model) normalised by the RMS noise (σ), with the mask shown by the dashed outline. Last column: Reconstructed source with Gouraud interpolation and caustics in black. The locations of I16 and H16 purported sub-haloes are shown with a triangle and a circle, respectively.

In the text
thumbnail Fig. 3.

Heat map showing the log Bayes factors of the models tested on the real band 7 data (y-axis model relative to x-axis model). a: Models tested on the real band 7 data. b: Models tested on the real data with H16 multipole coefficients (MP). The smooth models are fixed, and the regularisation is the only free parameter.

In the text
thumbnail Fig. 4.

Heat map showing the log Bayes factors of the models tested on mock band 7 data (y-axis relative to x-axis). a: Models tested on mock band 7 data where the ground truth is PL+MP. b: Models tested on mock band 7 data where the ground truth is PL+MP+PJ (p,m).

In the text
thumbnail Fig. 5.

Heat map showing the log Bayes factors of the models tested on the real band 7 data (y-axis model relative to x-axis model), using a regular source grid of 10 mas. a: Same as Fig. 3a, using a regular source grid. b: Same as Fig. 3b, using a regular source grid.

In the text
thumbnail Fig. 6.

Multipole coefficient mean and 1σ uncertainties for the real and mock data. Yellow and red circles show the posteriors for the mock data with ground truth PL+MP and PL+MP+PJ, respectively; cyan and black triangles show the posteriors for PL+MP and PL+MP+PJ models tested on the real data. The multipole coefficients are consistent within their 1σ uncertainties.

In the text
thumbnail Fig. A.1.

Two-dimensional posterior probability distributions of parameters of four models fit to the band 7 data of SDP 81. Contours show 1σ and 2σ volumes of the posterior.

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.