© ESO 2020

1. Introduction

Galaxy surveys, including the 2dF Galaxy Redshift Survey (Colless et al. 2001) and the Sloan Digital Sky Survey (SDSS; Zehavi et al. 2011), have shown that visible matter in our Universe is not uniformly distributed on intermediate scales ∼100 h⁻¹ Mpc. Instead, a web-like structure is observed with clusters of galaxies identifying the densest regions. N-body simulations predict the existence of these large-scale structures (e.g. Bond et al. 1996; Springel et al. 2005), suggesting a hierarchical structure formation for the cosmic web. We can classify the web (e.g. Eardley et al. 2015) into regions of clusters, filaments, sheets and voids. In this cosmic web, large under-dense regions (voids) are enclosed by anisotropically collapsed surface structures (sheets) and line structures (filaments) which intersect at the most over-dense isotropic regions (clusters). The Zel’dovich approximation predicts that ∼42% of the mass of the Universe is in a filament environment (Zel’dovich 1970), and this has been confirmed by simulations (Aragón-Calvo et al. 2010). However, as filament environments do not display a very high density contrast, this makes direct observations challenging.

One way to observe filaments is from the X-ray emission induced by the warm hot intergalactic medium (WHIM) with several inter-cluster filaments investigated in this way (Briel & Henry 1995; Kull & Böhringer 1999; Werner et al. 2008). There are also reported detections of filaments using overdensities of galaxies (Pimbblet & Drinkwater 2004; Ebeling et al. 2004). Recently, two independent studies (de Graaff et al. 2019; Tanimura et al. 2019) detected the Sunyaev–Zel’dovich (SZ) signal from the ionised gas in the cosmic web. They estimated the density of ionised gas to be ∼(28 ± 12)% of the total baryon density in the Universe, close to resolving the missing baryon problem (Bregman 2007).

In this paper, we investigate the use of weak gravitational lensing to detect filaments. Based on the distortion of light rays around massive objects, gravitational lensing probes the total mass traced by the large-scale structures and is, therefore, highly complementary to the SZ detection of the gas. Although Dietrich et al. (2012) made a direct weak lensing detection of a filament connecting two massive Abell clusters, the direct detection of typical individual filaments is limited by the low signal-to-noise measurement, and most studies stack large samples of candidate filaments and analyse the resulting average weak lensing signal instead (Mead et al. 2010; Clampitt et al. 2016; Epps & Hudson 2017; Kondo et al. 2019).

Clampitt et al. (2016; hereafter C16) determined the weak lensing signal at around 135 000 pairs of SDSS Luminous Red Galaxies (LRGs) with a projected separation 6 h⁻¹ Mpc ≤ R_sep ≤ 14 h⁻¹ Mpc and a redshift separation of Δz <  0.004. Using a “nulling” estimator that cancels the spherically symmetric contribution of the LRG haloes in the shear measurement, they reported a 4.5σ detection of the filament lensing signal. In another study, Epps & Hudson (2017) used ∼23 000 pairs of LRGs from the Baryon Oscillation Spectroscopic Survey (BOSS) LOWZ and CMASS samples as tracers of filaments. Using data from the Canada-France-Hawaii Telescope Lensing Survey (CFHTLenS), they performed a mass reconstruction of a set of stacked LRG pairs with a projected angular separation between 6 and 10 h⁻¹ Mpc and a redshift separation Δz <  0.003. After subtracting the signal from a mass reconstruction of a set of stacked LRG pairs with the same separation on the sky but a greater redshift separation (0.033 <  Δz <  0.04) such that haloes should not be physically connected, they reported a 5σ detection of a filament lensing signal. A more recent study (Kondo et al. 2019) used 70 210 pairs of LRGs from the CMASS sample with a projected separation between 6 and 14 h⁻¹ Mpc and a line-of-sight separation of less than 6 h⁻¹ Mpc. Using the Subaru Hyper Suprime-Cam (HSC) first-year galaxy shape catalogue and adopting the C16 nulling approach they reported 3.9σ detection of a filament signal.

We note that the methodology taken in these previous studies can be understood as a three-point galaxy–galaxy-shear correlation function conditioned on specific intervals of separation between lens galaxies (Schneider & Watts 2005). The full suite of the galaxy-galaxy-galaxy lensing (GGGL) statistics have been applied to the Red-Sequence Cluster Survey (Simon et al. 2008) and CFHTLenS (Simon et al. 2019) to measure the excess mass around galaxy pairs separated by ≲300 h⁻¹ kpc.

In this paper, we present the weak lensing signal measured between 11 706 LOWZ LRG pairs that have a separation of 3 − 5 h⁻¹ Mpc, combining three public weak lensing surveys; the KiDS+VIKING-450 survey (KV450; Hildebrandt et al. 2018; Wright et al. 2019), the Red Cluster Sequence Lensing Survey (RCSLenS; Hildebrandt et al. 2016), and the CFHTLenS (Heymans et al. 2012). We improve the nulling methodology described in C16 to deal with contamination from filament tracers and use a large suite of N-body simulations to validate our pipeline and compare our results. A standard ΛCDM cosmology has been adopted throughout this study to calculate distances with a matter density Ω_m = 0.3, energy density Ω_Λ = 0.7, effective number of neutrino species N_eff = 3.04, baryon density Ω_b = 0.0 and current Hubble constant H₀ = 100 h km s⁻¹ Mpc⁻¹ where h is the Hubble parameter h = 0.7.

This paper is structured as follows. In Sect. 2, we describe the survey data and simulations. Section 3 summarises the weak lensing formalism, the adopted filament model, and methodology. We show our results in Sect. 4 and conclude in Sect. 5. In Appendix A, we present a validation of our nulling technique. In Appendix B, we document the spherical rotation methodology that is required for high declination surveys.

2. Surveys and simulations

2.1. The lensing surveys

The properties of the three lensing surveys, KV450, RCSLenS, and CFHTLenS, are listed in Table 1. They share a similar data processing pipeline, where the shape measurement of galaxies was conducted using the lensfit model fitting code (Miller et al. 2013). This approach convolves the pixelised model Point-Spread-Function with an analytical surface brightness model consisting of bulge and disk components. It uses model fitting to estimate galaxy ellipticitices and with an associated inverse variance weight, w_s. The (reduced) shear (cf. Eq. (4)) is then given by the weighted average of ellipticities, . The observed shear is biased with respect to the true shear and is typically described by the linear bias model (Heymans et al. 2006) as

(1)

In all of these three surveys, the shear multiplicative bias terms were characterised as a function of the signal-to-noise ratio and size of the galaxies, thereby allowing us to calculate the bias for an arbitrary selection of galaxies. The correction for this multiplicative bias was carried out as it was in Velander et al. (2014).

Photometric redshifts, z_B, were estimated using the Bayesian photometric redshift algorithm (BPZ; Benítez 2000) as detailed in Hildebrandt et al. (2012). Wright et al. (2019) and Hildebrandt et al. (2018) show how the photometric redshifts distributions for KV450 are then calibrated using the “weighted direct calibration” method, with weights estimated using a deep spectroscopic training sample in 9-band ugriZYJHK_s magnitude space. No such calibration was performed for the 4-band (RCSLenS) or 5-band (CFHTLenS) surveys, and instead a probability distribution of true redshifts was estimated from the sum of the BPZ redshift probability distributions. This approach has been demonstrated to carry more systematic error (Choi et al. 2016; Hildebrandt et al. 2017). We discuss how we take this redshift uncertainty into account in our final analysis in Sect. 3.4.

2.2. The BOSS survey

We used the Baryon Oscillation Spectroscopic Survey (BOSS) galaxies from the 12th SDSS Data Release (Alam et al. 2015) to define candidate filaments. Among all three lensing surveys RCSLenS has the most SDSS overlap with almost double that of KV450 or CFHTLenS. Once a robust photometric redshift selection has been applied, however, RCSLenS has only 20%/30% the lensing source density in comparison to CFHTLenS/KV450, (see Table 1 for details).

Given the depth of the lensing surveys and the uncertainty in the high redshift tail of the redshift distribution for CFHTLenS and RCSLenS, we chose to limit our analysis to the LOWZ sample, selected based on colour and magnitude, using a redshift cut 0.15 <  z <  0.43 (Ross et al. 2012). We did not consider the higher redshift CMASS sample. The typical virial halo mass of LOWZ galaxies is ∼5.2 × 10¹³ h⁻¹ M_⊙ (Parejko et al. 2013). These haloes have a typical virial radius ∼1 h⁻¹ Mpc.

2.3. Simulations

We used the Scinet Light Cone Simulations (SLICS¹; Harnois-Déraps & van Waerbeke 2015; Harnois-Déraps et al. 2018) to test our methodology. This suite provided us with 819 independent light cones on a 10 × 10 deg² patch of the sky. Each light cone was constructed from the full non-linear evolution of 1536³ particles with m_p  =  2.88 × 10⁹ h⁻¹ M_⊙, within a 505³(h⁻¹ Mpc)³ cube. Particles were then projected onto mass sheets at 18 redshifts between 0 <  z <  3, and subsequently inspected to identify dark matter haloes. For each simulation, 100 deg² light cone mass sheets and haloes were extracted; the former were then ray-traced into lensing shear maps, while the latter were used to generate mock LOWZ galaxies with a halo occupation distribution, optimised so that the clustering of mock galaxies would be consistent with the LOWZ data (see Harnois-Déraps et al. 2018 for details). We find the mock sample variance to be in good agreement with the Jack-knife errors measured directly from the data. Source galaxy positions were drawn at random, with the shear assigned for a range of redshifts and high number density. We next randomly drew from the SLICS mock source galaxy sample so as to match the corresponding number density and redshift distribution of each of the three surveys, KV450, RCSLenS, and CFHTLenS. Intrinsic galaxy shapes were chosen to match the KiDS ellipticity dispersion (Hildebrandt et al. 2017), which is good description of the ellipticity dispersion that is also found in RCSLenS and CFHTLenS.

3. Summary of weak lensing formalism and methodology

In this section, we summarise the weak gravitational lensing theory following the more detailed derivations in Bartelmann & Schneider (2001). Assuming the thin lens approximation, a foreground object at a position θ has a 2D comoving surface mass density Σ(θ). The convergence is then defined as

(2)

where Σ_crit is the comoving critical surface density in a flat Universe given by

(3)

Here χ is the comoving distance and z_l,  z_s are the redshifts of the lens and source, respectively. Since we are interested in the large scale comoving surface density around filaments, this is the appropriate choice for Σ_crit (see Dvornik et al. 2018, for a discussion on the different definitions of Σ_crit).

The deflection potential, ψ(θ), is connected to the convergence, κ, via Poisson’s equation ∇²ψ = 2κ, and the complex shear is related to the second derivatives of the deflection potential via

(4)

where x₁, x₂ are the horizontal and vertical displacements on the projected sky.

For a filament aligned with the x₁ axis if we assume the deflection potential ψ and convergence κ are both invariant along the filament, then Eq. (4) immediately implies that partial derivatives with respect to the x₁-axis will be equal to zero. This leads to the approximation that for filaments, we should expect to measure γ₁ ≈ −κ and γ₂ ≈ 0. Motivated by the simulation results of Colberg et al. (2005) and Mead et al. (2010) considered the power law density profile around filaments and suggested the model for the convergence at a distance r from filament centre, measured perpendicular to the major filament axis (x₁-axis)

(5)

Here κ_c is the amplitude of the convergence at the filament centre (r = 0) and r_c is the half-maximum radius of the density profile.

3.1. Filament candidiate

Colberg et al. (2005) showed that cluster pairs separated by < 5 h⁻¹ Mpc are always connected by filamentary structures. We therefore select luminous red galaxy (LRG) pairs in the LOWZ catalogue with redshift separation δz <  0.002 and a projected separation 3 h⁻¹ Mpc ≤ R_sep ≤ 5 h⁻¹ Mpc as our candidate filaments² which we will refer to as our physical pairs (PP). Non-physical pairs (NP) are defined to have the same projected separation range but with large line-of-sight separations with 0.033 <  δz <  0.04 (corresponding to ∼100 h⁻¹ Mpc). With such a large physical separation we would not expect to detect a filament signal. The NP therefore provide an important null-test for our methodology.

Our candidates differ from the selection made by C16, Epps & Hudson (2017), de Graaff et al. (2019), Tanimura et al. (2019), and Kondo et al. (2019), who focused on separations of 6 − 10 h⁻¹ Mpc. Our choice maximises signal-to-noise, as shown in our analysis of numerical simulations in Sect. 4.1, but for completeness we also present an analysis of 6 − 10 h⁻¹ Mpc filaments in Sect. 5.

3.2. Stacking method

For each lens filament candidate at redshift z_l = z_f = (z_lens 1 + z_lens 2)/2, we measure the Σ_crit-weighted shear, E_f, on a grid (i, j) centred and oriented with the pair of LRGs, where

(6)

and the sum is taken over all sources³, s, with z_B >  z_l + 0.1 and

(7)

The complex ellipticity is the observed ellipticity of the source rotated into the reference frame where the filament lies along the x₁ axis. The Σ_crit weight converts from a shear estimate to an estimate of the surface mass density Σ in Eq. (2). The grid (i, j) has an extent [ − 2R_sep, 2R_sep]×[−2R_sep, 2R_sep] and 129² pixels, and the pair of lens galaxies that define the filament candidate are positioned to lie at the centre of the pixels at ( − 0.5R_sep, 0) and (0.5R_sep, 0).

We also construct a corresponding weight map for each filament candidate,

(8)

where the extra factor of provides optimal signal-to-noise weighting (Velander et al. 2011).

When rotating each filament pair into a common reference frame, we note that at high declination, the tangent plane method used in Epps & Hudson (2017) and the direct cartesian approximation in C16 results in non-uniform grid cells. As the nulling approach requires a flat geometry on the grid, we found that these approximations lead to a biased result. To solve this problem for high-declination patches, we use the spherical rotation method from de Graaff et al. (2019). This process is detailed in Appendix B, and illustrated in Fig. B.1, with the rotated shear map defined in Eq. (B.4).

As we have spectroscopic redshifts for the filaments but only photometric redshifts for the sources, the inverse critical surface mass density is calculated for each survey as

(9)

(10)

where is the probability distribution of the true redshift of the source galaxies that enter the measurement

(11)

For CFHTLenS and RCSLenS we use the per-source p_s(z_s|z_B) provided by each survey, even though this has been shown to introduce biases (Choi et al. 2016), which we account for in Sect. 3.4. For KV450 we use the weighted direct calibration method of Hildebrandt et al. (2017) to determine the source redshift distribution directly for an ensemble of sources. In practice we calculate in Eq. (10) at eight z_l values and interpolate to evaluate at each filament redshift. The effective n(z_s) for each survey, given by

(12)

is shown in Fig. 1, with CFHTLenS and KV450 providing a deeper source redshift than RCSLenS.

Fig. 1. The effective n(z) of all three lensing surveys, as defined in Eq. (12), when using LOWZ galaxies as lenses. Each curve is normalised such that ∫n(z) dz = 1.

We correct the measured shear signal E_f, with the signal measured around “random” filaments. This is now a standard procedure in galaxy–galaxy lensing studies (e.g. Mandelbaum et al. 2005) which removes any hidden systematics and reduces sampling variance noise. We create random filament samples for each survey patch, listed in Table 2, by randomly shifting filaments within the same patch while preserving their redshifts, position angles and number density in the patch. As we do not expect any physical signal from the random catalogue, we subtract any measured “random” signal from the data as follows:

(13)

Here, for a patch with N_p filament candidates is given by

(14)

and

(15)

where E_r and W_r are the weighted shear (Eq. (6)) and weight (Eq. (8)) but measured around a random pair. w_ran is the normalisation weight where N_R is the number of realisations which ensures N_R × N_p exceeds 100 000 in each patch⁴. This ensures that the random signal has low scatter so that we can take the mean as the random correction.

The random-corrected shear map and weight map W_f are then optimally combined over all filament candidates to determine the total weighted shear signal, 𝒯, over the full sample,

(16)

where w_f = w_lens 1 × w_lens 2 is the product of the recommended SDSS completeness weights for the LOWZ galaxies to account for any extra galaxy pairs in the case of fibre collision. We have also applied the multiplicative calibration correction at this stage which is given by Velander et al. (2014) as

(17)

The total weighted shear 𝒯(i, j) is a combination of both the shear contribution from the haloes surrounding the LOWZ LRGs and the contribution from any filament that connects them.

In order to isolate the filament we apply the “nulling” procedure described in Appendix A to get a final measurement of the shear contribution from the filament only, ℱ(r), as a function of the distance, r, from the central filament axis,

(18)

Here 𝒩 is the nulling operator given in Eq. (A.11) and the summation over pixels from an i_min to i_max runs along the filament from −0.438R_sep to 0.438R_sep. This value was found to minimise any residual contribution from the haloes positioned at ±0.5R_sep, that remains after a nulling analysis of the SLICS simulation (see Sect. 4.1). Our nulling operator 𝒩 combines the shear values measured at 8 different positions (including 4 positions from a reflection about the filament axis) which alternatingly rotate around the two haloes. In this way the isotropic contribution from the parent haloes sum to zero (i.e. “null”) and any anisotropic contribution in-between the two haloes can be recovered. In Appendix A, we provide a detailed proof and compare our nulling approach to that adopted in C16. Through tests on a fiducial model, we show that the C16 nulling approach produces a biased result on large scales.

In this derivation we have carried both components of the shear with ℱ = ℱ_γ1 + iℱ_γ2. Given the Mead et al. (2010) filament model, where κ = −γ₁, we expect ℱ_γ2 = 0. We will fit ℱ_γ1 using the two-parameter model in Eq. (5) with the amplitude parameter replaced by ℱ_c which is equivalent to Σ_critκ_c for a single lens–source pair.

3.3. Error estimation from SLICS

In order to estimate the error on the measured signal from observations, we use a large number of independent and representative lensing simulations from SLICS for each of the three surveys. SLICS allows us to correctly account for the sampling variance, which was found to be the dominant source of noise in Kondo et al. (2019). The source sample of galaxies differs for each filament pair owing to our source selection that z_B − z_f >  0.1. For KV450 we can apply this source selection criteria accurately as the SLICS simulations include mock KV450 photometric redshifts that re-produce the scatter, bias and catastrophic outlier populations found in the KV450 data (Harnois-Déraps et al. 2018). For CFHTLenS and RCSLenS, this information is not encoded. We therefore create mocks from SLICS, modelling source samples for four different filament bins with (z_min, z_max) = [0.1, 0.2],[0.2, 0.3],[0.3, 0.4], and [0.4, 0.5], respectively. For each filament bin, we calculate the source galaxy redshift distribution n(z), using Eq. (12) and the effective galaxy number density for sources with z_B − z_max >  0.1. We then populate 500 independent simulations using these distributions, and measure and combine the weighted shear and weight maps for each of the four filament bins. For KV450 we are able to verify that this binned approach is consistent to the unbinned methodology applied to the data using the KV450 SLICS simulations.

The covariance matrix of the signal is reweighted by the number of filament candidates for each survey, n_{fil, survey}, and estimated from the SLICS simulations as

(19)

where is the filament signal (Eq. (18)) averaged over all N_sim = 500 survey-specific SLICS simulations, and is the average number of filament candidates in these simulations. The covariance is then used to calculate the χ² when estimating parameters in the filament model as

(20)

where is our filament model defined in Eq. (5), calculated on a fine grid of parameters (ℱ_c, r_c). In analogy to Eq. (19), we also define a covariance for the ℱ_γ2 component which serves as a systematic null-test for our analysis.

We found that a simple bootstrap error analysis of the data, where the set of maps are resampled with repetitions before stacking, underestimates the true measurement error. This approach misses the sampling variance term which, like Kondo et al. (2019), we find is a significant component to the error for small-area surveys such as KV450 and CFHTLenS.

3.4. Accounting for uncertainty in the redshift distributions

As discussed in Sect. 2.1, the probability distribution from Eq. (11) have not been calibrated for RCSLenS and CFHTLenS. A systematic uncertainty on the resulting n(z_s) is thus expected. In order to take this into account, we use a nuisance parameter δz_s = 0.1 for RCSLenS and δz_s = 0.04 for CFHTLenS that captures the p(z_s) uncertainty determined by Choi et al. (2016). For KV450 we use δz_s = 0.025 following Wright et al. (2019). We shift the by ±δz_s in Eq. (10) to yield two new functions , and repeat the full measurement and error analysis using both the and calibration. We then estimate the likelihood of the data D given the model (ℱ_c, r_c), assuming a prior on the uncertainty δz_s which consists of 3 δ-functions and marginalised over as

(21)

A full marginalisation where many samples are taken at different redshift offsets, spanning the δz range, is unfortunately unfeasible given the complexity of the measurement pipeline.

4. Results

4.1. Filaments in SLICS

We validated our pipeline using the SLICS simulations of mock LOWZ lens galaxies and mock KiDS sources. We analysed both a noise-free catalogue and a catalogue with shape noise and mock photometric redshifts. The result is shown in Fig. 2 where the upper panels show measurements of ℱ_γ1, and the lower panels show measurements of ℱ_γ2. The left and right columns correspond to results from the 3 − 5 h⁻¹ Mpc and 6 − 10 h⁻¹ Mpc filament candidates, respectively.

Fig. 2. Weighted shear measured between LRG pairs in the SLICS simulation of a random 400 deg2 degree survey. Left: results from 3 to 5 h−1 Mpc filament candidates. Upper panel: average surface mass density ℱγ1, and the lower panel shows the measured cross-shear component. On both panels, the blue data shows the result for a KiDS-like survey depth and shape noise, and the orange data points show the measurement for a noise-free simulation with the errorbar given by the error on the mean of all 158 realisations. Right: set of results from 6 to 10 h−1 Mpc filament candidates.

For the noise-free simulations, we find that ℱ_γ2 is consistent with zero for both the 3 − 5 h⁻¹ Mpc and 6 − 10 h⁻¹ Mpc length filaments. This demonstrates that our nulling procedure correctly removes the contribution from the LRG haloes in the analysis. However, for the simulations with shape noise, both the ℱ_γ1 and ℱ_γ2 measurements are consistent with zero, which suggests that even though KV450 is deeper than the KiDS-450 data simulated in SLICS, we should not expect a significant detection from KiDS alone. As the variance in the noise-free simulation reflects the level of sample variance, we also report that, by measuring the noise level from these two sets of simulations, the sample variance is comparable with the shape noise.

Constraining the parameters of Mead et al. (2010) model with the noise-free SLICS results we find and . The and demonstrate that the model is a good fit to the data (ν = 13 degrees of freedom). We find that the surface mass density of the filament is a factor of 2.5 smaller for the 6 − 10 h⁻¹ Mpc filament and will therefore be more challenging to detect using gravitational lensing.

With the noise-free simulations we are able to analyse whether our signal depends on the redshift of the filament. We constrain the amplitude ℱ_c and the scale r_c parameters of the filament model for 4 redshift quantiles of the SLICS LOWZ filament samples using the same background sources, with the result shown in Fig. 3. The choice of parameterisation on the axes is motivated by the filament model equation as well as for visual simplicity. Here we see significant differences between the samples which could be caused by an evolution in the bias of the LOWZ-like galaxy samples in the SLICS mocks, or the evolution of the filament density field. Given the different source redshift distributions of the three lensing surveys (see Fig. 1) which makes the effective redshift of the average lens differ, this result suggests that we should not necessarily expect the results of these surveys to agree perfectly.

Fig. 3. Redshift evolution of the filament signal in the noise-free SLICS simulation. Left panel: ℱγ1 signal measured from filament candidates in SLICS for four redshift quantiles. Right panel: 68% and 95% confidence region of the model parameters (ℱc, rc) from the corresponding signal in the left panel, with the result from all samples combined shown in grey.

4.2. The detection of filaments with KV450, RCSLenS and CFHTLenS

Figure 4 presents our filament shear measurements and constraints on the two parameters of the filament model (Eq. (5)) for physical pairs, our filament candidates (upper panels), and non-physical pairs, our control sample (lower panels). The left panel shows the nulled ℱ(r) shear signal as a function of the distance from the centre of the filament measured in units of h M_⊙ pc⁻². The result is presented for each lensing survey. We also show the measurements from the three-surveys combined using inverse variance weighting. The blue data points are a measurement of ℱ_γ1, whereas the orange data points show the null-test ℱ_γ2. The shaded region corresponds to the statistical noise from our fiducial measurements, and the capped errorbars correspond to the systematic uncertainty captured by the photometric redshift bias nuisance parameter δz (see Sect. 3.4). The right panel shows the 68% and 95% confidence region of parameters ℱ_c and r_c in Mead et al. (2010) filament model. These estimated parameters can be compared to the best-fit parameter from the noise-free SLICS analysis in Sect. 4.1 which is also represented by the cross in the right panel. We note that the noise-free SLICS best-fit is consistent with all three surveys. We also present joint constraints from the combined signal using a block covariance

(22)

Fig. 4. The detection of the cosmic web between neighbouring luminous red galaxies as detected through the weak lensing of background galaxies from different lensing surveys. Left: x-axis is the distance measured perpendicular to the filament axis scaled such that 1 is equivalent to the projected separation between the pair of LRGs. The y-axis is the nulled shear signal where ℱγ1 (blue data points) measures the average surface mass density of the filament, and ℱγ2 is expected to be consistent with zero and hence serves as a null test. Lower right small panel: measurements from the three-surveys combined using inverse variance weighting. We note that this additional panel is purely for illustration, however, as our joint survey-constraints on the filament model are derived from a combination of the surveys on the likelihood-level. Right: estimated parameters in the filament density model Eq. (5) from the stacked signal for all surveys individually and their combination. Upper: results from 3 ∼ 5 h−1 Mpc Physical Pairs, our filament candidates; Lower: results from 3 ∼ 5 h−1 Mpc non-physical Pairs, our control sample.

and extended data vector ℱ_ext = (ℱ_CF, ℱ_RC, ℱ_KV)^T, providing an estimate of the average filament profile from all the filament candidates across the three surveys. Equation (22) assumes that the surveys are uncorrelated, which is a good approximation to make given the lack of overlap between the different surveys.

To quantify the significance of our measurements, we use the likelihood ratio test between the null hypothesis H₀ and the filament model H₁, where the likelihood ratio LR is

(23)

where B₀ and B₁ are the parameter space in each hypothesis, i.e., B₀ has no free parameters and B₁ = {ℱ_c, r_c}. By Wilks’ theorem (Wilks 1938; Williams 2001), the deviance defined as Dev = 2lnLR has an asymptotic chi-squared distribution with dim(B₁)−dim(B₀) = 2 degrees of freedom when H₀ is true. Estimating the maximum likelihoods from and using Eq. (21) and computing the deviance, we report the significance level for each individual survey as well as the combined analysis in Table 3. The reported suggests our model is a reasonable fit to the data in all cases even for KV450, where p(χ² <  7.40|ν = 13) = 0.12. We find the best-fit model parameters for ℱ_γ1 from all three surveys combined as ℱ_c = 10.5 ± 2.9 h M_⊙ pc⁻² and . We note that, the majority of the detection derives from CFHTLenS alone with a 3.1σ detection. Combining all three surveys we measure a 3.4σ detection of the filament weak lensing signal. CFHTLenS is the most constraining survey as it combines both survey depth with significant SDSS overlap. KV450, with roughly half the source density of CFHTLenS, and RCSLenS, with roughly 20% of the source density, will only start to add significant constraining power with the inclusion of additional overlapping SDSS area.

Our control sample of “non-physical pairs” (NP) were selected to be pairs of lens galaxies with projected separations 3 − 5 h⁻¹ Mpc, but distant in redshift space (0.033 <  Δz <  0.04). These non-physical pairs would not be connected by a filament, hence providing an important validation of our nulling approach to isolate the filament signal. We find that the measured signal for ℱ_γ1, 2 is consistent with zero for all surveys and the combined survey as shown in Table 3 and the lower panel of Fig. 4.

For consistency with other analyses in the literature we also analysed 27 880 filament pairs in LOWZ that have a physical separation of 6 − 10 h⁻¹ Mpc. In contrast to other studies, we do not detect a significant signal for these larger-separation filaments in any of the surveys individually. In combination we find a weak signal at 1.6 σ significance, with ℱ_c = 1.3 ± 0.6. The half-maximum radius of the density profile, r_c, is, however, unconstrained.

5. Conclusions and discussion

In this paper, we present a 3.4σ detection of filamentary structure connecting luminous red galaxies separated by 3 − 5 h⁻¹ Mpc. Through a series of null tests we have verified the robustness of this result. Our work extends the analysis presented in C16 by improving the methodology to null the weak lensing signal from the LRGs in order to isolate the weak lensing distortions induced by the filamentary structure alone. We note that this nulling method cannot distinguish between a pair of spherical haloes joined by a cylinder of matter, or two elliptical haloes which extend towards each other. Higuchi et al. (2014) shows that there is no hard line between a filament and its corresponding halo, with haloes typically extending along the filament. As we find a strong nulled signal on scales much larger than the typical (∼1 h⁻¹ Mpc) virial radius of the LRGs (Parejko et al. 2013), we would argue that it is unlikely to originate from two perfectly aligned haloes (see Xia et al. 2017). But nevertheless we prefer to refer to our detection as that of filamentary structure, rather than that of a filament per se.

Previous studies have focused on LRGs separated by 6 − 10 h⁻¹ Mpc. Using all three lensing surveys (KV450, CFHTLenS and RCSLenS), we do not detect a significant signal for these larger-separation filaments in either the surveys individually, or in combination. This is in contrast to the significant 5σ weak lensing detection of 6 − 10 h⁻¹ Mpc separation filaments reported by Epps & Hudson (2017) using the same CFHTLenS dataset that has been analysed in this study. We report that we are unable to reproduce their result, even when adopting the same methodology. In comparison to C16, we find no comparable detection to their reported ∼4σ detection. As we have shown that sampling variance makes a significant contribution to the overall error budget, the increase by a factor of 5 in the number of filaments studied by C16, in contrast to this analysis, is key to their detection, even though the C16 lensing source galaxy density is significantly shallower than the source densities of KV450, CFHTLenS and RCSLenS. Neither Epps & Hudson (2017) nor C16 constrain the parameters of the Mead et al. (2010) filament model, but in the case of C16 we can compare the amplitudes of the measured signals. When adopting the C16 nulling approach, we find the two shear measurements to be fully consistent. The mean amplitude of our measurement is about four times larger than the amplitude reported in Kondo et al. (2019) (see discussion in Appendix A). Given our error budget, however, the results are consistent.

Authors de Graaff et al. (2019) calculate the average density between galaxy pairs separated by 6 − 14 h⁻¹ Mpc using a CMB lensing convergence map, finding . This estimate assumes that the matter density follows a cylindrical filament model, with density

(24)

where ℓ defines the size of the filament in the line-of-sight direction, r_⊥ defines the distance perpendicular to the filament axis on the projected sky and σ is the intrinsic width. Integrating this density model over the line-of-sight, we can relate this to the surface mass density at the centre of the filament ℱ_c as

(25)

For the value σ = 1.5 h⁻¹ Mpc adopted by de Graaff et al. (2019), and which is our best-fit amplitude parameter ℱ_c, we find for the 3 − 5 h⁻¹ Mpc filament sample at the average filament redshift z = 0.299. For the 6 − 10 h⁻¹ Mpc filament sample we find , which is consistent with the de Graaff et al. (2019) result. Adapting the Mead et al. (2010) filament model, we can also integrate the model over the perpendicular distance and calculate the total mass enclosed between the two LRGs as

(26)

Taking the best-fit parameters ℱ_c and r_c for 3 − 5 h⁻¹ Mpc measurement, we find M_fil = 4.9 ± 2.0 × 10¹³ h⁻¹ M_⊙. It is worth noting that this estimate is based on the approximation that the deflection potential vanishes along the filament major axis. A more detailed analysis would attempt to obtain the excess mass map under the framework of galaxy-galaxy-galaxy lensing (Simon et al. 2008, 2019) but with a much larger separation of galaxy pairs.

Looking forward to upcoming deep weak lensing surveys, such as the European Space Agency’s Euclid mission⁵ and the Large Synoptic Survey Telescope (LSST⁶), and deep spectroscopic surveys such as the Dark Energy Spectroscopic Instrument (DESI⁷), the methodology that we present in this paper could be used to probe filamentary structure as a function of LRG mass and redshift. The combination of overlapping weak lensing surveys and spectroscopic surveys will provide the optimal datasets with which to fully explore the cosmic web.

Acknowledgments

We acknowledge support from the European Research Council under grant numbers 647112 (QX, CH, AA, MA, JHD, TT) and 770935 (HH, AW). CH also acknowledges the support from the Max Planck Society and the Alexander von Humboldt Foundation in the framework of the Max Planck-Humboldt Research Award endowed by the Federal Ministry of Education and Research. YC acknowledges the support of the Royal Society through the award of a University Research Fellowship and an Enhancement Award. We acknowledge support from the European Commission under a Marie-Sklodwoska-Curie European Fellowship under project numbers 656869 (JHD) and 797794 (TT). HH is supported by a Heisenberg grant of the Deutsche Forschungsgemeinschaft (Hi 1495/5-1). AK acknowledges support from Vici grant 639.043.512, financed by the Netherlands Organisation for Scientific Research (NWO). PS acknowledges support from the Deutsche Forschungsgemeinschaft in the framework of the TR33 “The Dark Universe”. Computations for the N-body simulations were performed in part on the Orcinus supercomputer at the WestGrid HPC consortium (www.westgrid.ca), in part on the GPC supercomputer at the SciNet HPC Consortium. SciNet is funded by: the Canada Foundation for Innovation under the auspices of Compute Canada; the Government of Ontario; Ontario Research Fund – Research Excellence; and the University of Toronto. Funding for SDSS-III has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, and the U.S. Department of Energy Office of Science. The SDSS-III web site is http://www.sdss3.org/. We thank SDSS-III for making their data products so easily accessible. SDSS-III is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS-III Collaboration including the University of Arizona, the Brazilian Participation Group, Brookhaven National Laboratory, Carnegie Mellon University, University of Florida, the French Participation Group, the German Participation Group, Harvard University, the Instituto de Astrofisica de Canarias, the Michigan State/Notre Dame/JINA Participation Group, Johns Hopkins University, Lawrence Berkeley National Laboratory, Max Planck Institute for Astrophysics, Max Planck Institute for Extraterrestrial Physics, New Mexico State University, New York University, Ohio State University, Pennsylvania State University, University of Portsmouth, Princeton University, the Spanish Participation Group, University of Tokyo, University of Utah, Vanderbilt University, University of Virginia, University of Washington, and Yale University. Author contributions: All authors contributed to the development and writing of this paper. The authorship list is given in two groups: the lead authors (QX, NR, CH), followed by an alphabetical group who contributed to either the scientific analysis or the data products.

References

Appendix A: Remarks on the nulling technique

In this appendix, we review the C16 nulling technique and develop an improved methodology to isolate the weak lensing signal from filaments. To explain the motivation behind nulling, we start with a single circularly symmetric halo positioned at the origin (0,0), for which the complex shear is given by

(A.1)

where is the mean convergence inside r. We define its counterpart γ^c as the complex shear at the same radial position, with a clockwise rotation of 90°, such that

(A.2)

The counterpart is therefore able to “null” the shear, as γ^c(r)+γ(r) = 0.

For a two-halo system as shown in Fig. A.1, the shear at each position (r, θ) is composed with the shear from halo h₁ and the shear from halo h₂. We can write this as

(A.3)

Fig. A.1. The bipolar configuration for two haloes centred symmetrically about origin O at h1 and h2.

where the co-ordinates (t, φ) are defined with halo h₁ at the origin, and the co-ordinates (s, ϕ) are defined with halo h₂ at the origin. Starting from position P₀ = (r₀, θ₀), shear is given by

(A.4)

A clockwise rotation around halo h₁ by 90° takes us to position P₁ at (r₁, θ₁). The shear here is given by

(A.5)

At position P₁, we see the shear contribution from halo h₁ is the counterpart to the shear contribution from halo h₁ at position P₀. We next rotate around halo h₂ by 90° to position P₂ = (r₂, θ₂), where the shear is given by

(A.6)

Similarly another 90° rotation about halo h₁ (see Fig. A.2) to position P₃ gives

(A.7)

Fig. A.2. Illustration of the nulling technique. Every two line segments with the same colour represents an anti-clockwise rotation with respect to one halo. Right: configuration described in C16. Left: our adopted “X19” configuration which starts below the horizontal axis.

We note that, after another 90° rotation about halo h₂, we come to position P₄ which is our starting point P₀. The sum of the shear from position P₀, P₁, P₂ and P₃ is given by

(A.8)

(A.9)

If we now add in a filament shear γ_f such that at each position γ = γ_f + γ_h1 + γ_h2, then

(A.10)

As we expect the filament shear profile to be symmetric about the horizontal axis, we also sum over the shear from positions and P₃′ that are reflections of P₀, P₁, P₂ and P₃ about the horizontal axis respectively, in order to get the average shear value at any distance away from filament axis. We therefore define the nulling operator 𝒩 as

(A.11)

It is interesting to note that the above equations also apply when two halos are of different masses given their circular symmetry.

In Fig. A.2 we show two configurations for nulling. C16 chose to start P₀ above the horizontal axis and our work starts P₀ below the horizontal axis. We investigate the difference between these two approaches by first constructing a noiseless shear map from two NFW halo profiles using Eq. (14) in Wright & Brainerd (2000) with M_vir = 10¹³ h⁻¹ M_⊙ and assuming a mass–concentration relation from Macciò et al. (2007). We assume both haloes are located at z = 0.3 with background sources at redshift 0.7, close to the mean value of KiDS. The resultant γ₁ map is shown in the top-left panel of Fig. A.3.

Fig. A.3. Left panel: shear maps from two NFW halos with or without fiducial filaments. Top row: γ1 map generated by two NFW haloes only, while the middle row shows the γ1 map with the addition of a fiducial filament profile. Lowest row: γ2 map generated by two NFW haloes. it Right panel: results of the nulling procedure corresponding to the shear map on the left. In the top row, we see under both the X19 configuration (solid line) and C16’s (dashed line), the resulting signal is consistent with zero. When adding a fiducial filament profile, in the middle row, we see our X19 configuration correctly recovers the input value whereas the C16 configuration is biased on large scales. In the lower row, we verify that both configurations null the γ2 signal from two NFW haloes.

For each pixel (x₁, x₂) on the map, we calculate 𝒩[γ_i(x₁, x₂)], and show in the top-right panel in Fig. A.3, the average of the sum along the horizontal axis, i.e., . In the lower-right panel in Fig. A.3, we show the resulting 𝒩[γ₂(x₁, x₂)]. As expected, the sum of nulling pixels are zeros (note the 10⁻¹⁶/10⁻¹⁸ on the y-axis). In the middle panel we repeat the analysis with the inclusion of a fiducial filament modelled using the power-law profile model in Eq. (5). The filament contribution to the shear is included in the 2D shear map as a power-law, symmetric about the x₁-axis (at x₂ = 0) with k_c = 0.02 and r_c = 0.4 h⁻¹ Mpc. These parameters were chosen to roughly replicate previous studies (Dolag et al. 2006), in order to model the true signal contrast in observations. If the nulling method is precise, we would expect the nulled filament signal measured from the map, to recover the input shear from the filament, , irrespective of the nulling method used. The middle-right panel shows that this is indeed the case, with our nulling method. The nulling method in C16, however, recovers a lower expectation value for the shear induced by the filament. This difference results from the C16 configuration. As shown in Eq. (A.10), for the C16 configuration, since P₀ and P₁ are both in the filament region, the nulling operator effectively mixes scales. This has an effect in producing a significant signal on large scales which does not reflect the underlying filament density profile. As P₁, P₂ and P₃ in our adopted configuration lies outside the bridge between two haloes, γ_f(r₁, θ₁)+γ_f(r₂, θ₂)+γ_f(r₃, θ₃) is negligible. This enables us to recover the density profile accurately. We note that Kondo et al. (2019) adopted the C16 estimator but in their Eq. (9) (the equivalent of our Eq. (A.11)), they included an additional factor of 4 in the denominator, which, in our test-case in Fig. A.3 would result in the underestimation of the filament signal by a factor of 4.

Appendix B: Remarks on the spherical rotation

Here we detail the spherical rotation used in Sect. 3.2 to project all filaments onto the same reference frame. As illustrated in Fig. B.1, we rotate all galaxies about a given axis so that the filament (the shorter arc connected by solid pink diamonds) is transformed to lie horizontally on the equator (hollow pink diamonds). To do this, we first transfer the right ascension and declination onto a 3D vector on a unit sphere, such that their positions are g₁ and g₂. The normal vector is defined as . Noting the rotation axis lies on the equator and is perpendicular to the normal vector, we can write down the rotation axis and angle β using components of , so that

(B.1)

(B.2)

Fig. B.1. Illustration of the spherical rotation for filament pairs (solid pink diamonds). Open symbols are the rotated positions.

We note that by defining g₁ to be always on the left of g₂, there is no ambiguity in the definition of the rotation axis, and the rotation angle will always lie between . Rodrigues’ rotation formula (Cheng & Gupta 1989) then allows us to rotate every point on the sky to the desired frame where the filament pair now lies on the equator, such that for each source galaxy at position g_s, the new position is at:

(B.3)

It is worth noting that because the shear was measured in each galaxy’s local (RA, Dec) coordinate frame, the angle of rotation is different for different source galaxies. The transformed shear map () for each galaxy is thus given by

(B.4)

where for each source galaxy, the rotation angle ϕ_s is defined by the angular change of its local coordinate frame, for example, the declination to the y-axis of the filament.

All Tables

All Figures

	Fig. 1. The effective n(z) of all three lensing surveys, as defined in Eq. (12), when using LOWZ galaxies as lenses. Each curve is normalised such that ∫n(z) dz = 1.
In the text

Fig. 2. Weighted shear measured between LRG pairs in the SLICS simulation of a random 400 deg2 degree survey. Left: results from 3 to 5 h−1 Mpc filament candidates. Upper panel: average surface mass density ℱγ1, and the lower panel shows the measured cross-shear component. On both panels, the blue data shows the result for a KiDS-like survey depth and shape noise, and the orange data points show the measurement for a noise-free simulation with the errorbar given by the error on the mean of all 158 realisations. Right: set of results from 6 to 10 h−1 Mpc filament candidates.

In the text

	Fig. 3. Redshift evolution of the filament signal in the noise-free SLICS simulation. Left panel: ℱγ1 signal measured from filament candidates in SLICS for four redshift quantiles. Right panel: 68% and 95% confidence region of the model parameters (ℱc, rc) from the corresponding signal in the left panel, with the result from all samples combined shown in grey.
In the text

Fig. 4. The detection of the cosmic web between neighbouring luminous red galaxies as detected through the weak lensing of background galaxies from different lensing surveys. Left: x-axis is the distance measured perpendicular to the filament axis scaled such that 1 is equivalent to the projected separation between the pair of LRGs. The y-axis is the nulled shear signal where ℱγ1 (blue data points) measures the average surface mass density of the filament, and ℱγ2 is expected to be consistent with zero and hence serves as a null test. Lower right small panel: measurements from the three-surveys combined using inverse variance weighting. We note that this additional panel is purely for illustration, however, as our joint survey-constraints on the filament model are derived from a combination of the surveys on the likelihood-level. Right: estimated parameters in the filament density model Eq. (5) from the stacked signal for all surveys individually and their combination. Upper: results from 3 ∼ 5 h−1 Mpc Physical Pairs, our filament candidates; Lower: results from 3 ∼ 5 h−1 Mpc non-physical Pairs, our control sample.

In the text

	Fig. A.1. The bipolar configuration for two haloes centred symmetrically about origin O at h1 and h2.
In the text

	Fig. A.2. Illustration of the nulling technique. Every two line segments with the same colour represents an anti-clockwise rotation with respect to one halo. Right: configuration described in C16. Left: our adopted “X19” configuration which starts below the horizontal axis.
In the text

Fig. A.3. Left panel: shear maps from two NFW halos with or without fiducial filaments. Top row: γ1 map generated by two NFW haloes only, while the middle row shows the γ1 map with the addition of a fiducial filament profile. Lowest row: γ2 map generated by two NFW haloes. it Right panel: results of the nulling procedure corresponding to the shear map on the left. In the top row, we see under both the X19 configuration (solid line) and C16’s (dashed line), the resulting signal is consistent with zero. When adding a fiducial filament profile, in the middle row, we see our X19 configuration correctly recovers the input value whereas the C16 configuration is biased on large scales. In the lower row, we verify that both configurations null the γ2 signal from two NFW haloes.

In the text

	Fig. B.1. Illustration of the spherical rotation for filament pairs (solid pink diamonds). Open symbols are the rotated positions.
In the text