© The Authors 2025

1. Introduction

Type Ia supernovae (SNe Ia) are excellent distance indicators in cosmology, and they are instrumental to the discovery of the accelerated expansion of the Universe (Riess et al. 1998; Perlmutter et al. 1999). SNe Ia are also crucial to precisely measuring dark energy and the Hubble constant (e.g. DES Collaboration 2024; Rubin et al. 2023; Brout et al. 2022a; Riess et al. 2022). In optical wavelengths, which is the regime where most constraints on cosmology from SNe Ia are obtained, standardisation of the SN Ia peak luminosity can reduce the scatter to ∼15%. The peak brightness is corrected for correlations with the light curve width and colour (e.g., Phillips 1993; Tripp 1998) and host galaxy properties (e.g., Kelly et al. 2010; Sullivan et al. 2010). The dependence of width- and colour-corrected luminosity on host galaxy stellar mass, commonly termed the ‘mass step’ is crucial for improving cosmological constraints. The origin of the mass step is poorly understood, but recent studies (e.g. Brout & Scolnic 2021) suggest this could be due to dust and/or intrinsic differences related to astrophysical properties, such as progenitor age (Rigault et al. 2020; Briday et al. 2022). As SN Ia cosmology is currently systematics limited, understanding the standardisation relations is crucial to constraining their cosmology. This is particularly important since several future Stage-IV dark energy missions are designed with a sizable component devoted to a high-redshift SN Ia survey (Hounsell et al. 2018; The LSST Dark Energy Science Collaboration 2018). At low redshift, a large sample of well-characterised SNe Ia has already been obtained by surveys such as the Zwicky Transient Facility (ZTF; Graham et al. 2019; Bellm et al. 2019; Dekany et al. 2020).

Type Ia supernovae siblings, i.e. multiple SNe Ia in the same parent galaxy (Brown 2015), are a powerful route to constraining the standardisation relations we have mentioned. Recently, SN Ia cosmological samples have been analysed using the SALT2 model (Guy et al. 2007, 2010), wherein the distance modulus μ is obtained by correcting the inferred apparent peak magnitude (m_B) for the light curve width (x₁) and colour (c) by the relation

$$\begin{array}{c}\hfill \mu =m_B+\alpha x_1-\beta c-M_B,\end{array}$$(1)

where α and β are derived from a simultaneous fit along with cosmology to minimize the scatter in the Hubble-Lemaitre diagram. In the SALT2 formalism, c is an observed colour, which can be viewed as a combination of the intrinsic colour and dust.

The parameter β, which is central to this work, is empirically derived and captures both the intrinsic and extrinsic colour-luminosity relations. In terms of the latter, β can be viewed as an analogue of the total-to-selective absorption ratio in the B-band, R_B, for a given dust law (Cardelli et al. 1989). A simultaneous cosmology fit using the largest compiled sample of SNe Ia inferred β = 3.04 ± 0.04 (Brout et al. 2022a), which is significantly lower than the R_B ∼ 4.1 seen in the Milky Way (MW; Cardelli et al. 1989; Fitzpatrick 1999). However, other analyses, such as the Union3 compilation (Rubin et al. 2023), have found a higher β compatible with the MW. In cosmological surveys of SNe Ia, it has been noticed that selection effects can lead to incompleteness in the distribution of SNe Ia properties due to correlations with the intrinsic dispersion. These effects will impact the standardisation relations, and they are corrected for using simulations (e.g., Kessler et al. 2019; Popovic et al. 2023). These simulations require detailed inputs of survey observations and the population models derived from the data (e.g. Scolnic & Kessler 2016). Studies exploring the colours of nearby SNe Ia (e,g, Nobili & Goobar 2008) and further expanding the wavelength coverage of the observations from UV to the near-IR (NIR; Burns et al. 2014; Amanullah et al. 2015) indicate a wide range of dust distributions in the interstellar medium (ISM) of SN Ia host galaxies to explain the observed colours. The procedure for the cosmological inference of α and β is convolved with effects such as K-corrections, selection effects, redshift uncertainties, and even MW extinction errors. It is, therefore, important to have independent methods for measuring the standardisation relations. Studies with cosmological samples have shown the likelihood of β values to be dependent on the host galaxy environment (Gonzá lez-Gaitán et al. 2021; Brout & Scolnic 2021), which is crucial for making an inference of cosmology based on current and future samples (e.g. see, Dhawan et al. 2024).

Owing to multiple SNe Ia exploding in the same galaxy, the inference obtained from sibling SNe Ia is insensitive to certain systematics, including cosmological model parameters, peculiar velocity corrections, and global host galaxy dependence. Therefore, it is a robust, independent test of the width-luminosity and colour-luminosity relations, as was demonstrated when constraining the colour-luminosity relation (β) from a single sibling pair in Biswas et al. (2022), where β is 3.5 ± 0.3. In the recent literature, it has been posited that SN Ia siblings could have a smaller dispersion in their luminosity compared to SNe Ia in different galaxies (Burns et al. 2020). This is also seen in the small distance dispersion for the three spectroscopically normal SNe Ia in NGC 1316 (Stritzinger et al. 2010), although the spectroscopically peculiar SN 2006mr has a distance modulus that differs by 0.6 mag from that of the other three. Other studies, however, have found no difference between the scatter in SN Ia siblings and non-sibling SNe Ia (Scolnic et al. 2020, 2022). Apart from understanding and improving the distance measurements for cosmology, comparing siblings also has interesting implications for SN Ia physics. Gall et al. (2018) analysed SN2007on and SN2011iv and found a difference of 14% and 9% in their distances from the optical and NIR, respectively. This was attributed to the differences in the progenitor systems, which are hypothesised to be due to different central densities of the primary white dwarf (e.g. Ashall et al. 2018). It is therefore interesting to study SN Ia siblings to both understand the luminosity corrections and test whether the absence of potential systematics in common can increase the precision of distance measurements. While we can collect a large sample of historical SN Ia sibling data (e.g. Anderson & Soto 2013; Kelsey 2023), as mentioned above, it is necessary to have sibling pairs observed on the same system to avoid systematic errors.

In this paper, we analyse a sample of SN Ia siblings homogeneously observed by ZTF. A large part of the sample of SN Ia siblings is derived from the second data release of SNe Ia observed by ZTF (ZTF DR2). We inferred the SALT2 parameters and subsequently the width-luminosity and colour-luminosity relations as presented in equation (1). With a sizable sample of siblings, we present both a cosmology-independent inference of α and β and an estimate of the observed scatter of the two parameters across the sample. Since all SNe Ia are on the same photometric system, the cross-calibration systematics are minimised, which have been a significant source of errors in cosmological studies (e.g. Scolnic et al. 2022; Brout et al. 2022b). We used the ZTF sample of SN Ia siblings to infer SN light curve parameters and simultaneously constrain the width-luminosity and colour-luminosity relations. We present the method we used in Section 2 and the results in Section 3, and we discuss our findings with respect to the literature, specifically, the cosmological sample of SNe Ia, in Section 4. Finally, we conclude in Section 5.

2. Data and methodology

Initial studies of multiple SNe in the same galaxy often focussed on a single pair or a small set of siblings (Hamuy et al. 1991; Stritzinger et al. 2010). With modern time-domain astronomy surveys having a long survey duration, it has been possible to assemble larger samples of SN Ia siblings (Scolnic et al. 2020, 2022; Burns et al. 2020). The first phase of ZTF (ZTF-I) comprises an optical imaging survey of the entire northern sky with a three-day cadence in the g- and r-bands and a ∼20.5 mag depth that operated between 2018 and 2020. It has been augmented to a two-day cadence since 2020 with its successor ZTF-II. This public g- + r-band survey is complemented with partnership surveys in the i-band and higher cadence observations. The unprecedented scanning speed and depth has made ZTF the ideal facility for discovering and characterising SN siblings (Biswas et al. 2022; Graham et al. 2022). Light curves for the objects in this paper were built using a variant on the standard IPAC forced photometry pipeline (Masci et al. 2019), with more details in associated papers (e.g. Smith et al. in prep.).

We constructed our sibling sample by starting with the sample of spectroscopically confirmed ZTF SNe Ia. We queried all transients classified as SNe Ia using fritz (van der Walt et al. 2019; Coughlin et al. 2023) within a 100 arcsecond radius from the SN Ia coordinates (at z ∼ 0.1, this corresponds to a physical separation of ∼35 kpc). We then saved the sample of pairs for which the second object is also associated with the same host galaxy. From this sample, we removed objects that showed continuous variability for more than 60 days before the date of maximum brightness in order to remove persistent transient sources such as active galactic nuclei (AGNs) and tidal disruption events (TDEs). Details of which sibling pairs passed the sample selection are presented in Section A.

For our analyses, we took a sample of 24 sibling pairs where both objects were spectroscopically classified as SN Ia (hereafter, spec) and a sample of 28 pairs with one object spectroscopically classified as SN Ia and one object photometrically classified as SN Ia (hereafter, photo-spec) that have multi-band light curves from ZTF. Below, we describe the process of determining the objects in the photo-spec sample without a classification as SNe Ia¹. While we do have two subsamples, we homogeneously analysed the entire sample with the same assumptions and selection cuts. However, for our analysis, we also made consistency checks between the spec and photo-spec subsamples and provided joint constraints. Most of the sibling pairs analysed in this work have at least one member in DR2 of ZTF SNe Ia ((Rigault et al. 2025) hereafter R24; Smith et al. in prep., hereafter S24), and a large fraction of them were classified with the Spectral Energy Distribution Machine (SEDM) on the Palomar 60-inch telescope (Blagorodnova et al. 2018; Rigault et al. 2019). We note that given the approximate rate of one SN Ia per galaxy per century, the total number of sibling pairs in our sample is consistent given the size of the entire DR2 sample is ∼3000 SNe Ia. Since previous studies with SN Ia siblings (e.g., Burns et al. 2020) suggest that cross-calibration systematics are a large error source in sibling analyses, we only constructed our sample from sibling pairs where both SNe are observed by ZTF. Our sample of siblings spans a large redshift range of 0.01 < z < 0.1, as shown in Figure 3.

Currently, the most widely used light curve fitting algorithm is the spectral adaptive light curve template SALT2 (Guy et al. 2010), which is based on the SALT method (Guy et al. 2005), and we used this algorithm in our analysis. The SALT2 model treats the colour entirely empirically and is used to find a global colour-luminosity relation, in contrast to models that distinguish between the intrinsic and dust components (e.g. SNooPy Burns et al. (2011), BayeSN (Mandel et al. 2022; Thorp et al. 2021)). We used the updated version of SALT2 presented in Taylor et al. (2021) and as implemented in sncosmo v2.1.0² (Barbary et al. 2016). We applied the model iteratively. In the first iteration, the fit is without the model covariance, and it is only used to guess the time of maximum. The second iteration fits the model only to the data between −10 and +40 days from the first guess time of maximum (see Rigault et al. in prep. for details on selection of the phase range) and with the model covariance to get all the SALT2 fit parameters simultaneously. In the fitting procedure, we corrected the SN fluxes for extinction due to the dust in the MW using extinction values derived in Schlafly & Finkbeiner (2011). We used the widely applied Galactic reddening law, proposed in Cardelli et al. (1989) and known as the ‘CCM’ law, to correct for MW extinction, with the canonical value for the total-to-selective absorption, R_V = R_B − 1 = 3.1. We also compared the inferred parameters with those from an older, more widely used version of SALT2 from Betoule et al. (2014), and we found consistent estimates of the inferred parameters. For the SNe in the sample without a spectroscopic classification, we fit template spectral energy distributions for SN Ib/c, IIn, and IIP, as provided very kindly by Peter Nugent³. Only the pairs where the SN without a spectroscopic classification also prefers a fit to an SN Ia template (by a Δχ² of at least five, though in most cases the fit is overwhelmingly prefers an SN Ia by a Δχ² of ∼50 or greater) are kept in the sample (Figs. 1 and 2).

Fig. 1. Zwicky Transient Facility RGB image of an example sibling pair from our sample: ZTF20abatows and ZTF20abcawtk. The crosses mark the position of the SNe Ia in the field. These siblings were closest in the time separation (∼5 days between peak for the two SNe) between their peaks and hence were detectable at the same time.

Fig. 2. Light curves in the ZTF g, r, i filters along with the SALT2 fits overplotted for the SN Ia pair from the spec subsample with the largest difference in x1, i.e. ZTF18abdmgab and ZTF20abqefja. As discussed in the text, the high Δx1 and low Δc are important for constraining the width-luminosity relation. The difference in x1 in the light curve shapes of the two SNe Ia can be seen as well as the time of the second maximum in the i-band and the r-band shoulder.

In this work, our aim is to infer the standardisation relations between the SN Ia luminosity and the light curve width and colour. The SALT2 model is typically used with a standardisation relation as parametrised in Tripp (1998):

$$\begin{array}{c}\hfill \mu =m_{\mathrm{B}}+\alpha x_1-\beta c-M_{\mathrm{B}}-{\delta }_{\mathrm{host}},\end{array}$$(2)

where the δ_host term corresponds to a ‘step-like’ correction either in global or local host galaxy properties, such as stellar mass and star formation rate (e.g. Briday et al. 2022; Scolnic et al. 2022). If the properties in consideration are global (e.g. MW extinction or peculiar velocity correction), then for SN siblings they would ‘l-out’ when inferring the standardisation relations since both SNe in each pair are in the same host galaxy. However, we note that several studies in the literature find a dependence of the SN Ia luminosity on local environmental properties, which would not cancel out in our analysis (Rigault et al. 2015; Roman et al. 2018; Rigault et al. 2020; Kelsey et al. 2021). We fit for α and β using the derived SALT2 parameters (instead of directly at the light curve level) since no degeneracy between α, β, and the SALT2 parameters has been observed (e.g., see Biswas et al. 2022). We note that no explicit K- or S-corrections are mentioned, as the SALT2 model is a spectral template, and hence, these corrections are taken into account when inferring the light curve parameters m_B, x₁, and c.

In our analyses, we used the Bayes theorem to get the posterior distributions on the parameters of interest, which are α, β, σ_int, and the spread in α and β. We fit for the α and β values, marginalising over the true x₁ and c values. If the observed stretch and colour differences are Δx₁^o, Δc^o, the true values can be written as Δx₁ = Δx₁^o + δx₁, Δc = Δc^o + δc, where δx₁ and δc are the deviations from the true differences. The difference in the observed distance modulus, for example, due to intrinsic differences in the SN siblings, is given by

$$\begin{array}{cc}\hfill \mathrm{\Delta }\mu & =\mathrm{\Delta }m+\alpha \mathrm{\Delta }{x}_1-\beta \mathrm{\Delta }c=\mathrm{\Delta }m+\alpha \mathrm{\Delta }{x}_1^o-\beta \mathrm{\Delta }{c}^o+\alpha \delta x_1-\beta \delta c\hfill \end{array}$$(3)

$$\begin{array}{cc}& \equiv \mathrm{\Delta }{\mu }^o+\alpha \delta x_1-\beta \delta c.\hfill \end{array}$$(4)

The uncertainty of Δμo, that is σΔμo, includes the uncertainties in the observed magnitudes, σmBi2, derived from the output covariance matrix of the SALT2 model fit; the intrinsic magnitude dispersion, σint2; and possible dispersions in the sample of α and β, σα, σβ:

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{\Delta }{\mu }^o}^2={\sigma }_{m_B^1}^2+{\sigma }_{m_B^2}^2+2{\sigma }_{\mathrm{int}}^2+{({\sigma }_{\alpha }\mathrm{\Delta }{x}_1^o)}^2+{({\sigma }_{\beta }\mathrm{\Delta }{c}^o)}^2.\end{array}$$(5)

For our analysis, we marginalised over the true values of x₁ and c, resulting in χ² and the following likelihood expression (see Section B for the derivation):

$$\begin{array}{c}\hfill {\chi }^2=-2ln(L)={\displaystyle \sum a-}\frac{b^2(c_{\alpha }+c_{\beta })}{(c_m{c}_{\alpha }+c_m{c}_{\beta }+c_{\alpha }{c}_{\beta })}+\\ \hfill ln\left[\frac{2\pi (c_m{c}_{\alpha }+c_m{c}_{\beta }+c_{\alpha }{c}_{\beta })}{c_m{c}_{\alpha }{c}_{\beta }}\right],\end{array}$$(6)

where

$$\begin{array}{cc}\hfill a& \equiv \frac{{(\mathrm{\Delta }{\mu }^o)}^2}{{\sigma }_{\mathrm{\Delta }{\mu }^o}^2}\hfill \end{array}$$(7)

$$\begin{array}{cc}\hfill b& \equiv \frac{\mathrm{\Delta }{\mu }^o}{{\sigma }_{\mathrm{\Delta }{\mu }^o}^2}\hfill \end{array}$$(8)

$$\begin{array}{cc}\hfill c_m& \equiv \frac{1}{{\sigma }_{\mathrm{\Delta }{\mu }^o}^2}\hfill \end{array}$$(9)

$$\begin{array}{cc}\hfill c_{\alpha }& \equiv \frac{1}{{\alpha }^2{\sigma }_{\delta x_1}^2}\hfill \end{array}$$(10)

$$\begin{array}{cc}\hfill c_{\beta }& \equiv \frac{1}{{\beta }^2{\sigma }_{\delta c^2}}.\hfill \end{array}$$(11)

In our inference, we fit α, β, and their dispersions as free parameters with uniform priors 0 < α < 1, 0 < β < 8, 0 < σ(α) < 5, 0 < σ(β) < 5, and 0 < σ_int < 2 mag. For the default analysis, we fit with five free parameters, including the intrinsic scatter. In alternate analysis cases, such as where the sample is fitted with a single α and β, we set σ_α, β = 0. We used PyMultiNest (Buchner et al. 2014), a python wrapper to MultiNest (Feroz et al. 2009), to derive the posterior distribution on the parameters. We used the optimal efficiency of the sampling for parameter inference and 1200 live points, i.e. points drawn uniformly from the prior for sampling. We used a tolerance of 0.05 in the error on the evidence as a criterion for the convergence.

3. Results

In this section, we present SN Ia light curve fit parameters and the inferred value of the luminosity-colour and luminosity-light curve width standardisation relations. We fit the SALT2 model to the light curves for our sample. To create the final sample used for computing the standardisation relations, we removed all objects that had an error on time of maximum σ(t₀) > 2 days and observations only in a single filter. Furthermore, we removed objects that lacked adequate sampling at early phases. To quantify this selection criterion, we used the best-sampled SNe, ZTF20acpmgdz and ZTF20achyvas, to perform a test of recovering the SALT2 parameters by downsampling the data and fitting in the absence of early time data. We found that for cases with data at least three days before maximum, we can recover the x₁ and c values from the full light curve; however, if the first observation is at a later epoch, the values are biased by > 2σ compared to the inference from the full light curve. We therefore only selected pairs where both SNe have at least one observation three days before to avoid any biases in the α and β measurements from biased x₁ and c inference. This left us with 12 spec and 13 photo-spec SN Ia pairs, making a total of 25 pairs.

The parameter distributions are shown in Figure 3 and reported in Table 1. The sibling pairs with the same ZTF name have an extremely small spatial separation but are separated in time. Since the aim is to constrain α and β and not cosmological parameters, we did not make selection cuts on the value of x₁ and c in order to allow for the full range of observed light curve widths and colours in our sample. For comparison, in Figure 3, we plot the complete parameter distribution of the ZTF DR2 sample as dashed lines (S24, R24). While the siblings do not extended to the highest redshifts in the DR2 distribution, they span the observed range of x₁ and c values of the entire DR2 sample. We note that for direct comparison, we plot the DR2 sample without the cosmological cuts on x₁ and c. While the c distribution has a high p-value (0.11) for a Kolmogorov-Smirnov test between the DR2 and the sibling samples, the x₁ distribution has a low p-value (0.014), suggesting that even though the siblings span the entire range of observed x₁ values in the DR2 sample, the distribution is not drawn from the same parent population. We note from Figure 3 that the mass distribution for the siblings is skewed towards higher values compared to the values for the DR2 distribution. This would be expected for a sibling sample since it is more likely that larger galaxies produce two SNe Ia. This may suggest that since more massive older galaxies typically host low x₁ events, we see on average more siblings that have a lower x₁ (Rigault et al. 2020; Nicolas et al. 2021).

Fig. 3. Parameters for the sibling SNe Ia in this study. The panels show the redshift (top left), SALT2 x1 (top right), c (bottom left), and host galaxy stellar mass (bottom right) distributions. We did not make any selection cuts on the values of z, x1, and c, unlike for the complete DR2 sample (Rigault et al. 2025). The equivalent distribution for the entire DR2 sample is overplotted as dashed lines Smith et al. in prep.

From Figure 4, we found that the difference in the peak magnitude is correlated more significantly with the difference in colour than the light curve width (the Pearson correlation coefficient is r = 0.42 compared to r = −0.045). Therefore, we expected stronger constraints on the β parameter; however, for our fiducial fit, we simultaneously inferred α and β.

Fig. 4. Difference in the inferred SALT2 mB versus the difference in the inferred x1 (top) and c (middle) for each sibling pair in the sample. The siblings with large differences in x1 (and similar values of c) predominantly constrain α precisely, whereas those with a large Δc constrain β. The colour bar shows the x1 for the wider SN (i.e. higher x1; top) and c for the redder (i.e. higher c; middle) SN Ia in the pair. For better visualisation, we plot the Δx1 versus the Δc in the bottom panel.

The sample has a wide distribution of angular separations for the siblings pairs. In both the spec and phot-spec samples, there are three sibling pairs where the separation is smaller than the pixel size of the camera, and hence, these are ‘same pixel’ siblings (similar to the pair presented in Biswas et al. 2022). This is interesting since the small separation would also indicate that the difference in the properties of the local environment of the SN is small.

We also fit with the dispersion in α and β to account for the scatter in the standardisation relation across the population. This is parameterised by two terms, σ(α) and σ(β), which multiply the Δx₁ and Δc terms and are added to the σ_m. In this formalism, we obtained α = 0.218 ± 0.055 and σ(α)≤0.195 and β = 3.084 ± 0.312 and σ(β)≤0.923, where σ_α, β are the dispersions in α and β. Below, we evaluate the constraints on the standardisation relations and their dispersion for the individual subsamples, spec and phot-spec, respectively. The constraints on the standardisation parameters are summarised in Table 2 When combining all the pairs to constrain α and β under the assumption of a single α and β, i.e. with the dispersion set to zero, we obtained α = 0.228 ± 0.030 and 3.162 ± 0.191 (Figure 5, black contours).

Fig. 5. Constraints on α, β, and the intrinsic dispersion (σint) for the complete (black), spectroscopic (violet), and photometric (brown)samples.

Fig. 6. Constraints on α and β from two sibling pairs. The β constraints are from the pair analysed in Biswas et al. (2022) (high Δc, low Δx1), and the α constraints are from the pair ZTF18abdmgab-ZTF20abqefja (high Δx1, low Δc). The figure illustrates the orthogonality in the constraints from high Δc (low Δx1) and high Δx1 (low Δc) sibling pairs. We chose this as a specific case to illustrate that a single pair will have a degeneracy in the α-β plane. However, this degeneracy can be broken by combining the constraints on α-β with constraints from a siblings pair that has orthogonal degeneracy in the two parameter plane. The combined constraints in black are for illustrative purposes; the final constraints on α and β from the combined sample are more stringent than presented here.

3.1. Spectroscopic SN Ia sample

The SALT2 fit parameters for the spec sample are summarised in Table 1. Unlike the complete ZTF DR2 sample (R24), we did not make selection cuts on the measured value of x₁ and c. This is because the objects with x₁ and c values outside the range of cosmological values allow for the possibility of a large Δx₁ and/or a large Δc, which is important for constraining α and β, respectively. Therefore, the sample has a greater range of observed properties than the complete DR2 sample. The diversity of the entire SN Ia sample from ZTF DR2 is discussed in a companion paper (Dimitriadis et al. 2025). The range of x₁ and c parameters provides a long baseline to fit for the standardisation parameters. For one of the pairs with high Δc (i.e. a difference > 0.1), ZTF20acehyxd+ZTF21abouuow, the colour excess attributed to extinction from MW dust (E(B − V)_MW) is also significant, i.e. 0.463 mag. We therefore tested the assumption of the MW R_V on the inferred β constraint from this SN Ia pair. We varied the MW R_V to the line of sight from 2.5 to 3.5 and found no significant shift in the inferred β value. Hence, for our analysis, we continued to adopt the fiducial MW R_V = 3.1.

Since our aim is to infer both the standardisation parameters and their dispersion in the sample, we fit the entire spec sample together with the likelihood expressed in equation (6). We found a value of α = 0.217 ± 0.061 and β = 3.084 ± 0.740. This would suggest a mean R_V ∼ 2.08 for the sample. Both the α and β dispersion have a high median, though it is consistent within 1σ with no dispersion.

We note that, as expected, the constraints on α are driven by the sibling pairs that have a high Δx₁ (and similar c) and the constraints on β by pairs that have a high Δc (and similar x₁). The similarity in the c values and large differences in x₁ allowed us to break the degeneracy between the α and β constraints since if Δx₁ and Δc were otherwise both large, there would be a strong correlation between the inferred α and β.

3.2. Photometric-spectroscopic SN Ia pairs subsample

Along with the sibling pairs of spectroscopically confirmed SNe Ia analysed in Section 3.1, ZTF has also discovered several pairs of SNe in the same galaxy where one SN in the pair is a spectroscopically confirmed SN Ia and the other is a likely SN Ia based on its light curve. We fit the photo-spec pairs with the same method as the spec pairs and report the SALT2 fit parameters values in Table 1. Similar to the spec sample in Section 3.1, we inferred α, β, and the dispersion. We found α = 0.186 ± 0.091 and β = 3.031 ± 0.501. Similar to the spectroscopic subsample, the median dispersion value is high. However, it is consistent with zero at 1σ. The constraints are shown as brown contours in Figure 5. We found that the parameter inference for the phot-spec sample is consistent with the spec sample, and hence, we could combine the samples to get the most precise constraints on α, β, and the dispersion.

Fig. 7. Constraints on α and β from sibling pairs in low (green) and high (brown) mass host galaxies split at the median log10(M*) = 10.57. The values are consistent at 1.5σ between the subsamples. For comparison, the constraints from the entire sample are plotted as well (black).

3.3. Host galaxy mass dependence

Recent studies have demonstrated that the reddening relations possibly depend on the properties of the host galaxy (e.g. Brout & Scolnic 2021; Gonzá lez-Gaitán et al. 2021). Here, we test whether there is a dependence of the inferred α and β and their dispersions on the host galaxy properties of the sibling pair. From the above results, we found that the spec and photo-spec samples yield consistent values. Therefore, we combined both subsamples for this analysis to gain more statistical precision.

We computed the stellar masses from the g − i colour and i-band absolute magnitude of the host galaxy using the relation provided in Taylor et al. (2011) given as

$$\begin{array}{c}\hfill {\log }_{10}\left(\frac{M}{M_{\odot }}\right)=1.15+0.7(m_g-m_i)-0.4M_i,\end{array}$$(12)

where m_g and m_i are MW-extinction corrected apparent magnitudes in the g- and i-bands and M_i is the absolute i-bandmagnitude.

To analyse the dependence, we split the sample into high and low host mass pairs, similar to analyses in the literature with Hubble residuals (e.g. Brout et al. 2022a; Johansson et al. 2021, and references therein). By dividing the sample into low and high mass bins at log(M_*/M_⊙) = 10, we found consistent results between the two subsamples. However, the statistics in the low mass bin are significantly smaller than for the high mass bin. We therefore split the combined sample at log₁₀(M_*/M_⊙) = 10.57 and found that the α and β values are consistent at 1.5σ.

3.4. Splitting α based on x₁

We tested the sibling sample for evidence of a difference in α between subsamples based on x₁. We divided the sample based on the x₁ of each SN Ia in a sibling pair to constrain α_low and α_high, i.e. a single pair can be fitted with a different α if the x₁ values for the individual SNe are on different sides of the break. In the literature, a broken stretch-luminosity relation fit has been used in the context of cosmological inference with SNe Ia for low-z and high-z samples (Scolnic et al. 2014; Rubin et al. 2015; Burns et al. 2018). We chose a break value of x₁ = −0.49 based on the analysis of the standardisation of the entire ZTF DR2 sample of SNe Ia in (Ginolin et al. 2025). Fitting this broken brightness-stretch relation, the low x₁ subsample has an α_low = 0.274 ± 0.045 and α_high = 0.133 ± 0.072. We perturbed the break value ranging from −1 to 0, including the median value of −0.27, and we found no significant differences in the inferred α values. We surmise that such a difference could arise due to strong selection effects at high z and intrinsic differences in the progenitor population.

4. Discussion

We constrained the standardisation relations of SNe Ia with the spec and photo-spec subsamples and a joint analysis with all the sibling pairs. Accurately constraining the standardisation relations is key to improving cosmological constraints with current and future SN Ia datasets (Brout et al. 2022a; The LSST Dark Energy Science Collaboration 2018).

One open question regarding the value of α and β is whether they have a unique value for all SNe Ia or whether there is scatter in the values across the populations (and specifically whether such a scatter is correlated to, e.g., host galaxy properties) (Brout & Scolnic 2021; Johansson et al. 2021; Wiseman et al. 2023). We note that in the sample of SN Ia siblings presented here, we can constrain the values of α = 0.228 ± 0.030 and β = 3.162 ± 0.191. The value of α is 2σ higher than the Union3 analysis (Rubin et al. 2023), 2.3 σ higher than the inference for the cosmological sample from the recent Dark Energy Survey results (Vincenzi et al. 2024; DES Collaboration 2024), and ∼3σ higher than the value from the Pantheon+ compilation (Brout et al. 2022a). The value for β is consistent (< 1σ difference) with the inference from the cosmological analysis. We inferred the scatter in both α and β by inferring the values along with the σ(α) and σ(β) terms, multiplying the x₁ and c difference in the error term while fitting for the parameters. For the total sample, we inferred a σ(α)≤0.195 and σ(β)≤0.923 at the 95% confidence level. For a robust test of the selection effects, we compared the α and β inference from the more nearby (z < 0.06) and distant (z ≥ 0.06) sibling pairs and found the α and β values are consistent for the two samples. We note that when allowing for the dispersion in the distributions, the median value of α and β for the fiducial case is consistent with the value from the fit to the cosmological SN Ia sample (Brout et al. 2022a). For the α value, this is due to the increase in the error on the median, while the central value does not shift. As a cosmology-independent method, the SN Ia siblings are a robust consistency check of the SN Ia standardisation relations. From Table 2, we observed this for the subsamples. While the central value of the dispersion can be high, it is still consistent with no dispersion in α and β at the ∼1σ level. For comparison, we also fit the individual SNe Ia in the siblings pairs with a cosmology-dependent method (although without any cuts on x₁ and c for the sample), i.e. from minimising the scatter of the Hubble residuals, as is done for cosmological analyses. We found α = 0.21 ± 0.03 and β = 2.84 ± 0.28 when fitting the Hubble residuals, which is consistent with the approach from fitting the sibling pairs in a cosmology-independent way (Table 3).

Compared to previous analyses inferring α and β from siblings in a cosmology-independent way (e.g. Biswas et al. 2022), this analysis constrains both α and β from the siblings alone. In Figure 6, we show how the degeneracy between α and β is broken by combining a sibling pair with high Δc and low Δx₁ from Biswas et al. (2022) with a pair with high Δx₁ and low Δc in our analysis, as opposed to only constraints on β from previous work. Moreover, the constraint from β has a 60% improvement in the uncertainty compared to previous studies, with a more conservative method to estimate uncertainties.

We therefore fit both the spec and phot-spec samples with only a single α and β for the entire population. We note that the individual subsamples each have β values consistent with the inference from the cosmological analyses (Brout et al. 2022a; Rubin et al. 2023; DES Collaboration 2024). However, we found a higher α value at the ∼2.5σ level. When comparing to the global α = 0.168 ± 0.009 inferred from Rubin et al. (2023), the difference is closer to ∼2σ, even though all the cosmological analyses take selection effects into account. An interesting test to analyse the cause of this difference in α would involve taking the input SN Ia population distributions from a volume-limited sample (Ginolin et al. 2025) and computing the selection effects.

We compared the α and β values we obtained with the cosmological value from Brout & Scolnic (2021). We simulated a sample similar to the sibling pairs we observed from the intrinsic colour-luminosity relation (β_int) and the dust properties (e.g. total to selective absorption, R_V) inferred in Brout & Scolnic (2021), convolving the expected scatter in both β_int and R_V. By inferring β from the two effects combined, which is what we are fitting with the sibling sample, we found the value we get is consistent with the β from cosmological analyses (e.g., Brout et al. 2022a; Rubin et al. 2023; DES Collaboration 2024). We performed a similar comparison for our sibling subsamples split by the host galaxy masses to the values from Brout & Scolnic (2021) for each mass bin. For the low-mass subsample, we found consistency with the cosmological analyses; however, the constraints have large errors, as the sample size is small. For the high-mass subsample, we found that the β corresponds to a larger R_V (taking β ∼ R_V + 1) compared to the mean for the high-mass subsample in Brout & Scolnic (2021) by ∼3σ. However, when convolving with the dispersion in the cosmological analyses, we found that the value from the siblings analysis is within range.

We also note that in all the cases studied here, we find that the intrinsic scatter in the sibling sample when doing a pairwise comparison is smaller than the typically observed scatter in the cosmological sample, with a median scatter of σ_int = 0.047 mag and a 95% confidence level upper limit of σ_int ≤ 0.088 mag. For our sample, the σ_int = 0.22 ± 0.04 mag when using the Hubble residuals without a pairwise comparison (i.e. with the same method as used for the complete sample in Rigault et al. (2025)). While on the higher end, this is statistically consistent with the scatter for the complete DR2 sample (Rigault et al. 2025). This has also been seen in the literature with the sample of 12 pairs analysed by Burns et al. (2020). The recent study of Dwomoh et al. (2024) analysed SN Ia siblings in the NIR and found evidence for residual intrinsic scatter. However, they suggest that it could possibly arise from a need for better observations and reduction in the NIR. With future surveys such as the Vera C. Rubin Observatory’s Legacy Survey of Space and Time (LSST), we expect to find ∼800 SN Ia siblings (Scolnic et al. 2020). Such a sample would also be crucial to confirming or refuting the ‘break’ in α between low and high x₁ SNe Ia by improving the constraints by a factor of approximately six. The stated improvements would make the future sibling constraints comparable to the current best constraints on the SN Ia standardisation parameters from cosmological compilations (Brout et al. 2022a; Rubin et al. 2023; DES Collaboration 2024). Given the rates of siblings from LSST, it is possible that a small subsample would also contain two or more SNe Ia, i.e. a SN Ia triplet (e.g. Ward et al. 2023), which can significantly help, depending on the shape and colour of the SNe, in constraining both α and β from a single host galaxy.

5. Conclusions

In this study, we have analysed a uniformly observed sample of sibling SNe Ia, i.e. multiple SNe in the same parent galaxy, from ZTF. This sample is the single largest of SN Ia sibling pairs observed with a single instrument and photometric system, and it allowed us to reduce the uncertainties on the inferred standardisation parameters, α and β. Our sample contains a total of 25 pairs, with 12 pairs having a spectroscopic classification for both SNe Ia and 13 pairs where one object has been spectroscopically classified as an SN Ia and the sibling is a photometricly classified SN Ia. Interestingly, three of the 25 pairs (and a further three that did not pass the quality cuts) were discovered to have a very small separation within the host galaxy and were observed within the same pixel on the detector.

We inferred α and β in a cosmology-independent way by comparing the standardisation of both siblings in the pair. This method does not require a computation of a cosmological distance or a correction of the observed redshift for peculiar motions in the local Universe since both of these quantities are the same for each sibling in the pair. For the fiducial analysis, we inferred α, β, and the spread in their values for the two subsamples. We find that the spec subsample indicates α = 0.217 ± 0.061 and β = 3.084 ± 0.740, and the photo-spec sample indicates α = 0.186 ± 0.091 and β = 3.031 ± 0.501. Both subsamples point to a median β value that is consistent with cosmological analysis, and this value points towards R_V ∼ β − 1 (assuming a single correction for intrinsic and extrinsic effects), which is significantly lower than the canonical MW value of 3.1. These results are consistent with the findings from the single sibling pair of Biswas et al. (2022). However, we note that with a large permissible dispersion in the β values, it is likely that an individual galaxy can have consistent dust properties with that of the MW.

While the fiducial analysis yields a high dispersion value for both α and β, the uncertainties on the dispersion parameters are large enough such that the samples could be consistent with having only a single α and β. We therefore constrained α and β using only a single linear relation and found α = 0.228 ± 0.030 and β = 3.162 ± 0.191. We also subdivided the sample based on host galaxy mass into the canonical low and high mass bins split at log(M_*/M_⊙) = 10.57. The α and β values for the subsamples are consistent and show no strong host galaxy dependence. Future surveys such as LSST are expected to find ∼800 SN Ia siblings, which will be an excellent sample to improve the uncertainties on α and β and understand the scatter in the width-luminosity and colour-luminosity relations (Scolnic et al. 2020).

Data availability

All data associated with this publication is made available via github at https://github.com/ZwickyTransientFacility/ztfcosmo.git as part of the ZTF second data release of Type Ia supernovae.

Acknowledgments

SD acknowledges support from the Marie Curie Individual Fellowship under grant ID 890695, UKRI Horizon Europe Underwriting EP/Z000475/1 and a Junior Research Fellowship at Lucy Cavendish College. This work has been supported by the research project grant “Understanding the Dynamic Universe” funded by the Knut and Alice Wallenberg Foundation under Dnr KAW 2018.0067. AG acknowledges support from Vetenskapsrådet, the Swedish Research Council, project 2020-03444. This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement n°759194 – USNAC) L.G. acknowledges financial support from the Spanish Ministerio de Ciencia e Innovación (MCIN), the Agencia Estatal de Investigación (AEI) 10.13039/501100011033, and the European Social Fund (ESF) “Investing in your future” under the 2019 Ramón y Cajal program RYC2019-027683-I and the PID2020-115253GA-I00 HOSTFLOWS project, from Centro Superior de Investigaciones Científicas (CSIC) under the PIE project 20215AT016, and the program Unidad de Excelencia María de Maeztu CEX2020-001058-M, and from the Departament de Recerca i Universitats de la Generalitat de Catalunya through the 2021-SGR-01270 grant. JHT and KM acknowledge support from EU H2020 ERC grant no. 758638. Based on observations obtained with the Samuel Oschin Telescope 48-inch and the 60-inch Telescope at the Palomar Observatory as part of the Zwicky Transient Facility project. ZTF is supported by the National Science Foundation under Grants No. AST-1440341 and AST-2034437 and a collaboration including partners Caltech, IPAC, the Weizmann Institute of Science, the Oskar Klein Center at Stockholm University, the University of Maryland, Deutsches Elektronen-Synchrotron and Humboldt University, the TANGO Consortium of Taiwan, the University of Wisconsin at Milwaukee, Trinity College Dublin, Lawrence Livermore National Laboratories, IN2P3, University of Warwick, Ruhr University Bochum, Northwestern University and former partners the University of Washington, Los Alamos National Laboratories, and Lawrence Berkeley National Laboratories. Operations are conducted by COO, IPAC, and UW. The Spectral Energy Distribution Machine is based upon work supported by the National Science Foundation under Grant No. 1106171. The ZTF forced-photometry service was funded under the Heising-Simons Foundation grant #12540303 (PI: Graham). The Gordon and Betty Moore Foundation, through both the Data-Driven Investigator Program and a dedicated grant, provided critical funding for SkyPortal.

References

Appendix A: Sample selection

In this section, we present the full list of siblings discovered in the ZTF data stream.

Appendix B: Derivation of the likelihood

We present the derivation of the likelihood expressed in equation (6). With the difference in the observed distance modulus,

$$\begin{array}{cc}\hfill \mathrm{\Delta }\mu & =\mathrm{\Delta }m+\alpha \mathrm{\Delta }{x}_1-\beta \mathrm{\Delta }c=\mathrm{\Delta }m+\alpha \mathrm{\Delta }{x}_1^o-\beta \mathrm{\Delta }{c}^o+\alpha \delta x_1-\beta \delta c\hfill \end{array}$$(B.1)

$$\begin{array}{cc}\hfill & \equiv \mathrm{\Delta }{\mu }^o+\alpha \delta x_1-\beta \delta c,\hfill \end{array}$$(B.1)

with uncertainty

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{\Delta }{\mu }^o}^2={\sigma }_{m_B^1}^2+{\sigma }_{m_B^2}^2+2{\sigma }_{\mathrm{int}}^2+{({\sigma }_{\alpha }\mathrm{\Delta }{x}_1^o)}^2+{({\sigma }_{\beta }\mathrm{\Delta }{c}^o)}^2,\end{array}$$(B.2)

where σmBi2 are the uncertainties in the observed magnitudes, σint2 the intrinsic magnitude dispersion, and σα, σβ possible dispersions in the sample of α and β, the likelihood is given by

$$\begin{array}{c}\hfill L=\mathrm{\Pi }\frac{1}{\sqrt{2\pi {\sigma }_{\mathrm{\Delta }{\mu }^o}^2}}exp[-\frac{1}{2}\frac{{(\mathrm{\Delta }{\mu }^o+\alpha \delta x_1-\beta \delta c)}^2}{{\sigma }_{\mathrm{int}}^2}]\frac{1}{\sqrt{2\pi {\sigma }_{\delta x_1}^2}}exp[-\frac{1}{2}\frac{{(\delta x_1)}^2}{{\sigma }_{\delta x_1}^2}]\frac{1}{\sqrt{2\pi {\sigma }_{\delta c}^2}}exp[-\frac{1}{2}\frac{{(\delta c)}^2}{{\sigma }_{\delta c}^2}].\end{array}$$(B.3)

In order to integrate the likelihood analytically over δx₁ and δc, we define k₁ ≡ αδx₁ and k₂ = βδc to get

$$\begin{array}{cc}\hfill L& =\mathrm{\Pi }\frac{1}{\sqrt{2\pi {\sigma }_{\mathrm{\Delta }{\mu }^o}^2}}exp[-\frac{1}{2}\frac{{(\mathrm{\Delta }{\mu }^o+k_1-k_2)}^2}{{\sigma }_{\mathrm{int}}^2}]\frac{1}{\sqrt{2\pi {\alpha }^2{\sigma }_{\delta x_1}^2}}exp[-\frac{1}{2}\frac{k_1^2}{{\alpha }^2{\sigma }_{\delta x_1}^2}]\frac{1}{\sqrt{2\pi {\beta }^2{\sigma }_{\delta c}^2}}exp[-\frac{1}{2}\frac{k_2^2}{{\beta }^2{\sigma }_{\delta c}^2}]\hfill \end{array}$$(B.4)

$$\begin{array}{cc}\hfill & =\mathrm{\Pi }\frac{1}{\sqrt{2\pi {\sigma }_{\mathrm{\Delta }{\mu }^o}^2}}\frac{1}{\sqrt{2\pi {\alpha }^2{\sigma }_{\delta x_1}^2}}\frac{1}{\sqrt{2\pi {\beta }^2{\sigma }_{\delta c}^2}}exp[-\frac{1}{2}(\frac{{(\mathrm{\Delta }{\mu }^o+k_1-k_2)}^2}{{\sigma }_{\mathrm{int}}^2}+\frac{k_1^2}{{\alpha }^2{\sigma }_{\delta x_1}^2}+\frac{k_2^2}{{\beta }^2{\sigma }_{\delta c}^2})],\hfill \end{array}$$(B.4)

where the α- and β-terms in the square roots in the denominator are important to re-normalise the likelihood correctly. We also introduce (not yet summed over all sibling pairs)

$$\begin{array}{cc}\hfill a& \equiv \frac{{(\mathrm{\Delta }{\mu }^o)}^2}{{\sigma }_{\mathrm{\Delta }{\mu }^o}^2}\hfill \end{array}$$(B.5)

$$\begin{array}{cc}\hfill b& \equiv \frac{\mathrm{\Delta }{\mu }^o}{{\sigma }_{\mathrm{\Delta }{\mu }^o}^2}\hfill \end{array}$$(B.5)

$$\begin{array}{cc}\hfill c_m& \equiv \frac{1}{{\sigma }_{\mathrm{\Delta }{\mu }^o}^2}\hfill \end{array}$$(B.5)

$$\begin{array}{cc}\hfill c_{\alpha }& \equiv \frac{1}{{\alpha }^2{\sigma }_{\delta x_1}^2}\hfill \end{array}$$(B.5)

$$\begin{array}{cc}\hfill c_{\beta }& \equiv \frac{1}{{\beta }^2{\sigma }_{\delta c^2}},\hfill \end{array}$$(B.5)

in terms of which

$$\begin{array}{cc}\hfill L& =\mathrm{\Pi }\frac{1}{\sqrt{2\pi /c_m}}\frac{1}{\sqrt{2\pi /c_{\alpha }}}\frac{1}{\sqrt{2\pi /c_{\beta }}}exp[-\frac{1}{2}(a+c_m{(k_1-k_2)}^2+2b(k_1-k_2)+c_{\alpha }{k}_1^2+c_{\beta }{k}_2^2)].\hfill \end{array}$$(B.6)

Finally, substituting q₁ = k₁ + (b − k₂c_m)/(c_m + c_α) and q₂ = k₂ − bc_α/(c_mc_α + c_mc_β + c_αc_β),

$$\begin{array}{cc}\hfill L& =\mathrm{\Pi }\frac{1}{\sqrt{2\pi /c_m}}\frac{1}{\sqrt{2\pi /c_{\alpha }}}\frac{1}{\sqrt{2\pi /c_{\beta }}}\hfill \end{array}$$(B.7)

$$\begin{array}{cc}\hfill & exp[-\frac{1}{2}(a-\frac{b^2(c_{\alpha }+c_{\beta })}{(c_m{c}_{\alpha }+c_m{c}_{\beta }+c_{\alpha }{c}_{\beta })}+q_1^2(c_m+c_{\alpha })+q_2^2\frac{(c_m{c}_{\alpha }+c_m{c}_{\beta }+c_{\alpha }{c}_{\beta })}{(c_m+c_{\alpha })})].\hfill \end{array}$$(B.7)

Integrating over δx₁ and δc, is equivalent to integrating over q₁ and q₂ from minus to plus infinity where

$$\begin{array}{cc}& {\displaystyle {\int }_{-\infty }^{\infty }}exp-\frac{q_1^2(c_m+c_{\alpha })}{2}d{q}_1=\sqrt{\frac{2\pi }{(c_m+c_{\alpha })}},\hfill \end{array}$$(B.8)

$$\begin{array}{cc}\hfill & {\displaystyle {\int }_{-\infty }^{\infty }}exp-\frac{q_2^2(c_m{c}_{\alpha }+c_m{c}_{\beta }+c_{\alpha }{c}_{\beta })}{2(c_m+c_{\alpha })}d{q}_2=\sqrt{\frac{2\pi (c_m+c_{\alpha })}{(c_m{c}_{\alpha }+c_m{c}_{\beta }+c_{\alpha }{c}_{\beta })}},\hfill \end{array}$$(B.8)

which substituted into the likelihood gives

$$\begin{array}{c}\hfill L=\mathrm{\Pi }\sqrt{\frac{c_m{c}_{\alpha }{c}_{\beta }}{2\pi (c_m{c}_{\alpha }+c_m{c}_{\beta }+c_{\alpha }{c}_{\beta })}}exp[-\frac{1}{2}(a-\frac{b^2(c_{\alpha }+c_{\beta })}{(c_m{c}_{\alpha }+c_m{c}_{\beta }+c_{\alpha }{c}_{\beta })})]\end{array}$$(B.9)

and

$$\begin{array}{c}\hfill {\chi }^2=-2ln(L)={\displaystyle \sum a}-\frac{b^2(c_{\alpha }+c_{\beta })}{(c_m{c}_{\alpha }+c_m{c}_{\beta }+c_{\alpha }{c}_{\beta })}+ln\left[\frac{2\pi (c_m{c}_{\alpha }+c_m{c}_{\beta }+c_{\alpha }{c}_{\beta })}{c_m{c}_{\alpha }{c}_{\beta }}\right].\end{array}$$(B.10)

All Tables

All Figures

	Fig. 1. Zwicky Transient Facility RGB image of an example sibling pair from our sample: ZTF20abatows and ZTF20abcawtk. The crosses mark the position of the SNe Ia in the field. These siblings were closest in the time separation (∼5 days between peak for the two SNe) between their peaks and hence were detectable at the same time.
In the text

Fig. 2. Light curves in the ZTF g, r, i filters along with the SALT2 fits overplotted for the SN Ia pair from the spec subsample with the largest difference in x1, i.e. ZTF18abdmgab and ZTF20abqefja. As discussed in the text, the high Δx1 and low Δc are important for constraining the width-luminosity relation. The difference in x1 in the light curve shapes of the two SNe Ia can be seen as well as the time of the second maximum in the i-band and the r-band shoulder.

In the text

Fig. 3. Parameters for the sibling SNe Ia in this study. The panels show the redshift (top left), SALT2 x1 (top right), c (bottom left), and host galaxy stellar mass (bottom right) distributions. We did not make any selection cuts on the values of z, x1, and c, unlike for the complete DR2 sample (Rigault et al. 2025). The equivalent distribution for the entire DR2 sample is overplotted as dashed lines Smith et al. in prep.

In the text

Fig. 4. Difference in the inferred SALT2 mB versus the difference in the inferred x1 (top) and c (middle) for each sibling pair in the sample. The siblings with large differences in x1 (and similar values of c) predominantly constrain α precisely, whereas those with a large Δc constrain β. The colour bar shows the x1 for the wider SN (i.e. higher x1; top) and c for the redder (i.e. higher c; middle) SN Ia in the pair. For better visualisation, we plot the Δx1 versus the Δc in the bottom panel.

In the text

	Fig. 5. Constraints on α, β, and the intrinsic dispersion (σint) for the complete (black), spectroscopic (violet), and photometric (brown)samples.
In the text

Fig. 6. Constraints on α and β from two sibling pairs. The β constraints are from the pair analysed in Biswas et al. (2022) (high Δc, low Δx1), and the α constraints are from the pair ZTF18abdmgab-ZTF20abqefja (high Δx1, low Δc). The figure illustrates the orthogonality in the constraints from high Δc (low Δx1) and high Δx1 (low Δc) sibling pairs. We chose this as a specific case to illustrate that a single pair will have a degeneracy in the α-β plane. However, this degeneracy can be broken by combining the constraints on α-β with constraints from a siblings pair that has orthogonal degeneracy in the two parameter plane. The combined constraints in black are for illustrative purposes; the final constraints on α and β from the combined sample are more stringent than presented here.

In the text

	Fig. 7. Constraints on α and β from sibling pairs in low (green) and high (brown) mass host galaxies split at the median log10(M*) = 10.57. The values are consistent at 1.5σ between the subsamples. For comparison, the constraints from the entire sample are plotted as well (black).
In the text