Towards solving model bias in cosmic shear
forward modeling

Benjamin Remy∗, François Lanusse and Jean-Luc Starck
Université Paris-Saclay, Université Paris Cité, CEA, CNRS, AIM, 91191, Gif-sur-Yvette, France

Abstract
As the volume and quality of modern galaxy surveys increase, so does the difficulty of measuring the cosmological signal imprinted in galaxy shapes. Weakgravitational lensing sourced by the most massive structures in the Universe generates a slight shearing of galaxy morphologies called cosmic shear, key probe forcosmological models. Modern techniques of shear estimation based on statistics
of ellipticity measurements suffer from the fact that the ellipticity is not a welldefined quantity for arbitrary galaxy light profiles, biasing the shear estimation.
We show that a hybrid physical and deep learning Hierarchical Bayesian Model,
where a generative model captures the galaxy morphology, enables us to recover
an unbiased estimate of the shear on realistic galaxies, thus solving the model bias.

1

Introduction

Over the next couple of years, a new generation of wide-field optical galaxy surveys will be coming
online, starting with the ESA Euclid mission [1] and the Vera C. Rubin Observatory Legacy Survey ofSpace and Time [2], soon followed by the Roman Space Telescope [3]. At the core of the cosmology
science that these surveys will enable is the measurement of the weak gravitational lensing effect, avery subtle shearing of the apparent images of distant galaxies in the presence of massive structures
along the line of sight. On cosmological scales, this effect, also known as cosmic shear, is aninvaluable probe of the total matter content of the Universe, and will play a key role in further
constraining and testing our cosmological model over the next decade.
Despite over 20 years of research in this field, unbiased measurements of gravitational shear fromobserved galaxy images remains one of the most challenging aspects of the scientific exploitation ofthese surveys. The prevailing strategy for measuring shear from data has been to build estimatorsof the ellipticity of galaxy light profiles, due to the fact that one can in theory relate how a galaxy
ellipticity e changes under a gravitational shear γ, as eobs ≃ ei + γ. Under the assumption thatgalaxies should have no preferential intrinsic ellipticity, a simple estimator for γ can therefore beconstructed from γ̂ = ⟨eobs ⟩ = ⟨ei ⟩ + γ. However, to this day, no shear measurement algorithm has
proven to be unbiased.
In fact, state-of-the-art approaches for shear measurement [4, 5] no longer attempt to build anunbiased shear estimator from images, but rather calibrate the shape measurement methods estimating
a multiplicative and an additive bias. There are many factors contributing to these biases, but afundamental issue comes from the fact that the ellipticity e is not a well defined quantity for arbitrary
galaxy light profiles.
In this paper we propose a new paradigm for cosmic shear estimation where we no longer rely on
ellipticities anymore. We propose to forward model the galaxy images at the pixel level, similarly to[6], but using non-parametric galaxy light profiles. The main idea is to train a deep generative model∗

Corresponding author: benjamin.remy@cea.fr

Machine Learning and the Physical Sciences workshop, NeurIPS 2022.

to learn galaxy morophologies, which can then be used as part of a hierarchical Bayesian model, ahybrid deep learning and physical model, simulating the observations and used to infer the shear.
The purpose of this paper is first to demonstrate that using parametric light profiles as part of theforward model induces a bias in the shear estimation, which we call a model bias, and then to showthat this bias can be removed by using a deep generative model for the morphology. We show in thispreliminary study, that using a hybrid deep learning-physical Bayesian forward model can not only
solve the model bias, but also provide uncertainty quantification to the shear estimation.

(a) Realistic galaxies fitted by a parametric profile(b) Realistic galaxies fitted by a deep morphology model

Figure 1: The first row of images correspond to observed galaxies. The second row corresponds tothe light profiles Maximum a Posteriori fits, convolved with their PSFs. The third row shows residual
images between the observations and the model fits. Column (a) is using Sérsic model as light profileand column (b) is using deepgal model as light profile.

2

Forward modeling weak lensing observations at the pixel level

From a Bayesian perspective, making the joint inference of cosmic shear and galaxy morphologies
means estimating the joint posterior distribution p(G, γ|D), where G = {zi }i=1...N is the set ofmorphology parameters and D = {di }i=1...N is the set of observed data. The Bayes identity tells us:p(G, γ|D) ∝ p(D|G, γ)p(G, γ),(1)where p(D|G, γ) is the data likelihood given by the forward model described bellow. Since theintrinsic morphology of galaxies are independent of the cosmic shear, then p(G, γ) = p(G)p(γ),
where p(G) and p(γ) are respectively the prior distributions of the morphological parameters and of
the cosmic shear.
We show the probabilistic graphical model in Figure 2. Nodes denote random variables and arrows
conditional relation between them. White nodes corresponds to random variable on which we dothe inference on, conditioned on the shaded nodes corresponding to observed variables. Dot nodes
correspond to fixed variables, e.g. neural network weights θ. In this paper we compare two forward
models: (1) using Sérsic law and (2) using the deepgal model to generate the galaxy light profile.
Once the galaxy light profiles are drawn, we shear images with a constant value γ shared by all thegalaxies. Each sheared galaxy i is then convolved by its Point Spread Function (PSF) Πi . Assuming
Gaussian noise corruption, the resulting image is therefore evaluated by a Gaussian likelhood of
2variance σpix.
Sérsic model The light profile of the Sérsic model isgiven by the following parametrization: I(r) = F ×
1
I0 exp(−bn [( rrhlr ) n − 1]), where rhlr is the half-lightradius, I0 is a normalization factor such that the integrationof I over all the space is F , which corresponds to the fluxof the galaxy. n is called the Sérsic index and bn is adimensionless constant depending on n. The galaxy isthen sheared with parameters e = (e1 , e2 ) corresponding
to its ellipticity. Finally, since the real galaxy may not
be perfectly centered in the postage stamp, we add two
shift parameters sx and sy . Therefore, each galaxy i hasa set of parameter zi = {n, rhlr , F, e1 , e2 , sx , sy }i . Prior
2

θ

zi

di

Πi

γ

σpix

galaxies i = 1 . . . N

Figure 2: Probabilistic Graphical Model.

distributions and parameters follow empirical distributions of the COSMOS shape catalog and aredescribed in subsection A.1.
Deepgal model deepgal is a VAE-based model, proposed in [7], generating complex galaxy light
profiles. In this model, the observations are generated following a random process mapping a latent
space representation zi ∈ R16 to the observed image di following the probabilistic model illustrated
in Figure 2. The neural network weights θ of the VAE are trained maximizing the Evidence Lower
Bound (ELBO) of the data. Moreover, another generative model was trained on top of the VAE tolearn latent space distribution using Masked Auto Regressive Flows [8]. This model was trained on
the HST/ACS COSMOS survey postage stamps.
In the weak lensing regime, the shear is very small, therefore we can assume that it does not change
the morphological statistics of the training sample. Thus, the model can be trained once, a priori, overa large catalog of galaxies and be used to perform the joint inference of the latent variables and thecosmic shear. For both light profile models, we denote the generated image I0 .
Cosmic shear In this model every galaxy image is sheared by a constant value γ = γ1 + iγ2 sharedby all the galaxies. Theimage shearing
by the following pixel
 transform, denoted I0 ⊗ γ, is1given′γ2 coordinate change: xy′ = S(γ) xy , such that S(γ) = √ 1 2 1+γ
γ2 1−γ1 .
1−|γ|Convolution with the PSF The Point Spread Function (PSF) model is assumed to be known andexact for each galaxy. We draw the PSF evaluated at the galaxy i position, denoted by Πi , in Fourier
space and then multiply it to the Fourier image of the light profile before returning in the real space,resulting in an image I, such that Ĩ = (Ĩ0 ⊗ γ) · Π̃.
Noise corruption We assume a known fixed gaussian noise identical for all the galaxies. For
simplicity, all of the galaxies in our observations are corrupted by a constant Gaussian noise of
2variance σpix
, as well as for our forward model. In a more realistic scenario where we would haveinformation about the noise correlation properties for each postage stamp, one could consider a noisemodel per galaxy. The likelihood is therefore given byp(D|G, γ) ∝


h

i

N di F † Ĩ0,i ⊗ γ · Π̃i , σ 2 , F † is the inverse Fourier operator. (2)Y

3Illustration on realistic data2

i=1,...,N

Sersic / Realistic galaxies
deepgal / Realistic galaxies
0.100 Sersic / Parametric galaxies
ground truth
ground truth
ground truth
0.075
0.050
0.025
0.000
0.025
0.050
0.075
0.1000.10
0.05 0.00
0.05
0.10 0.10
0.05 0.00
0.05
0.10 0.10
0.05 0.00
0.05
0.10
1

1

1

Figure 3: Shear marginal posterior contours against ground truth γ = (−0.03, 0.03). The left panel
shows the posterior using the Sérsic light profile over the parametric observations. The middle panel
shows the posterior using the Sérsic light profile over the realistic galaxies. The right panel shows theposterior using the deepgal model over the realistic galaxies.
3.1

Simulating observations

We designed a toy model where galaxies are assumed to be already detected, deblended and centered in
postage stamps, so that we only need to simulate the light profile, shearing transform and convolution
with the PSF. In order to demonstrate the model bias, we distinguish two sets of observations.
3(1) generated with parametric galaxy light profiles using a Sérsic model and (2) using realistic
galaxies based on the galaxy postage stamps from the COSMOS HST Advanced Camera for Surveys(ACS) field [9–11]. We used the postage stamps from [12] and Sérsic fit parameters compiled in[13]. Because the unconvolved galaxies are very noisy, we used the Auto-Encoder from [7] to get
denoised version of unconvolved galaxies, therefore this is not a demonstration on real galaxies yet.
However, the purpose of this paper is to demonstrate the model bias, so we find that realistic galaxy
morphologies that cannot be captured by analytic profiles are sufficient to illustrate the argument.
Each set of observation is a field of 256 galaxies, generated from the COSMOS catalog, with thefollowing cuts: n > 0.4, rhlr > 0.2 and 23.5 < i+ < 25.2 so rather faint galaxies. That way theentire galaxy is in the postage stamp, is not too peaky and not too resolved to be close to a Sérsic
model. Every postage stamps were corrupted by a Gaussian noise of standard deviation σpix = 0.01.
In order to simulate shearing, we removed any shear present in the selected sample from the COSMOS
catalog, by adopting the same strategy than in the GREAT3 Challenge [14], i.e. doubling each galaxy,
and rotating the doubled galaxies by 90°, as well as their PSF, to cancel the overall sample shear andthen apply the same constant shear to every galaxy γ = (−0.03, 0.03).
3.2

Approximate inference with mean-field VI

In order to estimate the posterior distribution p(G, γ|D), we relied on approximate inference using
mean-field Variational Inference (VI), where
Q Qwe approximate the posterior distribution using aparametric distribution qϕ (G, γ) = qϕγ (γ) i j qϕi,j (zi,j ), parameterized by ϕ. We used Normal
distributions qϕi,j to approximate the posteriror of parameters i of the j th galaxy, so the parameters
ϕi,j = (µi,j , σi,j ) correspond to means and standard deviations. Similarly, qϕγ is a 2d-diagonal
multivariate Gaussian. The parameters ϕi,j and ϕγ are then optimized maximizing the ELBO:ϕ∗ = arg max E(G,γ)∼qϕ [log p(G, γ, D) − log qϕ (G, γ)]ϕ(3)This remains a very simplistic approach, which does not assume any correlation between variables. A
more accurate model would condition the morphological parameters on the shear. In principle, themost accurate method would be to run a Markov Chain Monte Carlo, but several experiments showedthat sampling based-methods are too computationally expensive for such high-dimensional problem.
3.3

Identifying and solving the model bias

We ran three sets of experiments to demonstrate the model bias. First fitting the fully parametricmodel to the parametric observations. Then fitting the same parametric model to the realistic galaxies
and then fitting the hybrid model to the realistic galaxies. We used a Sérsic law as light profileas proposed in [6]. One could argue that there exist more complex parametric profiles such as thebulge+disk model which could reduce the aforementioned bias, but we propose here a model which
does not rely on ellipticities and in principle can capture any complexity of galaxy morphologies.
First, as illustrated in Figure 3 (a), the inferred shear using the parametric model on the parametricobservations is unbiased. This was indeed expected since the observations and the model are exactly
the same. Secondly, Figure 3 (b) shows that the parametric model induces a bias in the shear
estimation when the galaxy morphologies are no longer elliptical. Besides, one can observe thenon-zero residuals in the Figure 1 fit to real galaxies, illustrating the model discrepancy when applied
to complex morphologies. Finally, Figure 3 (c) shows that when we substitute the parametric light
profile model by the deepgal model, the shear estimation becomes unbiased. One can notice thatthe uncertainties using the generative model are smaller than for the Sérsic model, but these contours
are not representative of the error bars due to the mean-field approximation.

4

Discussion and conclusion

We have presented in this work the necessity of an accurate light profile model for forward modelingof the cosmic shear. We have shown that a hybrid physical and deep learning HBM has the potential
to solve the model bias, being able to capture the complex morphology of galaxies and simulating
physical processes such as shearing transform and convolution by the PSF. We identify many
directions of research towards a real cosmological analysis by forward modeling. On the one
4hand, many assumptions could be relaxed by building a more complex hierarchical model, e.g.
parametrizing the PSF model and jointly inferring it along with the shear. Concerning the detection
and deblending of galaxies, similar work have already been proposed, e.g. BLISS [15]. On the other
hand, accurate uncertainty quantification could be enabled by scaling MCMCs for HBMs or usingneural density estimator to capture the shape of the joint posterior distribution.

5

Acknowledgments

BR acknowledges support by the Centre National d’Etudes Spatiales (CNES). This work is supported
by a public grant overseen by the French National Research Agency (ANR) through the program
UDOPIA, project funded by the ANR-20-THIA-0013-01. This work was granted access to the HPC
resources of IDRIS under the allocation 2022-102038 made by GENCI.

6

Broader impact

The methodology described here finds many useful applications within and outside of cosmology.
Hierarchical Bayesian Models are very efficient tools for statistical inference but can suffer form
model misspecification. Combining deep generative models and physical models seems to be a
promising approach to solve this problem. We believe this work does not entail any negative
consequences or ethical issues.

References[1] Rene Laureijs, J Amiaux, S Arduini, J-L Augueres, J Brinchmann, R Cole, M Cropper, C Dabin,
L Duvet, A Ealet, et al. Euclid definition study report. arXiv preprint arXiv:1110.3193, 2011.
[2] Zeljko Ivezić et al. LSST: From Science Drivers to Reference Design and Anticipated Data
Products. The Astrophysical Journal, 873(2):111, 3 2019.
[3] D Spergel, N Gehrels, C Baltay, D Bennett, J Breckinridge, M Donahue, A Dressler, BS Gaudi,
T Greene, O Guyon, et al. Wide-field infrarred survey telescope-astrophysics focused telescope
assets wfirst-afta 2015 report. arXiv preprint arXiv:1503.03757, 2015.
[4] Eric Huff and Rachel Mandelbaum. Metacalibration: direct self-calibration of biases in shear
measurement. arXiv preprint arXiv:1702.02600, 2017.
[5] Erin S Sheldon and Eric M Huff. Practical weak-lensing shear measurement with metacalibration. The Astrophysical Journal, 841(1):24, 2017.
[6] Michael D Schneider, David W Hogg, Philip J Marshall, William A Dawson, Joshua Meyers,
Deborah J Bard, and Dustin Lang. Hierarchical probabilistic inference of cosmic shear. The
Astrophysical Journal, 807(1):87, 2015.
[7] François Lanusse, Rachel Mandelbaum, Siamak Ravanbakhsh, Chun-Liang Li, Peter Freeman,and Barnabás Póczos. Deep generative models for galaxy image simulations. Monthly Noticesof the Royal Astronomical Society, 504(4):5543–5555, 2021.
[8] George Papamakarios, Theo Pavlakou, and Iain Murray. Masked autoregressive flow for density
estimation. Advances in neural information processing systems, 30, 2017.
[9] Anton M Koekemoer, H Aussel, D Calzetti, P Capak, M Giavalisco, J-P Kneib, A Leauthaud,
O Le Fevre, HJ McCracken, R Massey, et al. The cosmos survey: Hubble space telescopeadvanced camera for surveys observations and data processing. The Astrophysical Journal
Supplement Series, 172(1):196, 2007.
[10] Nick Scoville, H Aussel, M Brusa, Peter Capak, C Marcella Carollo, M Elvis, M Giavalisco,
L Guzzo, G Hasinger, C Impey, et al. The cosmic evolution survey (cosmos): overview. The
Astrophysical Journal Supplement Series, 172(1):1, 2007.
5

[11] Nick Scoville, RG Abraham, H Aussel, JE Barnes, A Benson, AW Blain, D Calzetti, A Comastri,
P Capak, C Carilli, et al. Cosmos: Hubble space telescope observations. The AstrophysicalJournal Supplement Series, 172(1):38, 2007.
[12] Alexie Leauthaud, Richard Massey, Jean-Paul Kneib, Jason Rhodes, David E Johnston, Peter
Capak, Catherine Heymans, Richard S Ellis, Anton M Koekemoer, Oliver Le Fevre, et al. Weakgravitational lensing with cosmos: Galaxy selection and shape measurements. The AstrophysicalJournal Supplement Series, 172(1):219, 2007.
[13] Rachel Mandelbaum, Claire Lackner, Alexie Leauthaud, and Barnaby Rowe. Cosmos realgalaxy dataset. Zenodo, January 2012.
[14] Rachel Mandelbaum, Barnaby Rowe, James Bosch, Chihway Chang, Frederic Courbin, Mandeep Gill, Mike Jarvis, Arun Kannawadi, Tomasz Kacprzak, Claire Lackner, et al. The third
gravitational lensing accuracy testing (great3) challenge handbook. The Astrophysical Journal
Supplement Series, 212(1):5, 2014.
[15] Derek Hansen, Ismael Mendoza, Runjing Liu, Ziteng Pang, Zhe Zhao, Camille Avestruz, andJeffrey Regier. Scalable bayesian inference for detection and deblending in astronomical images.
arXiv preprint arXiv:2207.05642, 2022.

Checklist
1. For all authors...
(a) Do the main claims made in the abstract and introduction accurately reflect the paper’s
contributions and scope? [Yes](b) Did you describe the limitations of your work? [Yes](c) Did you discuss any potential negative societal impacts of your work? [Yes](d) Have you read the ethics review guidelines and ensured that your paper conforms tothem? [Yes]2. If you are including theoretical results...
(a) Did you state the full set of assumptions of all theoretical results? [N/A](b) Did you include complete proofs of all theoretical results? [N/A]
3. If you ran experiments...
(a) Did you include the code, data, and instructions needed to reproduce the main experimental results (either in the supplemental material or as a URL)? [No](b) Did you specify all the training details (e.g., data splits, hyperparameters, how they
were chosen)? [No] Because of a lack of space but notebooks reproducing the results
will be released once the paper is accepted. We mentioned the prior parameters but not
the optimization of hyperparameters.
(c) Did you report error bars (e.g., with respect to the random seed after running experiments multiple times)? [No] For this short paper, we show the experiments ran
once.
(d) Did you include the total amount of compute and the type of resources used (e.g.,
type of GPUs, internal cluster, or cloud provider)? [No] because of a lack of space.
Moreover, we did not intend to benchmark computation time here.
4. If you are using existing assets (e.g., code, data, models) or curating/releasing new assets...
(a) If your work uses existing assets, did you cite the creators? [Yes](b) Did you mention the license of the assets? [No](c) Did you include any new assets either in the supplemental material or as a URL? [No](d) Did you discuss whether and how consent was obtained from people whose data you’reusing/curating? [No](e) Did you discuss whether the data you are using/curating contains personally identifiable
information or offensive content? [No]5. If you used crowdsourcing or conducted research with human subjects...
6

(a) Did you include the full text of instructions given to participants and screenshots, if
applicable? [N/A](b) Did you describe any potential participant risks, with links to Institutional Review
Board (IRB) approvals, if applicable? [N/A](c) Did you include the estimated hourly wage paid to participants and the total amountspent on participant compensation? [N/A]A

Appendix

A.1

Prior distributions and parameters for Sérsic light profile

The following distributions were fitted from the COSMOS shape catalog [12, 13].
• log10 rhlr ∼ N (−0.68, 0.30)• log10 n ∼ N (0.1, 0.39)• log10 F ∼ N (−1.97, 0.53)• e1/2 ∼ N (0, 0.28)• γ1/2 ∼ N (0, 0.09)• sx/y ∼ N (0, 1)

7