Relative Intrinsic Scatter in Hierarchical Type Ia Supernova Sibling Analyses: Application to SNe 2021hpr, 1997bq, and 2008fv in NGC 3147

We present Young Supernova Experiment grizy photometry of SN 2021hpr, the third Type Ia supernova sibling to explode in the Cepheid calibrator galaxy, NGC 3147. Siblings are useful for improving SN-host distance estimates and investigating their contributions toward the SN Ia intrinsic scatter ( post-standardization residual scatter in distance estimates ) . We thus develop a principled Bayesian framework for analyzing SN Ia siblings. At its core is the cosmology-independent relative intrinsic scatter parameter, σ Rel : the dispersion of siblings distance estimates relative to one another within a galaxy. It quanti ﬁ es the contribution toward the total intrinsic scatter, σ 0 , from within-galaxy variations about the siblings ’ common properties. It also affects the combined distance uncertainty. We present analytic formulae for computing a σ Rel posterior from individual siblings distances ( estimated using any SN model ) . Applying a newly trained BAYESN model, we ﬁ t the light curves of each sibling in NGC 3147 individually, to yield consistent distance estimates. However, the wide σ Rel posterior means σ Rel ≈ σ 0 is not ruled out. We thus combine the distances by marginalizing over σ Rel with an informative prior: σ Rel ∼ U ( 0, σ 0 ) . Simultaneously ﬁ tting the trio ’ s light curves improves constraints on distance and each sibling ’ s individual dust parameters, compared to individual ﬁ ts. Higher correlation also tightens dust parameter constraints. Therefore, σ Rel marginalization yields robust estimates of siblings distances for cosmology, as well as dust parameters for sibling – host correlation studies. Incorporating NGC 3147 ʼ s Cepheid distance yields H 0 = 78.4 ± 6.5 km s − 1 Mpc − 1 . Our work motivates analyses of homogeneous siblings samples, to constrain σ Rel and its SN-model dependence.

SN Ia siblings -SNe Ia that occur in the same host galaxy -are valuable for investigating standardisation systematics (Elias et al. 1981;Hamuy et al. 1991;Stritzinger et al. 2011;Brown et al. 2015;Gall et al. 2018;Burns et al. 2020;Scolnic et al. 2020Scolnic et al. , 2021;;Graham et al. 2022;Kelsey 2023).Siblings share common properties, like distance, redshift, and host galaxy stellar mass.Consequently, the relative dispersion of siblings distance estimates has no contribution from a variation of these common properties.This mitigation of systematic uncertainties can yield insights into the origin of the total intrinsic scatter.Siblings are also useful for improving distance estimates to the host galaxy, by combining individual siblings distances.
In this work, we develop a principled Bayesian framework for robustly analyzing SN Ia siblings.We do this by introducing the 'relative intrinsic scatter' parameter, σ Rel .This is the intrinsic scatter of individual siblings distance estimates relative to one another within a galaxy.It is constrained without assuming any cosmology, and quantifies the contribution towards the total intrinsic scatter, σ 0 , from within-galaxy variations about the siblings' common properties.
We model σ Rel to analyze a unique trio of SN Ia siblings in NGC 3147 -the only Cepheid calibrator galaxy to host three SNe Ia.We advance on previous state-ofthe-art single-galaxy siblings studies.Like Hoogendam et al. (2022); Barna et al. (2023), we find strong consistency between individual siblings distance estimates, but we compute a σ Rel -posterior to assess the probability that σ Rel > 0. Further, we marginalize over σ Rel with an informative prior, σ Rel ∼ U (0, σ 0 ), to robustly quantify the uncertainty on a combined distance estimate.For the first time, we simultaneously fit siblings light curves to yield improved constraints on distance and host galaxy dust parameters, compared to individ-ual fits.Similar to Biswas et al. (2022), we study the SN Ia color-luminosity relation in a cosmology-independent fashion, but we constrain a physically-motivated common host galaxy dust law parameter R V .Like Gallego-Cano et al. (2022), we use a single siblings calibrator calibrator galaxy to estimate the Hubble constant.Our inference is dominated by the Cepheid-distance measurement error; nonetheless, for the first time, we robustly propagate the siblings' combined uncertainty into a cosmological inference, by marginalizing over σ Rel .
We present our correlated intrinsic scatter model in §2.We then apply these concepts to fit NGC 3147's siblings.The data and model are described in §3 and §4, respectively; this includes new Pan-STARRS-1 grizy photometry of SN 2021hpr from the Young Supernova Experiment (Jones et al. 2021), and a new 'W22' version of BayeSN (Mandel et al. 2022; Appendix A).We perform our analysis in §5, and conclude in §6.

Correlated Intrinsic Scatter Model
We decompose the distance estimate error terms originating from the total intrinsic scatter1 , δM s , into common and relative components.The achromatic magnitude offset common to all siblings in a galaxy is δM Common , and the relative achromatic offsets -specific to each SN sibling -are δM Rel,s .
The population distributions of each component are: δM Rel,s ∼ N (0, σ2 Rel ). (3) Their sum must have a dispersion equal to the total intrinsic scatter, σ 0 : The common intrinsic scatter, σ Common , is thus the contribution towards σ 0 from population variations of the siblings' common properties.Meanwhile, the relative intrinsic scatter, σ Rel , is the contribution towards σ 0 from within-galaxy variations about the siblings' common properties.This simple model can be applied in hierarchical Bayesian settings to robustly analyze siblings, and combine individual siblings distance estimates.
The uncertainty in a joint siblings distance estimate has a contribution, Var[δM ], from the intrinsic scatter components: where δM is the siblings' sample mean δM s .Therefore, if the siblings are assumed to be perfectly correlated, σ Rel = 0, the distance uncertainty contribution is σ 0 (it is maximised).But when the siblings are assumed to be perfectly uncorrelated, σ Rel = σ 0 , the distance uncertainty contribution is minimised: σ 0 / N Siblings .
Computing the precision-weighted average of individual siblings distance estimates is equivalent to adopting this uncorrelated assumption, σ Rel = σ 0 .Consequently, the joint distance uncertainty may be underestimated if in fact there is some correlation, i.e. σ Rel < σ 0 .

σ Rel Prior Knowledge
The size of the relative scatter is uncertain in the literature.On the one hand, Scolnic et al. (2020) and Scolnic et al. (2021) analyzed 8 and 12 siblings galaxies, respectively, using SALT2 (Guy et al. 2007(Guy et al. , 2010)), and results indicated σ Rel is large and comparable to the total intrinsic scatter, σ 0 .This implies the siblings distance estimates are highly uncorrelated.On the other hand, Burns et al. (2020) constrain the dispersion of SNooPy (Burns et al. 2011) distance differences between sibling pairs; after removing fast decliners, and SNe observed with the Neil Gehrels Swift Observatory, they use 11 siblings galaxies to estimate a 0.03 mag 95% upper bound on their dispersion hyperparameter.This is sub-dominant compared to the total intrinsic scatter in the Hubble diagram (typically ∼ 0.1 − 0.15 mag), which implies the siblings distance estimates are highly correlated.It is unclear then how σ Rel compares to σ 0 , and how this depends on the SN sample properties, and the SN model used to estimate distances.
Moreover, σ Rel is constrained without using 'redshiftbased distances': distances obtained by inputting estimates of cosmological redshift into a cosmological model.Therefore, σ Rel has no error contribution from estimating cosmological redshift, and assuming a cosmology; however, these two sources may contribute towards σ 0 if these additional errors -if any -are not already encompassed by the redshift-based distance uncertainties.
The contrast of σ Rel versus σ 0 thus indicates whether it is within-galaxy variations, or the population variation of the siblings' common properties, that dominates σ 0 .We conclude that 0 ≤ σ Rel ≤ σ 0 is expected, and σ Rel can be marginalized over with an informative flat prior: This is justified because the contributions towards σ Rel are a subset of those that contribute towards σ 0 .

Constraining σ Rel
We present analytic formulae for computing a σ Relposterior using individual siblings distance estimates.In a single galaxy, we adopt a simple normal-normal hierarchical Bayesian model, as described in Chapter 5.4 of Gelman et al. (2013).Re-writing in our notation, we have N Siblings individual distance estimates, μs , each with a measurement uncertainty, σfit, s , from fitting each set of siblings light curves, D s .They are estimates of the distance plus the common and relative achromatic offsets, µ + δM Common + δM Rel,s ; therefore, there is a dispersion of the latent variables, which comes from the relative scatter: σ Rel .This is written as: Using eq.5.21 from Gelman et al. (2013) to marginalize over the distance and common offset, µ + δM Common , and with a σ Rel -prior, p(σ Rel ), we can compute an unnormalised relative scatter posterior: where D = {D s } N Siblings s=1 . We note the single-galaxy σ Rel likelihoods are conditionally independent (Bishop 2006).Therefore, a multi-galaxy σ Rel -posterior is computed simply by multiplying the prior by the product of likelihoods in Eq. 8. We use the siblings in NGC 3147 to compute a single-galaxy σ Rel posterior in §5.1.

DATA
We apply these concepts to analyze light curves of NGC 3147's trio of SN Ia siblings.We present new Young Supernova Experiment Pan-STARRS-1 grizy photometry of SN 2021hpr in Table 1 (Chambers et al. 2016).First detected on 2 nd April 2021 (Itagaki 2021), SN 2021hpr is the third spectroscopically confirmed normal Type Ia supernova in NGC 3147 (z = 0.00938), which is a high-stellar-mass Cepheid calibrator galaxy (Rahman et al. 2012;Yim & van der Hulst 2016;Sorai et al. 2019).The data were reduced with Photpipe (Rest et al. 2005), which has been used for numerous photometric reductions of YSE data (e.g.Kilpatrick et al. 2021;Tinyanont et al. 2021;Dimitriadis et al. 2022;Gagliano et al. 2022;Jacobson-Galán et al. 2022a,b;Terreran et al. 2022).The YSE-PZ software was used for a transient discovery alert of SN 2021hpr, and for data collation, management, and visualisation (Coulter et al. 2022;Coulter et al. 2023).We present a spectrum of SN 2021hpr in Appendix B, which indicates it is a normal SN Ia.Recently, SN 2021hpr was studied us-ing optical photometry and spectral observations (Zhang et al. 2022;Barna et al. 2023;Lim et al. 2023).
Fig. 1 displays a Hubble Space Telescope image of NGC 3147, with the trio's locations marked.The image3 was constructed from stacked images produced by and retrieved from the Barbara A. Mikulski Archive for Space Telescopes, obtained from Programmes GO-15145 (PI Riess) and SNAP-16691 (PI Foley) between 29 th October 2017 and 30 th December 2021.The SN 2021hpr inset was created using the Trilogy software (Coe et al. 2012) using data from 30 th December 2021 obtained through Programme SNAP-16691.

MODELING
We use BayeSN to fit each sibling's light curves individually, or all the siblings' light curves simultaneously.BayeSN is an optical-to-near-infrared hierarchical Bayesian model of the SN Ia SEDs (Thorp et al. 2021;Mandel et al. 2022).It uniquely enables us to fit the trio's optical-to-near-infrared light curve data to coherently estimate distance and dust parameters.

New BayeSN Model: W22
The trio's data are in BV RI and grizy passbands, whereas previous iterations of BayeSN were trained either on BV RIY JH or griz, but not both.Therefore, we retrain a new and more robust BayeSN model (hereafter 'W22') simultaneously on optical-NIR BgV rizY JH (0.35-1.8µm) data.The training sample combines the Foundation DR1 (Foley et al. 2018;Jones et al. 2019) and Avelino et al. (2019) samples, so comprises 236 SNe Ia in host galaxies that are ≈ 70% highmass (log 10 (M * /M ⊙ ) > 10).This sample does not include the siblings trio analyzed in this work.Population hyperparameters learned in W22 model training include the global dust law shape, R V = 2.659, and the total intrinsic scatter, σ 0 = 0.094 mag.Appendix A provides further details on this new model.

Fitting Procedures
To correct for Milky Way extinction, we use the Fitzpatrick (1999) law, adopting R V ; MW = 3.1, and a reddening estimate E(B − V ) MW = 0.024 mag from Scolnic et al. (2021).To model host galaxy dust, we fit for each sibling's dust law shape parameter, R s V , using a flat prior: R s V ∼ U (1, 6).The lower bound is based on the Rayleigh scattering limit R V ≈ 1.2 (Draine 2003), and the upper bound is motivated by observational results of lines of sight in the Milky Way (Fitzpatrick 1999;Schlafly et al. 2016).
The time of B-band maximum brightness, T B; max , which defines rest-frame phase via: is fitted for in an individual fit, and thereafter frozen at the posterior mean time.Data outside the model phase range, −10 < t < 40 d, and data with SNR < 3, are removed from the fit.

Joint Siblings Fit Modeling Assumptions
In the new joint fits4 , we fit all light curves of siblings in a single galaxy simultaneously, while applying a common distance constraint; this yields a posterior on a single distance hyperparameter.
To jointly fit the siblings, we must make an assumption about their correlation.While the total

Assumption
σ Rel -prior a The contribution towards the joint siblings distance uncertainty originating from the total intrinsic scatter.b Marginalizing over σ Rel means a posterior on σ Rel is computed, which weights the common distance uncertainty accordingly.The uncertainty contribution given some value of σ Rel is shown in Eq. 5.
intrinsic scatter is learned in BayeSN model training, σ 0 = 0.094 mag, the σ Rel value is unknown.Therefore, we can fit under three δM modeling assumptions shown in Table .3. We can assume the siblings are perfectly uncorrelated, perfectly correlated, or fit for and marginalize over σ Rel while imposing an informative prior: σ Rel ∼ U (0, σ 0 ).Respectively, these are the δM -Uncorrelated, δM -Common, or δM -Mixed assumptions (Table .3).Our default assumption is to marginalize over σ Rel .

Individual Fits
We first fit each of NGC 3147's siblings individually with the new W22 BayeSN model.Fig. 2 shows the SN 2021hpr fit posteriors of (µ, A V , θ, R V ): the distance modulus, dust extinction, light curve shape, and dust law shape, respectively.Appendix B shows the SN 1997bq and SN 2008fv individual-fit posteriors.b The fitting uncertainties, or 'measurement errors', on the individual siblings distance estimates, computed using Eq.11.  1), and kernel density estimates of the posterior distributions.
Table 4 shows posterior summaries.The distance estimates are strongly consistent with one another, with a 0.060 mag sample standard deviation of distance point estimates (Fig. 3).This consistency indicates the data are reliable and the model is robust.However, this does not directly evidence that σ Rel < σ 0 .
Instead, the important information is captured by the posterior on σ Rel .We compute cosmology-independent σ Rel posteriors using the individual distance estimates, and Eqs.(8, 9).We test the sensitivity to the choice of σ Rel prior upper bound: (0.1, 0.15, 1) mag.To estimate distance measurement errors, or BayeSN 'fitting uncertainties', σfit,s , we use where, (σ µ,s , σ0,s ) are the posterior standard deviations of µ s and δM s , respectively5 .Fig. 3 shows the σ Rel -posteriors.Although they peak at zero, the 68% (95%) posterior upper bounds are strongly prior dependent, and extend out as far as σ Rel < 0.21 (0.67) mag.This shows, unsurprisingly, that the σ Rel constraints are weak, and the trio's data alone do not rule out that σ Rel > σ 0 .Therefore, the sample standard deviation of distance point estimates is a poor indicator of σ Rel .More siblings galaxies are required to tightly constrain σ Rel .
Nonetheless, we can robustly combine individual siblings distances by marginalizing over σ Rel , whilst imposing an informative prior, σ Rel ∼ U (0, σ 0 ).To perform this marginalization, we build three simple models in the probabilistic programming language Stan (Carpenter et al. 2017; Stan Development Team 2020), corresponding to the three δM modeling assumptions in   3).The distance uncertainty is minimised when σ Rel = σ0 is assumed, while the best estimate is obtained via σ Rel -marginalization.
Table 3.We run each fit for 100,000 samples, which reduces the Monte Carlo error to < 1 mmag.Results in Fig. 3 show the distance uncertainty is minimised, 0.087 mag, when adopting the δM -Uncorrelated assumption.Meanwhile, marginalizing over σ Rel returns a larger and more robust distance uncertainty: 0.105 mag.This methodology is inexpensive, and can be implemented into cosmological analyses to improve joint siblings distance estimates.

BayeSN Joint Fits
Next, we apply the new BayeSN model architecture ( §4.3), to jointly fit the siblings' light curves simultaneously, while imposing that the siblings have a common, but unknown, distance.Fig. 4 displays the joint-fit posterior from marginalizing over σ Rel .As expected, the constraints on distance improve compared to individual fits (Table 4).Moreover, the constraints on each sibling's individual parameters also improve by jointly fitting the siblings' light curves.
We visualise this in Fig. 5, for SN 2008fv.This SN has optical BV R data only, and so benefits greatly from the joint fit; in particular, the uncertainties in the individual dust parameters (A s V , R s V ) reduce by ≈ (52, 29)%.Further, the individual-fit R s V constraint peaks at the lower prior boundary R s V = 1, while the joint fit rules out low and unphysical values, to constrain R s V = 2.86 ± 0.75.The benefits of siblings data can be understood by drawing parallels with NIR data.Just as NIR data provide added leverage on the dust parameters, and hence the distance, so too siblings data provide added leverage on the distance, which improves constraints on the remaining parameters, like the dust parameters.
To more precisely constrain R V we assume it is the same for all three siblings.The common-R V constraint from marginalizing over σ Rel is R V = 2.62±0.67.This is consistent with the W22 global value, R V = 2.659, and the population mean, R V = 2.70 ± 0.25, constrained in Thorp et al. (2021).Fig. 6 shows the increase in precision on R V constraints, compared to individual fits.
We further fit under the other two intrinsic scatter modeling assumptions: δM -Uncorrelated, and δM -Common.Fig. 6 shows that the R V uncertainties behave in the opposite sense to the distance uncertainty: with larger σ Rel values, there is a larger dispersion of δM parameters, and more freedom in the fit, which leads to larger R V uncertainties.We check and assert that this trend also applies to all other chromatic parameters (like the dust extinction, A s V , and the light curve shape parameters, θ s ).Therefore, like the distance, the best chromatic parameter constraints are obtained by marginalizing over σ Rel .0.2 0.4 0.6 0.8 = 33.169±0.107 BayeSN Common-Fit 2021hpr [grizy] 1997bq [BVRI]

2008fv
[BVR]    show that the more correlated the intrinsic scatter components are, the less freedom there is in the fit, which results in tighter constraints on chromatic parameters.This is opposite to the behaviour of the distance hyperparameter.Marginalizing over σ Rel yields RV = 2.62 ± 0.67.

H 0 Constraints
Our work culminates in an estimate of the Hubble constant, obtained by simultaneously fitting the siblings trio's light curves, a Cepheid distance to NGC 3147 (Riess et al. 2022), and BayeSN distance estimates to 109 Hubble flow SNe Ia.This inference is motivated by the conclusion in Scolnic et al. (2020), which states: "Understanding the limited correlation between SNe in the same host galaxy will be important for SH0ES to properly propagate the combined uncertainty".While understanding the root cause of correlated intrinsic scatter is beyond the scope of this work, our new relative scatter model marginalizes over this correlation, leading to robust uncertainties.Moreover, correlated intrinsic scatter was neglected in the recent siblings-H 0 estimate in Gallego-Cano et al. (2022).This sub-section thus demonstrates how siblings can be robustly modeled to estimate cosmological parameters.
We simultaneously fit the siblings' light curves with the Cepheid and Hubble flow distances, and we marginalize over σ Rel using the informative prior: σ Rel ∼ U (0, σ 0 ).This extracts the most information out of the data via jointly fitting the siblings light curves, whilst accurately quantifying the distance uncertainty via σ Rel marginalization.Appendix C details our H 0 -inference methodologies further.
We estimate H 0 = 78.4± 6.5 km s −1 Mpc −1 .The posterior is displayed in Fig. 7.In Appendix C, we show this result is insensitive to the siblings' R s V modeling assumptions, and their intrinsic scatter modeling assumptions.Assuming the siblings are uncorrelated yields H 0 = 78.4± 6.2 km s −1 Mpc −1 , meaning ≈ 1.7 km s −1 Mpc −1 is added in quadrature to the H 0 uncertainty as result of marginalizing over σ Rel .Equivalently, this is a ≈ 2.2% increase in the H 0 uncertainty.
Our Hubble constant estimate is consistent with typical local Universe measurements, in particular H 0 = 73.04 ± 1.04 km s −1 Mpc −1 (Riess et al. 2022), and also the lower Planck value: H 0 = 67.4± 0.5 km s −1 Mpc −1 (Planck Collaboration et al. 2020).The statistical error dominates, and including more calibrator galaxies will lead to a more accurate and precise H 0 inference, that is less sensitive to Cepheid modeling choices (as seen in figure 18 in Riess et al. 2022; see companion paper: Dhawan et al. 2023).The siblings' correlation is a sub-dominant effect in this case study; nonetheless, we have demonstrated, for the first time, how to robustly propagate the siblings' combined uncertainty to infer cosmological parameters.

CONCLUSIONS
We have analyzed light curves of the SN Ia trio of siblings in the high-stellar-mass Cepheid calibrator galaxy NGC 3147.This includes PS1 grizy photometry of SN 2021hpr from the Young Supernova Experiment (Jones et al. 2021), presented   We marginalize over the total intrinsic scatter, σ0, the siblings' relative intrinsic scatter, σ Rel , and the peculiar velocity dispersion, σpec.Marginalizing over σ Rel while jointly fitting the siblings' light curves yields the best NGC 3147 distance constraint.We estimate H0 = 78.4±6.5 km s −1 Mpc −1 .
• The relative intrinsic scatter, σ Rel , is the intrinsic scatter of individual siblings distance estimates relative to another within a galaxy.It quantifies the contribution towards the total intrinsic scatter, σ 0 , from within-galaxy variations about the siblings' common properties, and is constrained without assuming any cosmology.Therefore, it is distinct from σ 0 , and we expect σ Rel ≤ σ 0 .
• We present analytic formulae for computing a σ Rel posterior from individual siblings distance estimates ( §2.3).Applying this to NGC 3147's siblings, the posterior is wide, meaning the sample standard deviation of distance point estimates is a poor indicator of σ Rel , particularly in this small sample size limit.
• Assuming σ Rel = σ 0 will underestimate the uncertainty on a joint siblings distance estimate if in fact σ Rel < σ 0 .An inexpensive way of constraining µ from individual siblings distance estimates is to marginalize over σ Rel with an informative prior: σ Rel ∼ U (0, σ 0 ).
• For the first time, we hierarchically fit siblings light curves simultaneously, to estimate a common, but unknown, distance hyperparameter.The BayeSN joint fit benefits from a hierarchical sharing of information, returning improved constraints on the common distance, as well as each sibling's individual chromatic parameters (e.g.light curve shape, host galaxy dust parameters).For example, SN 2008fv's individual-fit R s V constraint peaks at R s V = 1, but the joint fit naturally rules out low unphysical values to constrain R s V = 2.86 ± 0.75.• Chromatic parameters are also affected by σ Rel , and in the opposite sense to the distance, with larger σ Rel values returning larger chromatic parameter uncertainties.Therefore, σ Relmarginalization yields robust estimates of chromatic parameters, as well as the distance.We infer a common dust hyperparameter R V = 2.62 ± 0.67.
• We estimate the Hubble constant, H 0 = 78.4± 6.5 km s −1 Mpc −1 , by hierarchically fitting the siblings light curves with a Cepheid distance and Hubble flow distances, whilst marginalizing over (σ Rel , A s V , R s V ) in the calibrator galaxy, and (σ 0 , σ pec ) in the Hubble flow.Marginalizing over σ Rel adds ≈ 1.7 km s −1 Mpc −1 , or ≈ 2.2%, in quadrature to the H 0 uncertainty, compared to assuming σ Rel = σ 0 .This inference is the first cosmological analysis to robustly propagate the combined uncertainty from SN Ia siblings.
This work provides a robust and principled Bayesian framework for hierarchically analyzing SN Ia siblings.Consequently, we have shown the relative intrinsic scatter, σ Rel , is a key parameter in any siblings analysis.Marginalizing over σ Rel yields robust inferences of siblings distances for cosmology, and chromatic parameters for SN-host correlation studies.With more siblings galaxies, there is the potential to tightly constrain σ Rel , and compare it against σ 0 , to investigate the dominant systematics acting in BayeSN and other SN Ia models.
A multi-galaxy σ Rel inference may be sensitive to cross-calibration systematics originating from siblings observed on different photometric systems.While extensive work has gone in to building distance covariance matrices for SALT2 distances in the Pantheon+ sample (Scolnic et al. 2021;Brout et al. 2022), building analogous covariance matrices for BayeSN is beyond the scope of this work.A workaround will be to analyze a sample of SN Ia siblings observed on a single photometric system, thus mitigating cross-calibration systematics.The near-future sample of ZTF DR2 siblings (Dhawan et al. 2022;Graham et al. 2022) is expected to contain ≈ 30 SN Ia siblings galaxies, so will be an interesting test set for investigating σ Rel , and its SN-model dependence.
Beyond constraining σ Rel , siblings are useful for studying within-galaxy dispersions of other features, like (A s V , θ s , R s V ), and host properties (Scolnic et al. 2020(Scolnic et al. , 2021;;Biswas et al. 2022).Identifying and studying features that have within-galaxy dispersions which are smaller than their total dispersions may help to explain any correlations amongst siblings distance estimates.
Intermediate sized samples of SN Ia siblings have already been compiled, numbering ∼ 10 − 40 siblings (e.g.Burns et al. 2020;Scolnic et al. 2020Scolnic et al. , 2021;;Graham et al. 2022).The largest siblings sample to date was recently presented in Kelsey (2023), and comprises 158 galaxies with 327 SNe Ia.The potential for discovering more siblings is promising, with surveys such as the Dark Energy Survey (DES; Scolnic et al. 2020;Smith et al. 2020), the Legacy Survey of Space and Time (LSST; The LSST Dark Energy Science Collaboration et al. 2018;Ivezić et al. 2019;Scolnic et al. 2020), the Nancy Grace Roman Space Telescope (Roman;Hounsell et al. 2018;Rose et al. 2021), the Young Supernova Experiment (YSE; Jones et al. 2021) and the Zwicky Transient Facility (ZTF; Dhawan et al. 2022;Graham et al. 2022).Therefore, SN Ia siblings will be increasingly valuable for studying SN Ia standardisation systematics in cosmological analyses in the years ahead.illustrate the intrinsic effect of W1(t, λ), the first functional principal component (FPC).Dust extinction and residual functions are set to zero, and the coefficient on the first FPC, θ, is varied about the population mean; this variation gives rise to the Phillips (1993) 'broader-brighter' empirical relation in optical passbands, while also affecting the NIR amplitudes and times of the second peak.For visual clarity, the BgV rizY JH light curves have been vertically offset by (1.25, 0, -1, -2.25, -4, 2.5, 0.25, -2.75, -4) mag, respectively.(right panel) Rest-frame intrinsic color curves are synthesised using the new W22 BayeSN model, to illustrate the predicted level of residual intrinsic color variation after correction for dust extinction and light curve shape.
The training procedure, and all hyperpriors, are exactly as described in Mandel et al. (2022), with the exception of the priors on the intrinsic residual standard deviation vector, σ ϵ .Where Mandel et al. (2022) used an uninformative half-Cauchy prior with unit scale on each element of this vector, here, we place an uninformative half-normal prior with a scale of 0.15, i.e.P (σ ϵ,q ) = Half-N (σ ϵ,q |µ = 0, σ = 0.15) for the qth element.This provides better regularization in the z-band and NIR, and can be interpreted as assuming that the residual scatter in any band at any given time is < 0.3 mag with 95% prior probability.We verify the robustness of the W22 model by computing the RMS and σ −pv values in the Hubble diagram, from fits to the full sample of 236 SNe Ia7 .The (RMS, σ −pv ) values are (0.129, 0.113) mag.From W22 fits to each of the T21 and M20 samples, the values are (0.129, 0.117) mag and (0.130, 0.104) mag, respectively.This compares well to the Hubble diagram statistics obtained from fitting T21, (0.124, 0.112) mag, and M20, (0.137, 0.109) mag, to their respective training samples.When fitting W22 to the 40 SN Ia 'NIR@max' sub-sample from Mandel et al. (2022), the fit statistics are (0.093, 0.085) mag, as compared to the M20 values: (0.096, 0.083) mag.
The effect of the first functional principal component (FPC) is showcased in the left panels of Fig. 8, and the population variations in intrinsic color are showcased in the right panels of Fig. 8.The posterior means of the population hyperparameters learned in W22 model training (in addition to the FPCs and residual covariance matrix) are: the global dust law shape, R V = 2.659, the dust extinction hyperparameter, τ A = 0.252 mag, and the total intrinsic scatter, σ 0 = 0.094 mag.

B. ADDITIONAL SIBLINGS DATA & FITS
The SN 2021hpr spectrum was obtained with the Kast spectrograph on the Lick Shane telescope (Miller & Stone 1993) on 13 th April 2021, ≈ 4.2 days before B-band maximum brightness.The data were reduced with standard CCD processing and extractions using a custom data reduction pipeline8 which employs IRAF9 .We fit low-order polynomials to calibration-lamp spectra and perform small shifts based on sky lines in object frames to determine a wavelength solution.We employ custom Python spectrum and remove telluric lines (Silverman et al. 2012).We remove data and linearly interpolate the spectrum between 5650 Å and 5710 Å in the observer frame so as to remove a ghosting artifact.Fig. 10 shows individual BayeSN fits to the Pantheon+ (Scolnic et al. 2021) light curves of SNe 1997bq and 2008fv.Table 5 records posterior summaries from fixed-R V individual fits to the trio (R * V = 2.659).

C. H 0 METHODS & ADDITIONAL RESULTS
In the upper Hubble flow rung, we use distances to a high-stellar-mass subsample of 109 SNe Ia from the Foundation DR1 sample (Foley et al. 2018;Jones et al. 2019); this minimises the systematic contribution from the mass step towards the H 0 inference (because NGC 3147 is also high mass).We use W22 BayeSN R * V = 2.659 photometric distance estimates, and redshift-based distances obtained from estimates of cosmological redshift and adopting our fiducial cosmology from Riess et al. (2016): (Ω M , Ω Λ , H 0 ) = (0.28, 0.72, 73.24 km s −1 Mpc −1 ).We use redshifts from the NASA/IPAC Extragalactic Database (NED), corrected to the CMB frame using the flow model of Carrick et al. (2015) 10 .The NGC 3147 Cepheid distance from Riess et al. (2022) comes from their 'fit 1 Baseline' analysis that uses optical+NIR 'reddening-free' Wesenheit magnitudes.We use the distance derived without the inclusion of any SNe Ia in any host (as in table 6; Riess et al. 2022).We add 0.025 mag in quadrature to the Cepheid distance uncertainty to account for the geometric distance error.The literature Cepheid distance is μCeph NGC 3147 = 33.014± 0.167 mag We sample a 5 log 10 H 0 parameter, and transform the posterior samples to H 0 .We marginalize also over the SN Ia total intrinsic scatter, σ 0 , and the peculiar velocity dispersion, σ pec , which is primarily constrained by the set of Hubble flow BayeSN distance estimates, and the redshift-based distances.The calibrator galaxy's true-distance parameter, µ NGC 3147 , has an uninformative prior.We also marginalize over ∆M 0 , which is the difference between the calibrated absolute magnitude constant of the SNe, and the fiducial value fixed during model training; without the calibrator galaxy, this parameter is degenerate with H 0 .Our priors are: log 10 H 0 ∼ U (log 10 50, log 10 100) (C1) In Table 6, we show H 0 inferences are insensitive to various modeling assumptions in the calibrator galaxy.The default configuration is to simultaneously fit the siblings' light curves with the Cepheid and Hubble flow distances (see 'Siblings Light Curves' block in Table 6).To demonstrate that these methodologies for siblings can be applied to distances obtained with any SN model, we perform additional H 0 inferences using only one of the individual siblings distance estimates ( §4.2; 'Individual Fit Distance'), or all individual siblings distance estimates at the same time ( §4.3; 'All Individual Fit Distances').These 2-step inferences are faster to compute, so we test the sensitivity of H 0 estimates to the R s V assumptions, and the intrinsic scatter modeling assumptions, in the calibrator galaxy.D. SIMULATIONS We present here two sets of simulations to show how the siblings' intrinsic scatter modeling assumptions affect the common distance uncertainties to a single siblings galaxy, and in turn the H 0 uncertainties in the distance ladder.We show that adopting the δM -Uncorrelated assumption returns underestimated distance uncertainties, and hence H 0 uncertainties, compared to those obtained when the true σ Rel is known (provided the true σ Rel < σ 0 ).Firstly, we simulate three individual siblings distance estimates under a true distance µ = 30 mag, with individual fitting uncertainties σfit = 0.06 mag (the median in the W22 sample), a total intrinsic scatter σ 0 = 0.1 mag, and various choices of relative intrinsic scatter: σ Rel = (0, 0.25, 0.5, 0.75, 1) × σ 0 .We then fit for the common distance, with σ 0 = 0.1 mag, and under the three intrinsic scatter modeling assumptions: δM -Uncorrelated (σ Rel = σ 0 ), δM -Mixed (where σ Rel is marginalized over using a uniform prior σ Rel ∼ U (0, σ 0 )) and δM -Common (σ Rel = 0).We perform 100 simulations, then compare against δM -True inferences, when the true value of σ Rel is known in each fit.
Table 7 shows the δM -Uncorrelated uncertainties are always smaller than the δM -True uncertainties when σ Rel < σ 0 .Therefore, adopting the δM -Uncorrelated assumption will always underestimate the common distance uncertainty, unless the siblings are in the extreme regime of being perfectly uncorrelated: σ Rel = σ 0 .The distance uncertainties increase as the siblings are assumed to be more correlated, and the most conservative estimates are obtained by adopting the δM -Common assumption (i.e.σ Rel = 0).In the rows where σ Rel ≲ σ 0 /2, the δM -Mixed uncertainties are slightly underestimated, ∼ O(1 mmag), compared to the δM -True uncertainties, which is reflective of the weak constraining power on σ Rel in the small sample size limit.However, these underestimates are an order of magnitude smaller than the ≈ 0.02 − 0.04 mag underestimates from adopting the δM -Uncorrelated assumption.The δM -Mixed uncertainties match or exceed the δM -True uncertainties for σ Rel ≳ σ 0 /2.Therefore, the most important modeling choice in the small sample size limit, which has the largest impact on siblings inferences, is whether to adopt the δM -Uncorrelated assumption, or either of the δM -Mixed/Common assumptions.
We further explore the effect of relative intrinsic scatter on inferences of H 0 , by assuming calibrator galaxies in the 2nd rung of the distance ladder contain siblings.We simulate 40 calibrator galaxies with sibling pairs, and Cepheid True |/σµ, averaged over 100 simulations, which shows how well estimated the calibrator galaxy distance uncertainties are.For each galaxy, we take the posterior median distance estimate, and divide by the posterior standard deviation, σµ; the vertical lines show the 68% quantiles, which should be close to 1 if the uncertainties are well estimated.The uncertainties are overestimated when adopting the δM -Common assumption (where σ Rel = 0 is assumed), whereas the δM -Uncorrelated uncertainties are underestimated (where σ Rel = σ0 is assumed).The δM -Mixed results (where σ Rel is marginalized over) and the δM -True results (where σ Rel is known and fixed in the model fit) agree with one another; moreover, the 68% quantiles are strongly consistent with 1, indicating the distance uncertainties are well-calibrated.distances with 0.1 mag measurement errors.In the 3rd rung we simulate 100 single-SN Hubble flow galaxies in the redshift range z CMB ∼ U (0.01, 0.1), with σ pec = 150 km s −1 Mpc −1 , and heliocentric redshift measurement errors of 0.0005.The simulated true Hubble constant is H 0 = 70 km s −1 Mpc −1 .As above, we assess recovery of H 0 for different intrinsic scatter modeling assumptions.Table 8 shows the H 0 uncertainties are underestimated with the δM -Uncorrelated assumption, in all cases except for when the true σ Rel = σ 0 .With 40 sibling pairs, σ Rel is constrained well, to the extent that the δM -Mixed uncertainties agree strongly with the δM -True uncertainties for all σ Rel values.
To better intuit these simulations, we plot the distance residuals and their uncertainties in Fig. 11.The left panel shows an example simulation, where the true σ Rel = σ 0 /2.The right panel shows a KDE of the 40 distance residuals, normalized by their uncertainties, averaged over 100 simulations.With well-calibrated uncertainties, the distance residuals should be ≈ 1σ consistent with zero in 68% of simulations; however, the residuals are consistent with zero in less than 68% of simulations if the uncertainties are underestimated, and vice versa.The right panel shows the δM -Common, -Mixed, and -Uncorrelated assumptions lead to overestimated, well-calibrated, and underestimated uncertainties, respectively.Table 6.Hubble constant constraints, using the siblings trio, a Cepheid distance Riess et al. (2022), and a Hubble flow sample of 109 SNe Ia in high-stellar-mass host galaxies from Foundation DR1 (Foley et al. 2018;Jones et al. 2019) posterior summary from the δM -Mixed fit simulations.For the posteriors peaking at the lower or upper prior boundaries, the 68% and 95% quantiles, or the 32% and 5% quantiles, are quoted, respectively.b H0 Inference when σ Rel is known and fixed in the model.

Figure 1 .
Figure 1.Hubble Space Telescope image of NGC 3147 (z Helio = 0.00938), with the SN Ia trio's locations marked, and the length scale in the bottom left.Overlayed is a 15.8 ′′ × 15.8 ′′ (3 kpc on a side) inset region of SN 2021hpr from 30 th December 2021 (≈ 9 months after SN detection).

Figure 3 .
Figure 3. Analysis of individual siblings distance estimates.(top) BayeSN distance estimates from individual fits to the siblings trio; error bars are the distance measurement errors [bold], and the total errors that include σ0 = 0.094 mag [faint].(middle) Posteriors of σ Rel computed using analytic formulae in §2.3.Posteriors are weakly constraining, so more siblings galaxies are required to tightly constrain σ Rel .The 0.06 mag standard deviation is thus a poor estimator of σ Rel .(bottom) The individual distance estimates are combined under the three intrinsic scatter modeling assumptions (Table3).The distance uncertainty is minimised when σ Rel = σ0 is assumed, while the best estimate is obtained via σ Rel -marginalization.

Figure 4 .
Figure 4. Posterior from the joint fit to NGC 3147's siblings' light curves.Constraints on the distance hyperparameter [black],and each sibling's individual parameters[colors], improve compared to individual fits.This fit marginalizes over σ Rel with an informative prior, σ Rel ∼ U (0, σ0), and σ0 = 0.094 mag, to yield robust estimates of the distance and chromatic parameters.

Figure 6 .
Figure 6.Chromatic parameter uncertainties are affected in the opposite sense to the distance hyperparameter by different intrinsic scatter modeling assumptions (Table3).Overlay of distance and common-RV posteriors from individual fits [colors], and joint fits[black].(top panel) Assuming the siblings are uncorrelated minimises the distance estimate uncertainty, while marginalizing over σ Rel yields robust uncertainties.(bottom panel) The common-RV posteriors[black]   show that the more correlated the intrinsic scatter components are, the less freedom there is in the fit, which results in tighter constraints on chromatic parameters.This is opposite to the behaviour of the distance hyperparameter.Marginalizing over σ Rel yields RV = 2.62 ± 0.67.

Figure 7 .
Figure7.Hubble constant inference from simultaneously fitting NGC 3147's siblings' light curves, the Cepheid distance fromRiess et al. (2022), and BayeSN distance estimates to 109 SNe Ia in the Hubble flow.We marginalize over the total intrinsic scatter, σ0, the siblings' relative intrinsic scatter, σ Rel , and the peculiar velocity dispersion, σpec.Marginalizing over σ Rel while jointly fitting the siblings' light curves yields the best NGC 3147 distance constraint.We estimate H0 = 78.4±6.5 km s −1 Mpc −1 .
7. ACKNOWLEDGMENTS by the NSF through grant AST-2108676.N.E.acknowledges funding from the Center for Astrophysical Surveys Fellowship through the National Center for Supercomputing Applications at the University of Illinois at Urbana-Champaign.A.G. is supported by the National Science Foundation Graduate Research Fellowship Program under Grant No. DGE-1746047.A.G. further acknowledges funding from the Center for Astrophysical Surveys Fellowship at UIUC/NCSA and the Illinois Distinguished Fellowship.This work was supported by a VILLUM FONDEN Investigator grant (project number 16599 and 25501).This work made use of the Illinois Campus Cluster, a computing resource that is operated by the Illinois Campus Cluster Program (ICCP) in conjunction with the National Center for Supercomputing Applications (NCSA) and which is supported by funds from the University of Illinois at Urbana-Champaign.

FacilitiesFigure 8 .
Figure8.(left panel) Rest-frame light curves in BgV rizY JH passbands are synthesised using the new W22 BayeSN model, to illustrate the intrinsic effect of W1(t, λ), the first functional principal component (FPC).Dust extinction and residual functions are set to zero, and the coefficient on the first FPC, θ, is varied about the population mean; this variation gives rise to the Phillips (1993) 'broader-brighter' empirical relation in optical passbands, while also affecting the NIR amplitudes and times of the second peak.For visual clarity, the BgV rizY JH light curves have been vertically offset by (1.25, 0, -1, -2.25, -4, 2.5, 0.25, -2.75, -4) mag, respectively.(right panel) Rest-frame intrinsic color curves are synthesised using the new W22 BayeSN model, to illustrate the predicted level of residual intrinsic color variation after correction for dust extinction and light curve shape.Dust extinction and the light curve shape, θ, are set to zero.The solid black lines depict the mean intrinsic color curves, while the dashed lines show a 1 σ variation about the mean, computed by taking the standard deviation of many residual function realisations, each drawn from the intrinsic SED residual covariance matrix.For visual clarity, the colors are offset by(-0.75,0.25, 1.25, 2.5, 3.5, 0.25, 2.25, 3.25)  mag for (B − V , g − V , V − r, V − i, V − z, V − Y , V − J, V − H), respectively.

Figure 9 .
Figure 9. Spectrum of SN 2021hpr, observed ≈ 4.2 d before the B-band maximum, using the Kast spectrograph on the Lick Shane telescope, in the observer-frame wavelength range: 3254-10894 Å.

Figure 10 .
Figure10.Individual BayeSN fits to SNe 1997bq and 2008fv (as in Fig.2).RV posteriors are wide owing to the moderate extinction and lack of NIR data.

Figure 11 .
Figure 11.(left panel) Example simulation showing the distance residuals to 40 calibrator galaxies with sibling pairs, where the true σ Rel = σ0/2 = 0.05 mag.The distance residual, μCalib.− µ Sim True , compares the true distance in the simulation, µ Sim True , with the calibrated SN siblings common distance estimate to the calibrator galaxy, μCalib. .(right panel) We plot a KDE of the 40 distance residuals, |μ Calib.− µ SimTrue |/σµ, averaged over 100 simulations, which shows how well estimated the calibrator galaxy distance uncertainties are.For each galaxy, we take the posterior median distance estimate, and divide by the posterior standard deviation, σµ; the vertical lines show the 68% quantiles, which should be close to 1 if the uncertainties are well estimated.The uncertainties are overestimated when adopting the δM -Common assumption (where σ Rel = 0 is assumed), whereas the δM -Uncorrelated uncertainties are underestimated (where σ Rel = σ0 is assumed).The δM -Mixed results (where σ Rel is marginalized over) and the δM -True results (where σ Rel is known and fixed in the model fit) agree with one another; moreover, the 68% quantiles are strongly consistent with 1, indicating the distance uncertainties are well-calibrated.

Table 1 .
(Jones et al. 2021)y photometry, observed with Pan-STARRS-1 as part of the Young Supernova Experiment(Jones et al. 2021)a .The full data will be made available online upon publication.b Phase is computed with (TB; max, z Helio ) = (59321.46d, 0.00938) using Eq.10. c All photometry in magnitudes.

Table 2 .
Dataset summary of NGC 3147's SN Ia trio.Only the passbands that can be fitted with BayeSN are shown, which excludes the SN 1997bq U -band data.

Table 3 .
Intrinsic scatter modeling assumptions for hierarchically analyzing siblings.

Table 4 .
Posterior summaries of SN parameters from individual fits, and the joint common-µ fit, to the siblings' light curves.The 68 (95)% quantiles are quoted for posterior peaking at zero.The dust law shape priors are R s V ∼ U (1, 6).

Posterior from Individual Distance Estimates
Rel -Rel Priors:Common-Posteriors (

mag): Combination of Individual Distance Estimates
Figure 5. Contrast plot of SN 2008fv fit posteriors showcases the improvement in parameter estimates from the individual fit [pink] to the joint fit [black].For example, the uncertainties on the individual dust parameters, (A s V , R s V ), shrink by ≈ (52, 29)%.With only BV R data in the individual fit, R s V constraints are weak; however, the joint fit rules out low unphysical values R s V ≈ 1.

Table 3
). Overlay of distance and common-RV posteriors from individual fits [colors], and joint fits [black].(top panel) Assuming the siblings are uncorrelated minimises the distance estimate uncertainty, while marginalizing over σ Rel yields robust uncertainties.(bottom panel) The common-RV posteriors [black] in Table 1.We use a new publicly available W22 version of the BayeSN hierarchical optical-to-near-infrared SN Ia SED model ( §4.1, Appendix A), which we have retrained simultaneously on BgV rizY JH (0.35-1.8 µm) data of 236 SNe Ia.We summarise the key conclusions from this work:

Table 5 .
Posterior summaries of SN parameters from individual fits with R * V = 2.659 (i.e.R s V fixed at the W22 training value).The fitting uncertainties, or 'measurement errors', on the individual siblings distance estimates, computed using Eq.11.

Table 7 .
. Intrinsic scatter modeling assumptions in Table3b NGC 3147 is calibrated with an individual-fit distance estimate to one sibling.c NGC 3147 is calibrated using all three individual fit distance estimates.d NGC 3147 is calibrated by jointly fitting the light curves of the siblings trio, while simultaneously fitting the Cepheid and Hubble flow distances.e Priors on RV .f Our best H0 constraint.Common distance uncertainties from three simulated siblings distances, averaged across 100 simulations.The median and sample standard deviation of common-distance uncertainties across 100 simulations (details in text).b Uncertainty on the common-µ inference when the true σ Rel is known and fixed in the model. a