<figure><img src="image_1.png" /></figure>

# Void Profile from Planck Lensing Potential Map

Teeraparb Chantaval<fnref n="1" />, Utane Sawangwit<fnref n="2" />, and Benjamin D. Wandelt<fnref n="3" />,<fnref n="4" />,<fnref n="5" />

## Abstract

We use the lensing potential map from Planck CMB lensing reconstruction analysis and the “Public Cosmic Void Catalog” to measure the stacked void lensing potential. We have made an attempt to fit the HSW void profile parameters from the stacked lensing potential. In this profile, four parameters are needed to describe the shape of voids with different characteristic radii $R_v$. However, we have found that after reducing the background noise by subtracting the average background, there is a residue lensing power left in the data. The inclusion of the environment shifting parameter, $\gamma_v$, is necessary to get a better fit to the data with the residue lensing power. We divide the voids into two redshift bins: $\text{cmass1}$ ($0.45 < z < 0.5$) and $\text{cmass2}$ ($0.5 < z < 0.6$). Our best-fit parameters are $\alpha = 1.989 \pm 0.149$, $\beta = 12.61 \pm 0.56$, $\delta_c = -0.697 \pm 0.025$, $R_s/R_v = 1.039 \pm 0.030$, $\gamma_v = (-7.034 \pm 0.150) \times 10^{-2}$ for the $\text{cmass1}$ sample with 123 voids and $\alpha = 1.956 \pm 0.165$, $\beta = 12.91 \pm 0.60$, $\delta_c = -0.673 \pm 0.027$, $R_s/R_v = 1.115 \pm 0.032$, $\gamma_v = (-4.512 \pm 0.114) \times 10^{-2}$ for the $\text{cmass2}$ sample with 393 voids at 68% C.L. The addition of the environment shifting parameter is consistent with the conjecture that the Sloan Digital Sky Survey voids reside in an underdense region.

**Key words:** dark matter – gravitational lensing: weak – large-scale structure of universe

## 1. Introduction

In the standard cosmological model, the universe is homogeneous and isotropic on large scales. The seeds of present-day large-scale structure of the universe are formed from the highly Gaussian and nearly scale-invariant power spectrum of matter density (Hinshaw et al. 2013; Planck Collaboration et al. 2015c). However, on small scales, the hierarchical clustering of matter leads to formations of complex cosmic structure such as clusters of galaxies, walls, filaments, and voids (Boylan-Kolchin et al. 2009). Among all of the large-scale objects in the universe, cosmic voids, which are large underdensities in the matter distribution, occupy the vast majority of the universe and hence provide the largest volume-based test on theories of structure formation (Ceccarelli et al. 2006). Being interesting objects in their own right, they contain a wealth of information on the fundamental properties of the universe. For example, the low-density environment of voids is a perfect place to study galaxies, as the galaxies are expected not to be affected by the complex astrophysical processes that modify galaxies in high-density environments, and allows galaxies to evolve independently without environmental effects (Beygu et al. 2013; Penny et al. 2015). In addition, since voids occupy the cosmic volume where the matter density is lowest, the difference between dark energy and modified gravity models for cosmic acceleration could be distinguishable within cosmic voids (Clampitt et al. 2013; Barreira et al. 2015; Cai et al. 2015; Zivick et al. 2015).

The computational approach called N-body dark matter simulations is one of the best tools to empirically understand various void properties such as number functions (Sheth & van de Weygaert 2004; Jennings et al. 2013) and void ellipticity functions (Biswas et al. 2010). However, the definition of voids is rather vague, and various definitions exist in the literature;

some are more suitable for theoretical calculations, while others are more suitable for observations or N-body simulations. This variety of definitions renders comparison of theoretical predictions on void properties and observations difficult. ZOBOV (Zones Bordering On Voidness; Neyrinck 2008) and WVF (Watershed Void Finder; Platen et al. 2007) are two of the popular void-finding algorithms. Both methods are based on some tessellation methods and the watershed concept of defining voids. ZOBOV requires no free parameters or assumptions about the shape and is based on Voronoi tessellation. However, the ZOBOV voids are unsmooth and rather edgy. WVF also requires no free parameters and is based on a watershed transform. However, WVF uses several techniques to smooth the density field so that the WVF voids are not edgy. With void identification algorithms being progressively developed for galaxy redshift surveys such as the Sloan Digital Sky Survey (SDSS), cosmic voids are being continually found, amounting to releases of public void catalogs (Pan et al. 2012; Sutter et al. 2012b, 2014a; Nadathur 2016).

Recently, there has been an increasing amount of attention on voids as objects for various aspects of cosmological studies. The dynamic of voids and redshift-space distortion (Kaiser 1987) is one of the probes of the growth of large-scale structure. The amount of dark matter in the universe could be obtained from the peculiar velocity fields (Courtois et al. 2012). The relationship between their extent angular size and the distance along the line of sight, known as the Alcock-Paczyński test (Alcock & Paczynski 1979), is predicted to be a promising probe of dark energy by using a stacking method to obtain a statistically averaged shape of voids in 2D or 3D spaces. By stacking a large number of voids, one would expect the difference in radial and transverse direction to be directly related to the product of angular distance and the Hubble

<footnote n="1">¹ Laboratory of Cosmology and Gravity, The Institute for Fundamental Study “The Tah Poe Academia Institute”, Naresuan University, Phitsanulok, 65000, Thailand; teeraparbc@nu.ac.th</footnote>

<footnote n="2">² National Astronomical Research Institute of Thailand (NARIT), Chiang Mai, 50200, Thailand</footnote>

<footnote n="3">³ Sorbonne Université, UPMC Univ Paris 6 et CNRS, UMR 7095, Institut d'Astrophysique de Paris, 98 bis bd Arago, 75014 Paris, France</footnote>

<footnote n="4">⁴ Sorbonne Université, Institut Lagrange de Paris (ILP), 98 bis Boulevard Arago, 75014 Paris, France</footnote>

<footnote n="5">⁵ Departments of Physics and Astronomy, University of Illinois at Urbana-Champaign, Urbana, IL 61801, USA

*Received 2016 August 1; revised 2016 December 14; accepted 2016 December 21; published 2017 February 16*</footnote>
parameter (Lavaux & Wandelt 2012; Sutter et al. 2012a, 2014b). Another geometrical study, the evolution of the ellipticity of voids, could be used as a tool in practice to constrain the dark energy equation of state (Lee & Park 2009; Bos et al. 2012). The Integrated Sachs-Wolfe (ISW) effect (Sachs & Wolfe 1967) caused by the evolution of the gravitational potential within voids can also be detected (Granett et al. 2008; Ilić et al. 2013; Cai et al. 2014; Chen & Kantowski 2015; Hotchkiss et al. 2015; Planck Collaboration et al. 2015e).

The measurement of weak gravitational lensing probes the matter distribution by means of the deflection of light from the background sources. The trajectories of photons from background sources are bent toward gravitating matter due to the distortion of spacetime caused by gravitational lensing (Einstein 1936). The scenario is reversed when voids are acting as the sources of gravitational lenses instead of dark matter. The delensing effect of voids has been investigated and recently observed through the distortions of background galaxies by a stacking method that enhances the signal (Higuchi et al. 2013; Krause et al. 2013; Melchior et al. 2014; Clampitt & Jain 2015; Gruen et al. 2016).

The cosmic microwave background (CMB) radiation, which is the signal from the last scattering surface, exists as the ubiquitous background for gravitational lensing. The gravitational anti-lensing effect of voids has been recently investigated by Bolejko et al. (2013), Chen et al. (2015), and Das & Spergel (2009). The CMB signals that are lensed by multiple voids are also a promising tool to obtain good constraints on cosmological parameters (Chantavat et al. 2016). Planck (Planck Collaboration et al. 2015d) has released lensing potential maps from the CMB<fnref n="6" /> that utilized quadratic estimators that exploit the statistical anisotropy induced by lensing (Okamoto & Hu 2003). From a theoretical point of view, if the matter density distribution within a void, known as a void profile, is known, then the lensing potential could be computed and vice versa. The statistical average void density profile known as the universal void profile (hereafter HSW; Hamaus et al. 2014a) has been released and potentially could be exploited to predict the lensing effect of voids at various sizes and redshifts with only a few parameters (more in Section 3). Hence, from the lensing potential data from Planck, one could, in principle, derive the HSW parameters. This would be a good consistency check if the HSW void profile could be reverse-engineered from observables.

The goals of this article are (1) to extract the stacked lensing potential from Planck lensing data by cross-correlation with voids from SDSS data (Sutter et al. 2014b), (2) to compare and cross-examine the derived HSW void parameters from the Planck lensing potential map with other methods, and (3) to better understand the effect of gravitational lensing from voids in the Planck lensing data. We shall begin with a description of the void catalog from the SDSS galaxy redshift survey (Sutter et al. 2014b) and Planck lensing data (Planck Collaboration et al. 2015d) in Section 2. The extraction and cross-correlation methods are also discussed. In Section 3, we describe the HSW void profile parameterization, and a brief overview of the gravitational lensing effect with voids is also introduced. Our parameter estimation is described in Section 4, and the results are shown in Section 5. The discussions and conclusions are

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Lensing potential map, $\psi(\hat{n})$, constructed from Planck lensing convergence, $\kappa_{\ell m}$. The gray shaded area marks the rejected pixels excluded by the analysis mask.</figcaption></figure>

given in Section 6. Throughout this article, our fiducial cosmological parameters are $\Omega_M = 0.315$, $\Omega_\Lambda = 0.685$, $H_0 = 67.3 \text{ km s}^{-1} \text{ Mpc}^{-1}$, $w = -1$, and $\Omega_k = 0$, which is consistent with a flat $\Lambda$CDM cosmology from Planck 2013 + WMAP polarization maximum likelihood cosmological parameters (Planck Collaboration et al. 2014).

## 2. Data and Methodology

This study utilizes the lensing potential measured from the lensed CMB maps to extract the potential in the vicinity around cosmic voids. Here, we describe the data sets and method for measuring the stacked voids' lensing potential, which will be used to constrain the void density profile, as discussed further in Section 3.

### 2.1. Planck Lensing Potential Map

The descriptions of data products and an overview of the scientific results of the Planck full-mission data release are given by Planck Collaboration et al. (2015a). Apart from CMB temperature, polarization frequency maps, and foreground component maps, the science team also released the reconstructed lensing potential map of the CMB (Planck Collaboration et al. 2015d). They applied quadratic lensing estimators and procedures described by Okamoto & Hu (2003) to the foreground-cleaned CMB map. The foreground-cleaned map was constructed from all the frequency band maps using the SMICA procedure (Planck Collaboration et al. 2015b). The contaminated regions of the SMICA map were further removed with the Galaxy, point-source, and SMICA-specific temperature and polarization masks, which leave ≈67% of the sky for the minimum-variance lensing reconstruction analysis. Given the sensitivity and angular resolutions of the Planck mission, the analysis thus results in the most significant measurement of the CMB lensing potential map to date.

The online data provided by the Planck science team are in the standard HEALPix format (Górski et al. 2005). The spherical harmonics coefficients of lensing convergence, $\kappa_{\ell m}$, are given for multipoles up to $\ell_{\max} = 2048$ instead of the lensing potential, $\psi_{\ell m}$. We use the standard definition of lensing convergence to calculate $\psi_{\ell m}$ from

$$ \kappa_{\ell m} = \frac{\ell (\ell + 1)}{2} \psi_{\ell m} \quad (1) $$

We then use the HEALPix synfast package to synthesize the lensing potential map, $\psi(\hat{n})$. The generated map has a resolution of $N_{side} = 2048$ (i.e., ≈1.7 × 1.7 pixels). The

<footnote n="6">⁶ Can be downloaded from http://pla.esac.esa.int/.</footnote>
projected lensing potential map is shown in Figure 1. This reconstructed lensing potential map is used to cross-correlate with the cosmic void catalog, where the azimuthally averaged potential around each void is extracted, scaled, and then combined to constrain the void density profile. For our analysis, we shall use only those regions that pass the analysis mask (Galaxy + point-source + SMICA).

## 2.2. Cosmic Void Catalog

In this work, we use the “Public Cosmic Void Catalog” (Sutter et al. 2012b) constructed using a modified and extended version of the watershed algorithm ZOBOV, called “Void IDentification and Examination” (VIDE; Sutter et al. 2014b). We applied the code to the SDSS Data Release 7 (Abazajian et al. 2009) main galaxy sample and SDSS-III Baryon Oscillation Spectroscopic Survey (BOSS) Data Release 10 (Ahn et al. 2014) LOWZ and CMASS samples. This results in ≈1500 individually detected voids. The samples represent volume-limited catalogs of voids at various redshift bins.

To minimize any possible evolution of the profile parameters, we choose to work with voids in individual redshift bins and do not stack voids across the bins. In order to have as many voids in a single bin as possible, we chose a high-redshift bin to increase the volume. Here, we report on the analysis using the “dr10cmass1” (0.45 < z < 0.5, Δt<sub>age</sub> = 0.4 Gyr, N<sub>void</sub> = 229) and “dr10cmass2” samples (0.5 < z < 0.6, Δt<sub>age</sub> = 0.7 Gyr, N<sub>void</sub> = 696), hereafter called “cmass1” and “cmass2,” respectively.

We therefore stacked the Planck lensing potential map (Section 2.1) around the central positions of cosmic voids taken from the Sutter et al. (2014b) “dr10cmass1” and “dr10cmass2” samples as two separated measurements. However, during our analysis we found that considerable numbers (106 and 303 for cmass1 and cmass2, respectively) of the voids are located in the negative lensing potential regions, which could be due to a number of reasons and warrants further investigations. We therefore restrict our analysis to 123 and 393 voids from samples cmass1 and cmass2, respectively. The radius and redshift distributions of the subsamples that we used do not show any significant difference from those of the excluded low-signal subsamples (see Figure 2).

## 2.3. Stacking Analysis

The stacking analysis is performed on the Planck lensing potential map (See Section 2.1) around each of the voids' centers. Since we need to compare our measurements to the theoretical prediction of the size-independent potential, $\tilde{\psi}(r/R_V)$ (see Equation (11) in Section 3.2), but the potential map is a 2D projection on the surface of a sphere, we therefore bin up the lensing potential according to physical separation and not angular separation. Therefore, for each void, we use the comoving angular diameter distance, $D_A(z)$, to scale the pixels of the lensing maps surrounding the vicinity of each void according to its redshift z and radius $R_V$. For the ith ($r/R_V$) bin, its corresponding angular bin is given by

$$ \theta_i(z, R_V) = R_V / D_A(z) \times (r/R_V)_i. \quad (2) $$

The lensing potential value in each $r/R_V$ bin around a void is azimuthally averaged and is then background subtracted by the measurement at the largest separation, 10 $R_V$. Next, we scale the amplitude of the overall lensing potential measured from each

void according to Equation (12) to account for their different redshifts and radii before we combine them together. Each lensing potential is scaled to its respective median redshift and radius,

$$ \tilde{\psi}(r/R_V)_i = \psi(r; z, R_V)_i \times \frac{R_{V,\text{med}}^3 (1+z_{\text{med}})^3 D_+(z_{\text{med}}) D_A(z)}{R_V^3 (1+z)^3 D_+(z) D_A(z_{\text{med}})}, \quad (3) $$

where $R_{V,med}$ and $z_{med}$ are the median void radius and redshift, respectively, of the bins given in Figure 2. The average lensing potential is then calculated from the scaled measurements of 123 and 393 voids for cmass1 and cmass2 samples, respectively. The uncertainties of our measurement are then calculated using the jackknife resampling technique. The void sample is separated into 12 subsamples according to their Galactic latitudes and longitudes. We then remeasure the lensing potential 12 times, each jth time leaving out one subsample, and the jackknife error is given by

$$ \Sigma_{ij}^{\text{JK}}(r/R_V) = \frac{11}{12} \sum_{i,j}^{12} (\tilde{\psi}_i(r/R_V) - \tilde{\psi}(r/R_V)) \times (\tilde{\psi}_j(r/R_V) - \tilde{\psi}(r/R_V)), \quad (4) $$

where $\tilde{\psi}(r/R_V)$ is the averaged lensing potential from 12 jackknife subsamples. The measured lensing potential and the estimated uncertainties are shown in Figure 3. Using the same jackknife sub-sampling method, we also estimate the correlations between measurements of different bins (off-diagonal elements). The estimated covariance matrices are then used in our fitting procedure as described in Section 4.

# 3. Theory

In this section, we describe all the relevant theories in this analysis, such as the parameterization of the void density profile and the theory of gravitational lensing potential applicable to voids.

## 3.1. Void Density Profile Parameterization

In general, voids will be observed with various shapes and orientations in the field of view. However, the averaged void density profile will be spherically symmetric and is well fitted by the universal void density profile (Hamaus et al. 2014a). The profile is given by

$$ \frac{\delta\rho_V(r)}{\bar{\rho}_M} = \delta_c \frac{1 - (r/R_S)^\alpha}{1 + (r/R_V)^\beta} + \gamma_V, \quad (5) $$

where $\bar{\rho}_M$ is the mean cosmic matter density and $\delta\rho_V$ is the density deviation for the mean density. $R_V$ is the characteristic void radius, and $R_S$ is a scale radius where $\rho_V(R_S) = \bar{\rho}_M \cdot \alpha$, $\beta$, and $\delta_c$ are the shape parameters. $\gamma_V$ is an environment shifting parameter that was not included in the original model. However, the benefit of the inclusion of the environment parameter is twofold. First, the parameter takes into account any systematic uncertainties that may occur in the data extraction process. The parameter $\gamma_V$ is considered as a nuisance parameter, which could be marginalized later. Second, there is a tendency that voids in the SDSS catalogs may reside
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Redshift (top) and radius (bottom) distributions of the void samples used in our analyses from the “dr10cmass1” (black solid lines) and “dr10cmass2” (blue dashed lines) samples, where their medians are $z_{\text{cmass1}} = 0.47$, $R_{V,\text{cmass1}} = 28.26 h^{-1}$ Mpc and $z_{\text{cmass2}} = 0.53$, $R_{V,\text{cmass2}} = 31.68 h^{-1}$ Mpc, respectively. For comparison, the redshift and radius distributions of voids rejected for the stacking analysis, as a result of a very weak lensing signal or having negative lensing potential, are also shown. We do not observe any significant deviation in their redshift and radius distributions.</figcaption></figure>

in an underdense region of the universe (Hamaus et al. 2014b). The constant shifting parameter will take the locality of the environment of voids into account. We shall take the average void profile as our estimate of the void profile in the analysis. Since $R_V$ has been given from the data, we shall use $R_S/R_V$ as one of the fitting parameters. Hence, our vector in the parameter space is $\mathbf{X} = \{\alpha, \beta, \delta_c, R_S/R_V, \gamma_V\}$.

## 3.2. Void Gravitational Lensing Potential

We encourage readers to consult Bartelmann & Schneider (2001) for a general review of gravitational weak lensing. The gravitational potential at redshift $z$ is given by

$$ \nabla^2 \Psi_N = 4\pi G \bar{\rho}_{M0} (1+z) D_+(z) \delta_M(z=0), \quad (6) $$

where $D_+(z)$ is the growth function normalized to unity at $z=0$. $\rho_{M0}$ is the matter density at the present epoch. The lensing potential is defined as the integral over the line-of-sight direction $\hat{n}$,

$$ \psi(\hat{n}) = -\frac{2}{c^2} \int d\chi \nabla_{\perp} \Psi_N(\chi \hat{n}), \quad (7) $$

<figure><img src="image_4.png" /><figcaption>**Figure 3.** Lensing potential from Planck data where voids are detected with SDSS data.</figcaption></figure>

where $\chi$ is the comoving distance. $\nabla_{\perp}$ is the transverse derivative

$$ \nabla_{\perp} \equiv \hat{\theta} \frac{\partial}{\partial \theta} + \frac{\hat{\phi}}{\sin \theta} \frac{\partial}{\partial \phi}. \quad (8) $$

The excess surface density is given by the line-of-sight integral

$$ \delta \tilde{\sigma}_V(b) = \int_{-\infty}^{\infty} dx \delta\rho_V(\sqrt{x^2 + b^2}), \quad (9) $$

where $x = r/R_V$ and $b$ is the scaled impact parameter. It is convenient to define scale-invariant quantities that are a function of the scaled radius $x$. In this article, all the scale-invariant quantities are denoted by the tilde symbol.

From Equation (7), the lensing potential is given by

$$ \tilde{\psi}(b) = \int_{0}^{b} \frac{dy}{y} \int_{0}^{y} dx \, x \, \delta\tilde{\sigma}(x). \quad (10) $$

The lensing potential will be given by

$$ \psi(r; R_V, z) = S(R_V, z) \times \tilde{\psi}(r/R_V), \quad (11) $$

where the scaling factor $S(R_V, z)$ is

$$ S(R_V, z) = \frac{16\pi G}{c^2} \bar{\rho}_{m0} \left( \frac{R_V}{\text{Mpc h}^{-1}} \right)^3 \times \frac{(1+z)^3 D_+(z)}{(D_A(z)/\text{Mpc h}^{-1})}, \quad (12) $$

where $D_A(z)$ is the comoving angular diameter distance.

## 4. Parameter Estimation

In order to find the best-fitting set of parameters for the lensing potential, we shall adopt the maximum likelihood estimator method as a fitting criterion (Hald 1999). Using the
log-likelihood function of the form,

$$ \log \mathcal{L}[\tilde{\psi}; \mathbf{X}] = -\frac{1}{2} \sum_{i,j}^{M} (\tilde{\psi}(x_i; \mathbf{X}) - \tilde{\psi}_i^D)(\Sigma^{\mathrm{JK}})_{ij}^{-1} \\ \times (\tilde{\psi}(x_j; \mathbf{X}) - \tilde{\psi}_j^D), \quad (13) $$

where $\log \mathcal{L}[\tilde{\psi}; \mathbf{X}]$ is the log-likelihood functional of the lensing potential $\tilde{\psi}(x_i; \mathbf{X})$ and $\tilde{\psi}_i^D$ is the lensing potential from the data as in Equation (3). $\Sigma^{\mathrm{JK}}$ is the data covariance matrix in Equation (4). The summation is running over all the data points $\mathbf{M}$. The best-fit parameters will be the set of parameters $\mathbf{X}_{\mathrm{best}}$ that maximizes the functional.

In order to explore the parameter space effectively, a Markov chain Monte Carlo (MCMC) method with the Metropolis-Hastings algorithm is implemented (Metropolis et al. 1953; Hastings 1970). The Metropolis-Hasting algorithm will advance the state with an acceptance probability from a parameter state $\mathbf{X}$ to a parameter state $\mathbf{Y}$ given by

$$ \alpha(\mathbf{X}, \mathbf{Y}) = \min \left\{ 1, \frac{\pi(\mathbf{Y}) q(\mathbf{Y}, \mathbf{X})}{\pi(\mathbf{X}) q(\mathbf{X}, \mathbf{Y})} \right\}, \quad (14) $$

where $q(\mathbf{Y}, \mathbf{X})$ is the proposal probability distribution from $\mathbf{X}$ to $\mathbf{Y}$. We shall take a multivariate normal distribution as our density proposal distribution,

$$ q(\mathbf{Y}, \mathbf{X}) = \frac{1}{((2\pi)^N \det(\Sigma))^{1/2}} \\ \times \exp \left\{ -\frac{1}{2}(\mathbf{Y} - \mathbf{X})^\top \Sigma^{-1} (\mathbf{Y} - \mathbf{X}) \right\}, \\ \equiv \mathcal{N}(\mathbf{Y}; \mathbf{X}, \Sigma), \quad (15) $$

where $\Sigma$ is the covariance matrix. Since our proposal distribution is symmetric, our acceptance probability is $\alpha = \min\{1, \pi(\mathbf{Y})/\pi(\mathbf{X})\}$. $\pi(\mathbf{X})$ is the weighing distribution, which shall be taken as the likelihood function, i.e., the exponential of Equation (13).

The MCMC will sample the parameter space giving a chain, $\{\mathbf{X}_i, i = 0, 1, \dots\}$, for the $i$th iteration. The chain will continue until an equilibrium state is reached, and the MCMC will sample the underlying posterior distribution. For an MCMC with a fixed variance of the proposal distribution $\Sigma$, this could lead to a situation where the acceptance rate is either too small or too large. This could result in a final posterior distribution being localized or a slow convergence rate, respectively. For a better learning performance, an adaptive MCMC algorithm (Andrieu & Thoms 2008) shall be implemented. The algorithm has the advantage of adjusting the variance according to the acceptance probability by introducing the scaling factor $\lambda^j$ for the marginal variance $[\Sigma]_{i,j}$. In addition, the mean value of the posterior distribution

$$ \mu \equiv \langle \mathbf{X} \rangle \quad (16) $$

and the covariance matrix

$$ \Sigma \equiv \langle (\mathbf{X} - \mu)^\top \cdot (\mathbf{X} - \mu) \rangle \quad (17) $$

are updated for each iteration. This is achieved by a comparison of $\alpha(\mathbf{X}, \mathbf{Y})$ with a preferred acceptance rate $\alpha_*$. If $\alpha(\mathbf{X}, \mathbf{Y}) < \alpha_*$ for most transition attempts, then $\Sigma$ should be increased. If, on the other hand, $\alpha(\mathbf{X}, \mathbf{Y}) > \alpha_*$ for most transition attempts, then $\Sigma$ should be decreased. However, the

acceptance rate will be compared with the preferred acceptance rate component-wise to allow the case where $\Sigma$ should decrease in one direction in the parameter space and increase in the other direction. Hence, our proposal distribution will be given by $\mathcal{N}(\mathbf{Y}; \mathbf{X}, \Lambda^{1/2} \Sigma \Lambda^{1/2})$, where $\Lambda$ is a scaling matrix,

$$ \Lambda = \mathrm{diag}(\lambda^1, \dots, \lambda^N). \quad (18) $$

The algorithm is explained in detail in the following, where $\hat{e}_k$ is a unit vector with zeros everywhere except the $k$th component and $\zeta_i$ is a nonincreasing function of $i$, the iteration number:

(a) Initialize $\mathbf{X}_0, \mu_0, \Sigma_0$, and $\lambda_0^1, \dots, \lambda_0^N$ for $i = 0$.

(b) Iterate $i + 1$ using the following procedure:

1. For a given $\mu_i, \Sigma_i$ and $\lambda_i^1, \dots, \lambda_i^N$ sample $\Delta\mathbf{X}_i$ from the distribution $\mathcal{N}(\mathbf{X}_i, \Lambda_i^{1/2} \Sigma_i \Lambda_i^{1/2})$.

2. Propose a new state $\mathbf{Y}_{i+1} = \mathbf{X}_i + \Delta\mathbf{X}_i$ and transverse to the new state with probability $\alpha(\mathbf{X}_i, \mathbf{Y}_{i+1})$; otherwise, $\mathbf{X}_{i+1} = \mathbf{X}_i$.

3. Update the scaling factor for $j = 1, \dots, N$,

$$ \log(\lambda_{i+1}^j) = \log(\lambda_i^j) + \zeta_{i+1} [\alpha(\mathbf{X}_i, \mathbf{X}_i + \Delta\mathbf{X}_i(j)\hat{e}_j) - \alpha_*]. \quad (19) $$

4. Update mean and covariance matrix,

$$ \mu_{i+1} = \mu_i + \zeta_{i+1} (\mathbf{X}_{i+1} - \mu_i) \quad (20) $$

and

$$ \Sigma_{i+1} = \Sigma_i + \zeta_{i+1} [(\mathbf{X}_{i+1} - \mu_i)^\top \times (\mathbf{X}_{i+1} - \mu_i) - \Sigma_i]. \quad (21) $$

(c) Repeat the procedure until an equilibrium state is achieved.

There are no constraints on the functional form of $\zeta_i$ as long as it is nonincreasing (Roberts & Rosenthal 2007). We shall set $\zeta_i$ as

$$ \zeta_i = \frac{1}{10 i^{1/3}}. \quad (22) $$

## 5. Results

We divide our lensing potential data into three data sets: Data-I, where we include the lensing potential from the center to 1 $R_V$, and Data-II and Data-III, where we include the lensing potential to 2 $R_V$ and 3 $R_V$, respectively. Data-I will help in exploring the interior structure of voids, while Data-II and Data-III will get an overall fit to the void profile. With the inclusion of an environment parameter from the HSW void profile, we shall refer to the model without the parameter, $\gamma_V$, as Model A and the model with the environment parameter as Model B. The initial position in the parameter space for each chain will be randomly selected from the parameter range, as shown in Equation (23):

$$ \begin{align} \alpha &\in [0.00, 4.50], & \beta &\in [0.00, 24.00], \nonumber \\ \delta_c &\in [-0.80, -0.20], & R_S/R_V &\in [0.80, 1.60], \nonumber \\ \gamma_V &\in [-0.30, 0.30]. & & & \tag{23} \end{align} $$

The range for each parameter is large enough to ensure that our fitting parameters will fall in agreement with other works (Hamaus et al. 2014a, 2015). In the adaptive MCMC described in Section 4, each data set was run for a total number of 100
<table><caption>**Table 1**
Fitting Parameters for Model A and Model B</caption><thead><tr><th rowspan="3">Parameters</th><th colspan="6">cmass1 Sample</th></tr><tr><th colspan="2">Data-I</th><th colspan="2">Data-II</th><th colspan="2">Data-III</th></tr><tr><th>Model A<br>(without γ<sub>V</sub>)</th><th>Model B<br>(with γ<sub>V</sub>)</th><th>Model A<br>(without γ<sub>V</sub>)</th><th>Model B<br>(with γ<sub>V</sub>)</th><th>Model A<br>(without γ<sub>V</sub>)</th><th>Model B<br>(with γ<sub>V</sub>)</th></tr></thead><tbody><tr><td>α</td><td>2.471 ± 0.127</td><td>2.016 ± 0.105</td><td>2.612 ± 0.116</td><td>1.923 ± 0.130</td><td>2.642 ± 0.107</td><td>1.989 ± 0.149</td></tr><tr><td>β</td><td>12.35 ± 0.62</td><td>11.59 ± 0.45</td><td>13.53 ± 0.49</td><td>11.68 ± 0.48</td><td>13.92 ± 0.51</td><td>12.61 ± 0.56</td></tr><tr><td>δ<sub>c</sub></td><td>-0.739 ± 0.024</td><td>-0.699 ± 0.018</td><td>-0.766 ± 0.022</td><td>-0.689 ± 0.024</td><td>-0.786 ± 0.020</td><td>-0.697 ± 0.025</td></tr><tr><td>R<sub>S</sub>/R<sub>V</sub></td><td>1.148 ± 0.025</td><td>0.969 ± 0.028</td><td>1.199 ± 0.021</td><td>1.006 ± 0.034</td><td>1.236 ± 0.021</td><td>1.039 ± 0.030</td></tr><tr><td>γ<sub>V</sub> (×10<sup>-2</sup>)</td><td>N/A</td><td>-9.679 ± 0.663</td><td>N/A</td><td>-8.505 ± 0.298</td><td>N/A</td><td>-7.034 ± 0.150</td></tr></tbody><thead><tr><th rowspan="3">Parameters</th><th colspan="6">cmass2 Sample</th></tr><tr><th colspan="2">Data-I</th><th colspan="2">Data-II</th><th colspan="2">Data-III</th></tr><tr><th>Model A<br>(without γ<sub>V</sub>)</th><th>Model B<br>(with γ<sub>V</sub>)</th><th>Model A<br>(without γ<sub>V</sub>)</th><th>Model B<br>(with γ<sub>V</sub>)</th><th>Model A<br>(without γ<sub>V</sub>)</th><th>Model B<br>(with γ<sub>V</sub>)</th></tr></thead><tbody><tr><td>α</td><td>2.580 ± 0.129</td><td>1.935 ± 0.119</td><td>2.687 ± 0.116</td><td>1.617 ± 0.126</td><td>2.959 ± 0.122</td><td>1.956 ± 0.165</td></tr><tr><td>β</td><td>12.22 ± 0.52</td><td>11.43 ± 0.49</td><td>13.85 ± 0.50</td><td>12.32 ± 0.52</td><td>13.71 ± 0.53</td><td>12.91 ± 0.60</td></tr><tr><td>δ<sub>c</sub></td><td>-0.726 ± 0.021</td><td>-0.668 ± 0.019</td><td>-0.784 ± 0.022</td><td>-0.657 ± 0.021</td><td>-0.817 ± 0.020</td><td>-0.673 ± 0.027</td></tr><tr><td>R<sub>S</sub>/R<sub>V</sub></td><td>1.172 ± 0.023</td><td>0.991 ± 0.031</td><td>1.260 ± 0.017</td><td>1.002 ± 0.037</td><td>1.276 ± 0.016</td><td>1.115 ± 0.032</td></tr><tr><td>γ<sub>V</sub> (×10<sup>-2</sup>)</td><td>N/A</td><td>-5.716 ± 0.693</td><td>N/A</td><td>-6.747 ± 0.239</td><td>N/A</td><td>-4.512 ± 0.114</td></tr></tbody></table>

Note. All the uncertainties are $1\sigma$ for the different data sets.

chains with 20,000 sampling points per chain. As for the burn-in period, we abandon the initial 15,000 points and take only the last 5000 points for our parameter estimation. For Data-I our fitting parameters are all the void parameters without the environment parameter, i.e., $γ_V = 0$. This is due to the fact that the fit to Data-I with only void parameters is already good without the environment parameter. Adding the environment parameter would produce an unnecessary overfit to the data. However, the void parameters with the environment parameter are used to fit Data-II and Data-III as the merit of an additional parameter to the model overcomes the penalty of adding more parameters (Occam's razor). The interpretation of the environment parameter will be further discussed in Section 6.

We take the average mean and covariance matrix as follows:

$$ \mu_{\text{ave}} = \frac{1}{N_c} \sum_{i=1}^{N_c} \mu_i \quad (24) $$

and

$$ \Sigma_{\text{ave}}^{-1} = \frac{1}{N_c} \sum_{i=1}^{N_c} \Sigma_i^{-1}, \quad (25) $$

where $N_c$ is the number of chains and $\Sigma_i$ is the final covariance matrix for the ith chain. We shall use the Gaumixmod algorithm (Press et al. 2007) to find the best-fit $\mu_i$ and $\Sigma_i$ for each chain. Our best-fitting parameters for Model A and Model B are shown in Table 1. The resulting void profile and lensing potential from the best-fit parameters in each data set are shown in Figure 4. For comparison, the best-fit parameters for Model A and Model B for each data set are shown in Figure 5.

## 6. Discussions and Conclusions

From the void profiles shown in the left panel of Figure 4 and the void parameters in Table 1, if we calculate the excess mass $\delta m = 4\pi \int_{0}^{\infty} \delta V(x) * x^2 dx$, we could see that the excess

mass lies within the range (−2.5, −1.5) for Model A and (−1.5, −0.5) for Model B. The range of value indicates that voids found in SDSS data are mostly undercompensated, having less density than the average. This indicates the fact that voids in underdense regions are easily found with high significant levels in SDSS data since large voids are usually found in underdense regions. We also notice that without $γ_V$ the excess mass is higher in the negative value.

In Table 1, we show the best-fit parameters within $1\sigma$ uncertainties. The parameters from different data sets all agree within the uncertainty ranges. The values of the environment parameter $γ_V$ are ~−0.08 for cmass1 and ~−0.05 for cmass2 with a high significant detection level from zero. The negative value of $γ_V$ indicates that on average the voids in SDSS data are in a locally underdense region, which is also confirmed by the integrated mass discussed previously. This effect could be seen from the lensing potential shown in Figure 3. Since the deflection angle is equal to the gradient of the lensing potential ($\alpha \equiv \nabla_{\perp} \psi$), the nonzero gradient of the lensing potential at large radius indicates that there is some residue lensing power. The lensing power mainly comes from the deficit (or excess) mass from the voids (or clusters). The constancy of the gradient translates into the constancy of the excess mass. Our MCMC has shown that this residue power is caused by the mass deficit by having a negative value of $γ_V$. However, this comes with the caveat, as stated in Planck Collaboration et al. (2015d), that the lensing potential has a very red spectrum and when cutting the map into small regions it could cause a leakage issue. This is the reason the team chose to release the lensing convergence rather than the potential. And for our purposes we need to convert it back to the lensing potential (Section 2.1).

To justify the necessity of the environment parameter, the best-fit lensing potentials from the models with and without the environment parameter are shown in Figure 5. The likelihood ratios between the models with and without the environment parameter $L_w/L_{w/o}$ are 1.01, 1.36, and 77.4 for Data-I, Data-II,
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Best-fit void profiles and lensing potentials for each data set. The blue, red, and green bands are for Data-I, Data-II, and Data-III, respectively. All the error bars are $1\sigma$ for the cmass1 sample.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 5.** Best-fit lensing potential for Data-I, Data-II, and Data-III (left to right). The solid lines are void parameters with the environment parameter, while the dashed lines are void parameters without the environment parameter. Even though all the models with the environment parameter fit better with the data, the model without the environment parameter also fits well for Data-I and Data-II. All the plots are subjected to the condition $\psi(r = 3R_v) = 0$. The error bars are $1\sigma$.</figcaption></figure>

and Data-III, respectively, for cmass1. Without taking the penalty of an additional parameter in $\mathcal{L}_w$ into account, we could see that the environment parameter is preferred in Data-III, while it is less preferable for Data-I and Data-II. The best-fit parameters for all samples and models are shown in Table 1. From Table 1, the significant detection levels for the environment parameter for Data-I, Data-II, and Data-III are approximately 15$\sigma$, 30$\sigma$, and 47$\sigma$ for Data-I, Data-II, and Data-III, respectively, for the cmass1 sample. The general results hold similarly for the cmass2 sample.

To analyze the effect of $\gamma_V$ on the other parameters, we shall take the best-fit parameters from Data-III for the cmass1 sample from Table 1 for a comparison. The effects of how the HSW parameters alter the shape of the void profile are explicitly shown in Figure 8 in Barreira et al. (2015). In both models, the values of $\beta$ between the two models are not much different given the uncertainties in the values—the values differ by 2.3$\sigma$. The similarity in the values of $\beta$ indicates that the extensions of the compensation region are similar. However, the values of $\alpha$,

$\delta_c$, and $R_S/R_V$ are significantly different by 4.4$\sigma$, 3.6$\sigma$, and 6.6$\sigma$, respectively. $\alpha$ describes the slope of the underdense region, $\delta_c$ describes the depth of the void profile, and $R_S/R_V$ is the zero-crossing radius. In general, Model B (with $\gamma_V$) gives a shallower void profile than Model A (without $\gamma_V$) by having smaller values of $\alpha$ and $\delta_c$. This indicates that $\gamma_V$ has a direct degenerate effect with both $\alpha$ and $\delta_c$—lowering the mean density could be compensated by having a shallower profile.

Attempts to recover or fit the HSW void parameters are found in Nadathur et al. (2014, hereafter N14) and Hamaus et al. (2015, hereafter H15). In N14, the stacked void profiles for different radii and redshifts are compared between a mock luminous red galaxy (LRG) catalog from the Jubilee simulation and the SDSS LRG and Main Galaxy samples. They have found that the void profiles from the simulations and the SDSS galaxy samples matched. H15 investigated the redshift-space distortions between pairs of galaxies from a mock galaxy catalog with HOD parameters from Zheng et al. (2007) and Manera et al. (2013) and stacked voids in redshift space. The
inference on the HSW void parameters was made by assuming the Gaussian streaming model, where the distribution of the pairwise line-of-sight velocities is assumed to be Gaussian. The HSW parameters in N14 and H15 are $(\alpha, \beta, \delta_c, R_S/R_v) = (1.57, 5.72, -0.69, 0.81)$ and $(\alpha, \beta, \delta_c) = (0.96 \pm 0.14, 8.84 \pm 1.16, -0.912 \pm 0.052)$, respectively. The main differences in the parameter constraints from our work and theirs are from the different methodology used in deriving the parameters and the inclusion of $\gamma_v$. We also notice that the value of $\delta_c$ from N14 (observationally derived) and our value, especially for Data-I, are similar, while the value of $\delta_c$ from H15 (mock catalog derived) is remarkably different. This may indicate a systematic bias between the SDSS samples and the mock catalog.

To summarize this work, we cross-correlated the Planck lensing map with 516 voids found in the SDSS data and stacked them to obtain the stacked lensing potential from voids. From the stacked void lensing potential, we recover the HSW void parameter from three different data sets: Data-I to $1R_v$, Data-II to $2R_v$, and Data-III to $3R_v$. We have found that it is necessary to include the environment parameter in the void profile to obtain a good fit to the data. The environment parameter has a physical interpretation that voids found in the SDSS data are mostly undercompensated voids and reside within an underdense region. The effects of the deficit mass are shown in Figure 3, where the gradient of the lensing potential is constant at large distances.

We would like to thank Paul Sutter for his useful comments. T.C. acknowledges the support from the National Astronomical Research Institute of Thailand (NARIT). This work is supported by a NARIT research grant and its High Performance Computer Facility. This work has been done within the Labex ILP (reference ANR-10-LABX-63) part of the Idex SUPER, and received financial support managed by the Agence Nationale de la Recherche, as part of the programme Investissements d'avenir under the reference ANR-11-IDEX-0004-02.

## References

* Abazajian, K. N., Adelman-McCarthy, J. K., Agüeros, M. A., et al. 2009, ApJS, 182, 543

* Ahn, C. P., Alexandroff, R., Allende Prieto, C., et al. 2014, ApJS, 211, 17

* Alcock, C., & Paczynski, B. 1979, Natur, 281, 358

* Andrieu, C., & Thoms, J. 2008, Statistics and Computing, 18, 343

* Barreira, A., Cautun, M., Li, B., Baugh, C. M., & Pascoli, S. 2015, JCAP, 8, 028

* Bartelmann, M., & Schneider, P. 2001, PhR, 340, 291

* Beygu, B., Kreckel, K., van de Weygaert, R., van der Hulst, J. M., & van Gorkom, J. H. 2013, AJ, 145, 120

* Biswas, R., Alizadeh, E., & Wandelt, B. D. 2010, PhRvD, 82, 023002

* Bolejko, K., Clarkson, C., Maartens, R., et al. 2013, PhRvL, 110, 021302

* Bos, E. G. P., van de Weygaert, R., Dolag, K., & Pettorino, V. 2012, MNRAS, 426, 440

* Boylan-Kolchin, M., Springel, V., White, S. D. M., Jenkins, A., & Lemson, G. 2009, MNRAS, 398, 1150

* Cai, Y.-C., Neyrinck, M. C., Szapudi, I., Cole, S., & Frenk, C. S. 2014, ApJ, 786, 110

* Cai, Y.-C., Padilla, N., & Li, B. 2015, MNRAS, 451, 1036

* Ceccarelli, L., Padilla, N. D., Valotto, C., & Lambas, D. G. 2006, MNRAS, 373, 1440

* Chantavat, T., Sawangwit, U., Sutter, P. M., & Wandelt, B. D. 2016, PhRvD, 93, 043523

* Chen, B., & Kantowski, R. 2015, PhRvD, 91, 083014

* Chen, B., Kantowski, R., & Dai, X. 2015, ApJ, 804, 130

* Clampitt, J., Cai, Y.-C., & Li, B. 2013, MNRAS, 431, 749

* Clampitt, J., & Jain, B. 2015, MNRAS, 454, 3357

* Courtois, H. M., Hoffman, Y., Tully, R. B., & Gottlöber, S. 2012, ApJ, 744, 43

* Das, S., & Spergel, D. N. 2009, PhRvD, 79, 043007

* Einstein, A. 1936, Sci, 84, 506

* Górski, K. M., Hivon, W., Banday, A. J., et al. 2005, ApJ, 622, 759

* Granett, B. R., Neyrinck, M. C., & Szapudi, I. 2008, ApJL, 683, L99

* Gruen, D., Friedrich, O., Amara, A., et al. 2016, MNRAS, 455, 3367

* Hald, A. 1999, Statist. Sci., 14, 214

* Hamaus, N., Sutter, P. M., Lavaux, G., & Wandelt, B. D. 2015, JCAP, 11, 036

* Hamaus, N., Sutter, P. M., & Wandelt, B. D. 2014a, PhRvL, 112, 251302

* Hamaus, N., Wandelt, B. D., Sutter, P. M., Lavaux, G., & Warren, M. S. 2014b, PhRvL, 112, 041304

* Hastings, W. K. 1970, Biometrika, 57, 97

* Higuchi, Y., Oguri, M., & Hamana, T. 2013, MNRAS, 432, 1021

* Hinshaw, G., Larson, D., Komatsu, E., et al. 2013, ApJS, 208, 19

* Hotchkiss, S., Nadathur, S., Gottlöber, S., et al. 2015, MNRAS, 446, 1321

* Ilić, S., Langer, M., & Douspis, M. 2013, A&A, 556, A51

* Jennings, E., Li, Y., & Hu, W. 2013, MNRAS, 434, 2167

* Kaiser, N. 1987, MNRAS, 227, 1

* Krause, E., Chang, T.-C., Doré, O., & Umetsu, K. 2013, ApJL, 762, L20

* Lavaux, G., & Wandelt, B. D. 2012, ApJ, 754, 109

* Lee, J., & Park, D. 2009, ApJL, 696, L10

* Manera, M., Scoccimarro, R., Percival, W. J., et al. 2013, MNRAS, 428, 1036

* Melchior, P., Sutter, P. M., Sheldon, E. S., Krause, E., & Wandelt, B. D. 2014, MNRAS, 440, 2922

* Metropolis, N., Rosenbluth, A. W., Rosenbluth, M. N., Teller, A. H., & Teller, E. 1953, JChPh, 21, 1087

* Nadathur, S. 2016, MNRAS, 461, 358

* Nadathur, S., Hotchkiss, S., Diego, J. M., et al. 2014, arXiv:1412.8372

* Neyrinck, M. C. 2008, MNRAS, 386, 2101

* Okamoto, T., & Hu, W. 2003, PhRvD, 67, 083002

* Pan, D. C., Vogeley, M. S., Hoyle, F., Choi, Y.-Y., & Park, C. 2012, MNRAS, 421, 926

* Penny, S. J., Brown, M. J. I., Pimbblet, K. A., et al. 2015, MNRAS, 453, 3519

* Planck Collaboration, Adam, R., Ade, P. A. R., et al. 2015a, A&A, 594, A1

* Planck Collaboration, Adam, R., Ade, P. A. R., et al. 2015b, A&A, 594, 9

* Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2014, A&A, 571, A16

* Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2015c, arXiv:1502.01589

* Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2015d, arXiv:1502.01591

* Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2015e, arXiv:1502.01595

* Platen, E., van de Weygaert, R., & Jones, B. J. T. 2007, MNRAS, 380, 551

* Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P. 2007, Numerical Recipes: The Art of Scientific Computing (3rd ed.; New York: Cambridge Univ. Press)

* Roberts, G. O., & Rosenthal, J. S. 2007, J. Appl. Probab., 44, 458

* Sachs, R. K., & Wolfe, A. M. 1967, ApJ, 147, 73

* Sheth, R. K., & van de Weygaert, R. 2004, MNRAS, 350, 517

* Sutter, P. M., Lavaux, G., Wandelt, B. D., et al. 2014a, MNRAS, 442, 3127

* Sutter, P. M., Lavaux, G., Wandelt, B. D., & Weinberg, D. H. 2012a, ApJ, 761, 187

* Sutter, P. M., Lavaux, G., Wandelt, B. D., & Weinberg, D. H. 2012b, ApJ, 761, 44

* Sutter, P. M., Pisani, A., Wandelt, B. D., & Weinberg, D. H. 2014b, MNRAS, 443, 2983

* Zheng, Z., Coil, A. L., & Zehavi, I. 2007, ApJ, 667, 760

* Zivick, P., Sutter, P. M., Wandelt, B. D., Li, B., & Lam, T. Y. 2015, MNRAS, 451, 4215