Family dispute: do Type IIP supernova siblings agree on their distance?

Context. Type II supernovae provide a direct way to estimate distances through the expanding photosphere method, which is independent of the cosmic distance ladder. A recently introduced Gaussian process-based method allows for a fast and precise modelling of spectral time series, which puts accurate and computationally cheap Type II-based absolute distance determinations within reach. Aims. The goal of the paper is to assess the internal consistency of this new modelling technique coupled with the distance estimation empirically


Introduction
Sibling supernovae are transients that exploded in the same host galaxy.As they are located at essentially the same distance from us, they allow us to test distance estimation methods and investigate their systematics empirically.Such consistency checks were recently carried out in the literature for pairs of Type Ia supernovae (see e.g.Burns et al. 2020;Scolnic et al. 2020) as well as for a pair of a Type Ia and a Type IIP supernova (Graham et al. 2022), yielding good matches.However, no similar tests have yet been conducted on Type II supernovae specifically.
In the past, Type II supernova siblings have been analysed with different goals.For example, Poznanski et al. (2009) used such siblings to check the colour term in their standardizable candle model under the assumption that they share the same distance.More recently, Vinkó et al. (2012) and Szalai et al. (2019) used a Type IIP-IIb and a IIP-IIP pair respectively to constrain the distances to their hosts better than with a single transient.In contrast, here we perform the first dedicated distance consistency check for Type IIP supernova siblings.
Type II supernovae (SNe II) correspond to the final gravitational collapse of massive (≥8 M ), red or blue supergiant stars, which is supported by pre-explosion images (Smartt 2009) and theoretical models (Tinyanont et al. 2022).These supernovae are historically subdivided into two main classes, Type IIP (plateau) and Type IIL (linear), based on their light curves (e.g.Patat et al. 1994).However, there are several indications that these objects can be rather explained as a continuous distribution, as opposed to distinct classes (Anderson et al. 2014;Sanders et al. 2015;Galbany et al. 2016;Morozova et al. 2017;Pessi et al. 2019).
Due to their high intrinsic luminosity and the fact that these supernovae are the most frequent stellar explosions in the Universe (Li et al. 2011), Type II SNe make for excellent distance indicators.To date, mainly two methods have been used to estimate the distances to Type II SNe: the expanding photosphere method (EPM, Kirshner & Kwan 1974), which is a geometric technique relating the photospheric radius of the SN to its angular size, and the standardised candle method (SCM, Hamuy & Pinto 2002), which is based on an empirical relation between the photospheric expansion velocity and the plateau luminosity of the supernovae.However, only the EPM provides a distance estimation that does not require any external calibration, as opposed to the SCM and other distance ladder formalism methods.It is independent of any other distance measurements and hence of the cosmic distance ladder.Owing to these advantages, it has already been applied to various sets of SNe to derive the Hubble constant (see e.g.Schmidt et al. 1992).The classical EPM analysis is, however, prone to several uncertainties as, for example, shown by Jones et al. (2009): the results can be subject to systematic differences depending on which atmospheric model -and thus, which dilution factors are assumed (e.g.Eastman et al. 1996or Dessart & Hillier 2005a) -or which photometric passband set is used for the photometry.To avoid such issues, as pointed out by Dessart & Hillier (2005a), it is necessary to carry out the EPM based on the radiative transferbased modelling of the spectra and estimate the input parameters through these models.We call this augmented version the tailored-EPM analysis, which bears more similarities to the Spectral-fitting Expanding Atmosphere Method (SEAM) introduced by Baron et al. (2004).This step not only allows for a more precise estimation of physical parameters but also avoids the detour of choosing the dilution factors for the EPM.
Although the number of available spectroscopic observations has grown significantly over the past years, the spectral modelling remains a time-consuming and laborious process.To change this situation, Vogl et al. (2020) developed an emulator based on spectra calculated with the Tardis radiative transfer code (Kerzendorf & Sim 2014), which allows for a maximum likelihood parameter estimation and modelling of the spectral time series several orders of magnitude faster than the conventional methods.To showcase the code, Vogl et al. (2020) performed the tailored-EPM analysis for SNe 1999em and 2005cs, showing that a few per cent precision in the derived distance can be achieved.
Here, we attempt to further test the method and investigate its internal consistency empirically by applying it to sibling supernovae.In their case, the maximum possible separation between the siblings is set by the line-of-sight extension of the thick disk of their host galaxies.For face-on galaxies (such as all the hosts in our sample), this is at most about 10 kpc (Gilmore & Reid 1983).At a distance of ∼6 Mpc, for the closest host galaxy in our sample, this corresponds only to a maximum relative error of 0.2%.Hence, the EPM should yield the same distance for the sibling pairs within the uncertainties.As a result, siblings provide a simple empirical way to assess the consistency and robustness of the algorithm: they allow us to test whether we find the same distances for the pairs even though the underlying conditions vary (such as the metallicity, the mass of the progenitor, the amount of ejecta-CSM interaction, the reddening, etc.), as well as the overall data quality and the level of calibration.Such distance comparisons also allow us to assess whether the inferred uncertainties are reasonable.
The paper is structured as follows.In Sect.2, we give a brief overview of the data collected for this study.Section 3 describes the calibration steps we applied to the data to achieve a similar level of calibration for all objects, then provides the background of the emulator-based modelling and gives an outline of the distance estimation.Section 4 shows the results of the modelling for the individual host galaxies, while Sect. 5 discusses these results and details the consistency check of the method.In Sect.6 we present our summary and conclusions.

Data
To obtain a set of objects compatible with the goals of this study, we filtered the catalogue of known supernovae through the Open Supernova Catalog1 (OSC, Guillochon et al. 2017).In addition to looking for supernovae that exploded in the same host, we set further constraints on the individual objects to make sure that our method could be applied to them: the objects should have at least one observation in their early photospheric phase (more precisely, in the epoch range of 10 to 35 days with respect to the time of explosion; see Sect.3.2) to be compatible with our radiative-transfer modelling, a well-covered light curve during these epochs (optimally, in multiple bands for calibration purposes; see Sect.3.1), and a well-constrained time of explosion from either non-detections or from fitting of the rise of the light curve.
By filtering the catalogue, we found that SN IIP pairs in four host galaxies met all the conditions that we described above: NGC 772, NGC 922, NGC 4303 (M 61), and NGC 6946 (Fig. 1).To retrieve the data for the supernovae in these hosts, we made use of the OSC and WISeREP2 (Yaron & Gal-Yam 2012).The properties of the final dataset are summarised in Table 1.

EPM and spectral modelling inputs
Before we could fit the spectral time series and perform the EPM analysis, we needed to obtain the necessary input data: an estimate of the time of explosion, t 0 , photometry interpolated to the spectral epochs, and flux-calibrated spectra.We describe how we inferred t 0 from a parametric fit in Sect.3.1.1.Section 3.1.2explains the Gaussian process (GP) interpolation of the light curves, and Sect.3.1.3details how we used the interpolated magnitudes to recalibrate the spectra.

Time of explosion
The time of explosion is a crucial parameter in the EPM.It sets the size of the photosphere and, thus, the model luminosity.Any error in t 0 causes an error in the distance.In cases where the distance to a supernova is estimated using a single epoch, the uncertainty in the time of the explosion, t 0 , translates directly into the distance uncertainty: based on linear relative error propagation, in accordance with the equations shown in Sect.3.3.The parameters with the ∆ denote the uncertainty of the given values.Assuming a 10% error on the Θ/v measurement for a single spectrum at the epoch of 20 days, a t 0 uncertainty of ±2 days would yield ∼14% error on the final distance, which is too large for our purposes.
Analyzing multiple spectral epochs helps in limiting the final uncertainties of the fit parameters, partially owing to the additional constraining effect exerted on t 0 by the EPM regression.For example, having four high-quality observations in the epoch range of 10 to 25 days, along with the above t 0 uncertainty of 2 days, would yield an EPM distance error of about ∼11%.While the improvement in precision is significant compared to the above case, bringing it to the required levels would require a Table 1.Summary table of the compiled SN sample, along with important properties of the host, namely, its heliocentric redshift (as adopted from the OSC) and reddening caused by the dust in the Milky Way (Schlafly & Finkbeiner 2011).

Host
z helio E(B − V) factor of a few more spectral epochs, which is rarely available.This demonstrates that even with a sufficiently large number of observations, it is crucial to obtain a well-constrained prior estimate on t 0 for the distance determination.
We can use the constraints on t 0 from the early light curve to significantly reduce these uncertainties.Also, if the time of explosion is known to a high level of precision (to an uncertainty of less than a day ideally) and independently of the EPM analysis, even a single spectral epoch is enough to obtain a meaningful distance.Observationally, the common approach for estimating the time of explosion has been to take the midpoint between the first detection and last non-detection (e.g.Gutiérrez et al. 2017).
A129, page 3 of 18 This is not accurate enough for our analysis.The method does not take into account the depth of the non-detections in comparison to the first detection; this can bias the estimated time of explosion and, thus, the distance.Also, the approach neglects the information from data on the rise of the light curve.
To minimise these possible biases and improve the precision, we estimated the time of explosion for each SN through the fitting of their early light curves (the initial plateau, i.e. the first 10−40 days, depending on the individual SNe).Following Ofek et al. (2014) and Rubin et al. (2016), we fit the flux f in band W with the following: where t is the time, t 0 is the time of explosion, f m,w is the peak flux, and t e,W is the characteristic rise time in the particular band.Whenever the data allowed, we carried out this fitting for multiple photometric bands simultaneously to increase the accuracy of the method.In this joint fit, each of the different bands had its own f m,W and t e,W values, but the fits were connected through the global t 0 value, which was the same for all bands.To improve the fitting, we introduced two additional constraints: (i) we took the depth of the non-detections into account by placing a constraint stating that the model flux should not exceed the limiting flux (but the light curves were allowed to extend to times before the non-detection), and (ii) we required that the characteristic time scale, t e , for the light curve rise increases with wavelength for a given supernova, as it was found previously (see, e.g.González-Gaitán et al. 2015).In this way, even the t e values, which otherwise correspond to their own band, were constrained globally in the joint fit.
To fit the model to the set of light curves in different bands simultaneously, we applied the UltraNest3 package (Buchner 2021), which allows for Bayesian inference on complex, arbitrarily defined likelihoods and derives the necessary posterior probability distributions based on the nested sampling Monte Carlo algorithm MLFriends (Buchner 2016(Buchner , 2019)).We assumed a flat prior for each of the input parameters and used a Gaussian likelihood for the parameter inference.To account for a possible underestimation of the photometric errors, which can influence the inferred parameters, we extended the likelihood with an additional term corresponding to the error inflation as described in Hogg et al. (2010).Using this procedure, we obtained a t 0 posterior distribution for each supernova.We then used the maximum of this distribution as our single "best estimate" value for t 0 to set the phase of the spectral observations for the fitting.The full distribution is used as a prior for the EPM analysis, as described in Sect.3.3.This unique approach for obtaining t 0 in parallel to the EPM significantly enhances the precision of the distance measurements, as will be shown in Sect. 4.

Interpolated photometry
The epochs of spectral observations are not necessarily covered by individual photometric datapoints.To obtain photometry at all phases, we fit the light curves using GPs by applying the george4 python package (Ambikasaran et al. 2015).In brief, GPs present a useful way to interpolate the light curves by providing a non-parametric way of fitting, while taking into account the uncertainties of the datapoints.Following several works that already used GPs for the modelling of supernova light curves (Inserra et al. 2018;Yao et al. 2020;Kangas et al. 2022), we adopted covariance functions from the Matérn family for our calculations (Rasmussen & Williams 2006).We chose a smoothness parameter of 3/2 for our analysis.We attempted to keep the length scale of the best-fit GP curve high enough to retain its robustness against the scatter present in the datapoints.To achieve this, we omitted the rise of the light curve and the drop from the plateau from the fits.C.1 lists the magnitudes at the spectral epochs.The interpolated magnitudes are direct inputs for the EPM, as described in Sect.3.3.

Flux-calibrated spectra
A reliable flux calibration is crucial for the spectral modelling and the determination of the host extinction.We therefore applied a linear flux correction to avoid possible biases in the analysis.This correction was based on the photometry: for every epoch, we calculated a set of synthetic magnitudes using the transmission curves from Bessell & Murphy (2012, for BVRI magnitudes) or Dekany et al. (2020, in case of ZTF photometry), and we then compared these values to the corresponding interpolated magnitudes.The correction curve was then calculated by fitting the linear trend present in the ratios of synthetic and interpolated fluxes against the effective wavelengths of the passbands.Finally, the spectra were multiplied with the obtained correction trends.In some cases, additional more complex flux calibration steps were required and these are described separately for each individual supernovae.

Spectral modelling
To fit the individual spectra, we applied the methodology introduced by Vogl et al. (2020).The fitting method is based on a spectral emulator that predicts synthetic spectra for a set of supernova parameters.The training set of synthetic spectra was calculated with a modified version of the Monte Carlo radiative transfer code Tardis (Kerzendorf & Sim 2014;Vogl et al. 2019), which is described in detail in Vogl et al. (2020).To allow for a fast and reliable interpolation of model spectra for a given set of physical parameters, along with correct absolute magnitudes, an emulator was built on this training set with the following procedure.First, the dimensionality of the spectra was greatly reduced by the use of principal component analysis (PCA) then GPs were trained on the physical input parameters (photospheric temperature and velocity, T ph and v ph , metallicity Z, time since explosion t exp , and the exponent of the density profile n = −d ln ρ/d ln r) to interpolate spectra and to predict absolute magnitudes.As described in Vogl et al. (2020), the emulator predicts spectra that in almost all cases match the simulations with a precision of better than 99%, as measured by the mean fractional error.The original version of the emulator covered the early photospheric phase, from 6.5 to 22.5 days after the explosion (Vogl et al. 2020).This was extended towards later epochs up until 38 days, and additional high-temperature models with an NLTE (non-local thermal equilibrium) treatment of He were also added (Vogl 2020, and in prep.;Vasylyev et al. 2022).The 38day upper end is set by the modelling limitations, since after this epoch, time-dependent effects in the excitation and ionization A129, page 4 of 18 Notes.Depending on the spectral epoch, one of these two simulation sets was chosen for building the emulator.
balance tend to become important5 .The physical range covered by the emulator is summarised in Table 2.
To infer the physical parameters, we used a Gaussian likelihood: where, where θ SN = (v ph , T ph , Z, t exp , n) denotes the set of physical parameters, f obs λ and f λ denote the observed and the reddened emulated spectra, respectively, N is the number of spectral bins and C is the bin-to-bin covariance matrix (see e.g.Czekala et al. 2015).The matrix C is important for the inference.It should account for uncertainties in the data and the model, and capture the correlations of these uncertainties across wavelength; otherwise, we would significantly underestimate the parameter uncertainties (see, Czekala et al. 2015).Constructing a matrix with these properties is a challenging, unsolved problem.
We thus used a simple diagonal matrix with constant values as in Vogl et al. (2020).Since we cannot infer reasonable uncertainties under these circumstances, we performed only maximum likelihood estimation for the parameters.Throughout the fitting, the t epoch of each spectrum was fixed.We treated the reddening towards the supernova separately from the other physical parameters: instead of directly including it in the likelihood function, we set up a grid of reddened spectra corresponding to various E(B − V) values and then evaluated the likelihood for each separately.Apart from reducing computational time, this allowed us to quantify the uncertainties caused by the reddening using the treatment described below.We applied the reddening correction according to the Cardelli et al. (1989) law with R V = 3.1.For the lower limit of the E(B − V) grid, we always assumed the Galactic colour excess towards the supernova, which was determined based on the Schlafly & Finkbeiner (2011) dust maps.The best-fit E(B − V) was chosen as the average of the E(B − V) values that resulted in the lowest χ 2 for the individual spectra.
Apart from calculating the best-fit physical parameters, we also evaluated the angular diameter values for every E(B − V) gridpoint to allow for the distance uncertainties resulting from the unknown amount of host reddening to be quantified, as described in the next section.0 .30 0 .3 1 0 .3 2 0 .3 3 0 .3 4 0 .3 5 0 .36 0 .37 0 .38 0 .3 9 0 .40 0 .4 1 0 .4 2 0 .4 3 0 .4 4 0 .4 5 0 .46 0 .47 0 .48  2011) dust map, which sets the lower limit for the KDE.The green and grey curves show the Gaussian kernels for a single observation obtained by our setting and Silverman's rule, respectively.A normalisation procedure was applied to the KDE histogram to obtain a better comparison.For more details on this supernova, see Sect.4.4.

Distance determination
To infer the distances, D, for the supernovae, we used a variant of the tailored-EPM method (Dessart & Hillier 2006;Dessart et al. 2008;Vogl et al. 2020).As a first step, the photospheric angular diameter of the supernova (Θ = R ph /D, where R ph denotes the radius of the photosphere) had to be inferred for each of the spectral epochs.The predicted apparent magnitude for a passband S depends on Θ as follows: Here, Σ * denotes the set of physical parameters corresponding to the best fit, M ph S is the absolute magnitude predicted by the radiative transfer model at the position of the photosphere, and A S denotes the broadband dust extinction in the bandpass; M S denotes the absolute magnitude as defined by the distance modulus formula.The absolute magnitudes are transformed with above formula so that the variations in the size will have no effect.With this definition, we determined the best-fitting angular diameter Θ * by minimising the square difference between the observed m obs S and model magnitudes: Finally, assuming homologous expansion R ph = v ph t (which is well motivated by observations and models for normal SNe IIP; see e.g.Woosley 1988;Dessart & Hillier 2005b), we determined the distance to the supernova and its time of explosion through a Bayesian linear fit to the ratios of the angular diameters and the photospheric velocities (Θ/v ph ) against time, t.In the fit, we assumed Gaussian uncertainties for Θ/v ph of 10% of the measured values for a given colour excess, as in Dessart & Hillier (2006), Dessart et al. (2008), Vogl et al. (2020).We used a flat prior for the distance, whereas for the time of explosion, we used a normalised histogram of the t 0 posterior from the early light A129, page 5 of 18 curve fit as the prior.However, instead of applying the standard χ 2 -based likelihood, we used a different approach to account for the correlated errors caused by the reddening.This is important since the reddening can affect the final EPM distance, as E(B − V) influences the Θ measurements through the de-reddening of the observed magnitudes and changes in the best-fit parameters, such as the photospheric temperature.Hence, the uncertainty of the reddening introduces a correlated error to the final Θ/v ph measurements, which transitions over to the distance obtained from the EPM fit.Depending on the exact model, higher extinction usually leads to shorter inferred distances.To take this into account, we extended the likelihood with the uncertainty of the reddening by applying kernel density estimation.
Kernel density estimation (KDE) is an effective tool for approximating the underlying probability density distribution of an observable using only a number of realizations.We use the scipy KDE implementation assuming Gaussian kernels (Virtanen et al. 2020) to estimate an underlying distribution for the reddening values obtained from the individual spectral fits.The distributions calculated by KDE depend on the bandwidth of the individual Gaussian kernels.We set the bandwidth to a constant value of 0.025.This choice was based on empirical comparison with Silverman's rule (Silverman 1986), which is regularly used for KDE: for multiple epochs, the bandwidth estimate obtained using Silverman's rule is in agreement with our preset value, thus, the resulting distributions match in the two cases.However, Silverman's rule cannot be used for a single epoch or when a sole value of E(B − V) is favoured by all fits; in these cases, our bandwidth choice still provided a realistic KDE.Finally, we set the lower limit of the estimated distribution to the Galactic reddening based on the Schlafly & Finkbeiner (2011) map, to exclude non-physical cases.One of the obtained KDE distributions can be seen in Fig. 2.
To incorporate the correlated uncertainty in the fit, we first drew a large number of reddening samples using the obtained KDE.Then, a respective sample of Θ/v ph values was generated based on the Θ/v ph − E(B − V) linear interpolation and by adding a random offset to each value assuming 10% uncertainty.This sample not only contained the Θ/v ph values for the relevant reddening values only, but it also carried information about the correlated errors present in them.To represent these distributions and set up the final EPM regression, we applied a multi-dimensional Gaussian KDE on the set of Θ/v ph values and then used the corresponding probability density function as the likelihood for the fitting.We then evaluated this likelihood using UltraNest: at each step, the sampler drew an assumed distance and time of explosion value, for which the model (Θ/v ph ) * were calculated and then the log probability density function was evaluated.As a result, we obtained a posterior distribution for the distance that includes the correlated error introduced by the uncertainty in the reddening.

Results
In this section, we present the distances obtained for the individual host galaxies based on their two sibling supernovae.

M61
The spiral galaxy M61 (NGC 4303, Fig. 1) hosted two Type II supernovae, SN 2008in and 2020jfo.This pair of transients has not been investigated together in the literature yet, since the second supernova occurred only recently.
The first supernova, SN 2008in was discovered on JD 2454827.29 by Nakano et al. (2008).As M61 was monitored by the Robotic Optical Transient Search Experiment (ROTSE) starting from the week before the discovery (Roy et al. 2011), the fitting of the early light curve allowed for the estimation of the time of explosion independent of the EPM, from which we obtained t 0 = JD 2454824.51+0.19  −0.14 (Fig. 3).This supernova was found to belong to the peculiar, subluminous group of SNe IIP (Roy et al. 2011).The spectral sequence for our analysis consisted of spectra obtained by Roy et al. (2011) complemented with more recently published spectra from Hicken et al. (2017).
SN 2020jfo was discovered on JD 2458975.70 by the Zwicky Transient Facility (ZTF, Schulze et al. 2020).Owing to the relatively short cadence of observations of ZTF (Bellm et al. 2019), the time of explosion, t 0 , is not only constrained by a nondetection 4 days pre-discovery, but the supernova was discovered on the rise; hence, it was possible to estimate t 0 precisely by the fitting of the early light curve (Fig. 3).We obtained a time of explosion of t 0 = JD2458975.37+0.10  −0.10 .The spectral time series of this object was presented by Sollerman et al. (2021) and Teja et al. (2022).For our study, we used seven early time spectra (t < 20 days).
The spectral time series and the emulator fits are shown in Fig. 4. For the calculation of the best-fitting models the telluric regions (as marked on the figure) were masked.After fitting, we performed the EPM analysis on both supernovae (Fig. 5).The light curve fits used for the flux calibration of the spectra A129, page 6 of 18 4000 5000 6000 7000 8000 9000 [Å]

NGC 772
The spiral galaxy NGC 772 (Fig. 1) is unique in the sense that it not only had been host to two Type II supernovae, but both objects were also observable simultaneously, as they exploded within one and a half months of one another.Both supernovae were followed up by the Carnegie Type II Supernova Program (CATS) and their spectral sequences were previously analysed by Jones et al. (2009).The first supernova, SN 2003hl, was discovered on JD 2452872.0 by Moore et al. (2003).By applying the method described in Sect.3.1 on the early light curve, along with the unfiltered pre-discovery KAIT (Katzman Automatic Imaging Telescope, Filippenko et al. 2001) non-detection on JD 2452863.0,we estimated the time of explosion to be JD 2452864.62 +1.18  −1.15 (Fig. 6).We assigned the KAIT non-detection to the V-band, given that the colour of the supernova is close to zero owing to the combined effect of the very blue spectral energy distribution in such early phases and the reddening, which places the effective wavelength close to the V-band even for a red-sensitive CCD.Only one spectrum was obtained for this supernova in the temporal range covered by our emulator.Nevertheless, using the time of explosion, we could still derive an approximate EPM distance.
The second supernova, SN 2003iq was discovered by Llapasset et al. (2003) on JD 2452921.5 during the photometric follow-up observations of SN 2003hl.Owing to the monitoring of the host, a pre-discovery image was taken on JD 2452918.5 by Llapasset et al. (2003) shortly before the first detection, which already constrained the explosion date of SN 2003iq to a range of 3 days.By fitting the early light curve, we estimated its time of explosion to be JD 2452919.71+0.59  −0.65 (Fig. 6).Apart from the more precisely known t 0 , this supernova has a more thorough spectral record during the photospheric phase than its sibling, as four spectra were acquired during its early evolution.
The fitted spectral time series of the supernovae are displayed in Fig. 7.The EPM regression derived from the obtained A129, page 8 of 18  Jones et al. 2009).However, the differences between the estimates can be explained by the significantly higher colour excess we obtained for SN 2003iq (with our best estimate being E(B − V) = 0.14 mag) and by the fact that we could not make use of more than one spectrum for SN 2003hl, due to the limitations in our modelling approach.On the other hand, by comparing our estimates with previous SCM distances, we find our solution for SN 2003iq is consistent with the previous result of Poznanski et al. (2009;D = 26.6 ± 1.25 Mpc), while our distance for SN 2003hl is not in tension with the SCM estimate of Olivares (2010; D = 25.6 ± 3.30 Mpc).

NGC 922
The peculiar SBcd type galaxy NGC 922 (Fig. 1) hosted two Type II supernovae six years apart : SN 2002gw and SN 2008ho. First, SN 2002gw was discovered on JD 2452560.8 (Monard 2002).Although the latest non-detection occurred too far from the discovery to be useful for constraining the explosion date, the fitting of the early light curve (including the unfiltered CCD observations obtained by Itagaki et al. 2002 andMonard 2002) gave an estimate of JD 2452556.58+0.97  −1.40 (Fig. 9), which is consistent with the value obtained through performing the EPM regression in Jones et al. (2009).In total, three optical spectra were acquired for this supernova by the CATS program (Hamuy et al. 2006) in the epoch range that could be used for our EPM analysis.
The second supernova, SN 2008ho, was discovered on JD 2454796.61 by Pignata et al. (2008).The last pre-discovery non-detection occurred on JD 2454787.77,which provides a weak constraint on the time of explosion: we found t 0 = JD 2454789.63 +3.63  −2.54 from the light curve fit (Fig. 9).In total, two spectra were obtained for this supernova through the Carnegie Supernova Project (Hamuy et al. 2006).The spectra were recalibrated using the BV observations obtained through the CSP project (Anderson et al., in prep., priv. comm.).For the EPM, we made use of the BV photometry from the same campaign.
The fitted spectral sequences of SN 2002gw and SN 2008ho are shown in Fig. 10.From the EPM regression we obtained a distance for SN 2008ho of D = 40.02± 5.07 Mpc (Fig. 11).On the other hand, we found D = 43.46 ± 3.77 Mpc for BV bands only and D = 43.85± 3.78 Mpc for the full BV I set for SN 2002gw (see Fig. 12 for the comparison).Our distance estimates for SN 2002gw fall between the values derived by Jones et al. (2009) for the different dilution factors (D = 37.4 ± 4.9 Mpc and D = 63.9 ± 17.0 Mpc) and agree well with the previous SCM (D = 48.1 ± 6.2 Mpc, Olivares 2010) and photospheric magnitude method estimates (D = 45.10 ± 3.11 Mpc, Rodríguez et al. 2014).Although the difference is only ∼1%, we used the BV instead of the BV I distance for SN 2002gw for the plots and the distance consistency check (Sect.5) to make the analysis of the siblings as similar as possible.We suspect that the majority of the offset between the distances of the two SNe can be attributed to the relatively large uncertainties on the times of explosion.Nevertheless, the two distances are fully consistent within 1σ.

NGC 6946
NGC 6946 (Fig. 1) is a bright face-on SABcd type galaxy, which produced several Type II supernovae: SNe 1948B, 1980K, 2002hh, 2004et, and 2017eaw.Due to its proximity, the distance of this galaxy could be estimated through the tip of the red giant branch (TRGB) method (Anand et al. 2018) and the planetary nebulae luminosity function (PNLF) relation (Herrmann et al. 2008).Three of the aforementioned supernovae (SNe 2002hh, 2004et and 2017eaw) were observed spectroscopically during the early photospheric phase, which makes it possible to measure their distances through the tailored-EPM.However, only SNe 2004et and 2017eaw were optimal for our purposes: even though the spectral time series and the photometric coverage would have allowed the analysis of SN 2002hh, it was found that this supernova exhibited a very strong, two-component reddening (one component arising from the joint effect of the interstellar dust in the host and the Milky Way, the other from large quantities of local dust, see Barlow et al. 2005;Pozzo et al. 2006 andWelch et al. 2007).With the currently available data and methodology, it is not possible to obtain an accurate extinction correction.Considering the very high reddening, any error in our extinction correction would significantly impact the inferred distance.Hence, we chose to exclude SN 2002hh from the analysis.
SN 2004et was discovered on JD 2453273 by Zwitter et al. (2004), and extensively followed up on by the 2-m Himalayan Chandra Telescope (HCT) of the Indian Astronomical Observatory (IAO) and the 3-m Shane telescope at the Lick Observatory.The resulting spectral times series have been previously studied by Sahu et al. (2006).However, the spectra taken with the HCT were subject to calibration issues, which could not be corrected with our standard linear flux re-calibration.Since these issues can influence the fitting significantly, we attempted A129, page 9 of 18 to correct them empirically, using spectra taken with the Keck and Lick telescopes.This procedure is detailed in Appendix A.
We completed this spectral time series with three additional early-time spectra obtained at the David Dunlap Observatory (Takáts & Vinkó 2012, courtesy of József Vinkó).However, since these spectra covered only a very narrow wavelength range, we could not apply the re-calibration procedure from Sect.3.1.3.The pre-discovery non-detection of the supernova and the subsequent early photometry from Li et al. (2004) has allowed for an accurate determination of the time of explosion through the fitting of the early light curve (Fig. 13), which yielded t 0 = JD 2453271.19+0.16 −0.17 .Along with the tight con-straints on the time of explosion, the six spectra obtained in the first month after explosion also allowed for a precise EPM analysis.SN 2017eaw was discovered by Wiggins (2017) on JD 2457885.78, and then extensively followed up on with spectroscopic observations at the Las Cumbres Observatory 1 m telescope and at the McDonald Observatory with the Low-Resolution Spectrograph 2 mounted on the 10 m Hobby-Eberly Telescope.The observations are described in Szalai et al. (2019).The supernova was discovered early, on the rise, and a non-detection is also available from a pre-discovery observation of the host, which yielded a time of explosion of t 0 = A129, page 10 of 18 JD 2457886.01 +0.25 −0.20 .In terms of spectroscopy, this supernova had the best-sampled and most homogenous spectral sequence in our sample, resulting in eight well-calibrated spectra.
The calibrated and fitted spectral sequences are shown in Fig. 14.For SN 2004et, the EPM regression yielded a distance of D = 5.99 ± 0.36 Mpc, while for SN 2017eaw the resulting value was D = 6.44 ± 0.29 Mpc (Fig. 15).The distance for SN 2004et is 20% higher than the classical EPM value inferred by Takáts & Vinkó (2012) and Szalai et al. (2019), while for SN 2017eaw the value we obtained is about 5-20% lower than the value from Szalai et al. (2019; depending on the assumed reddening and the chosen epoch range of spectra of the reference work, see Szalai et al. 2019 for details).We note that our distance estimates are in much better agreement with one another, which highlights the advantage of the tailored-EPM over the classical EPM approach.The distances are also consistent with the various independent distances obtained for NGC 6946: apart from the previous EPM-based values, our estimates are in agreement with SCM distances of Poznanski et al. (2009) andde Jaeger et al. (2017) (both of which yielded D = 6.69 ± 0.30 Mpc for SN 2004et) and the PNLF distance (Herrmann et al. 2008, D = 6.1 ± 0.6 Mpc), but is slightly lower than the latest TRGB estimate (D = 6.95 ± 0.20 Mpc, Anand et al. 2022, from the Extragalactic Distance Database 6 ).

Discussion
As shown in the previous sections, the fitting procedure not only yields distances with a claimed uncertainty of ∼10% or better (Table 3), it does so by requiring only a limited amount of modelling choices.Due to various uncertainties in the observations and the modelling, these distances can be slightly different for the supernova siblings.The main question, in this case, is whether the estimated distances of the siblings are consistent with one another.Assuming an average galaxy, a reliable upper limit for the thickness of the disc is 10 kpc, while the width of the disc is on the order of 100 kpc (Gilmore & Reid 1983;Zanisi et al. 2020).Consequently, in an ideal case, assuming face-on galaxies, the distances inferred for the supernova siblings should match to an uncertainty of ±0.01 Mpc, which corresponds to only ∼0.2% even for the most nearby pair.For more inclined galaxies, the uncertainties can be one order of magnitude larger; however, since the closest hosts we investigated are very likely to be low-inclination galaxies (M 61 and NGC 6946) and the ones farther away are not viewed edge-on either (NGC 722 and NGC 922), the offset between the siblings should remain sub-1% relative to the distance of the galaxy.Hence, the distance estimation to siblings provides an empirical test to assess the effect of these uncertainties and systematics (e.g. the CSM interaction, or data calibration issues, among others).It allows us to test whether our results are not only precise but accurate as well.
To test the consistency of the obtained distances, we performed a Bayesian model comparison.We adopted the methodology used by Wong et al. (2020) for assessing whether pairs of strongly lensed quasars favour a single global set of cosmological parameters or individual cosmologies for each lens that are inconsistent with one another.We applied this method to the distance estimates of the siblings: thus, we asked the question of whether the measured distance posteriors of a pair have more likely been generated from a single underlying distance, D gal , or from two distinct distances, D ind,1 and D ind,2 (and would thus prove inconsistent).The probability ratio of the two scenarios is called the Bayes factor and can be calculated as: where d i and d j denote the distance posteriors of the individual supernovae in each pair.We chose a uniform prior P(D) around the average distance of the host (based on previous measurements, as quoted on NASA/IPAC Extragalactic Database, NED7 ), ranging from half to twice that distance.The quoted Bayes factors can thus be understood as lower limits since they scale linearly with the width of the distance prior.If the obtained Bayes factor exceeds unity, this can be interpreted as a sign of consistency between the two posterior distributions in accordance with the table of Kass & Raftery (1995).The higher the F value, the stronger the consistency between the siblings.The distance posteriors from the EPM analysis are shown in Fig. 16, while the obtained Bayes factors are listed in Table 4.The consistency of the distance posteriors for the various supernovae is good, with Bayes factors being close to or over 5 for two cases, showing remarkable agreement.In the case of NGC 772 and NGC 922 the values of the Bayes factors are lower, which might be attributed to the lack of a well-sampled spectral 4000 5000 6000 7000 8000 9000  sequence for SN 2003hl and the lack of well-constrained explosion date for SNe 2002gw and 2008ho.However, the remaining cases show that the fitting method and the analysis result in distances that are not only consistent for supernovae in the same host but are also precise to a degree of a few percent.The analysis also showed that some of the greatest limiting factors of the method are: an insufficient amount of spectra, poor quality and calibration-level of the obtained spectral times series, and weak constraints on the time of explosion.If all the conditions for the above were met, the EPM analysis would be able to yield highly precise and consistent values for the supernova pairs, as in the case of SNe 2004et and 2017eaw or in the case of SNe 2008in and 2020jfo.
We offer the caveat that some of our results may be influenced by the choice of R V = 3.1.If the true R V is different, this could, in principle, cause notable offsets between the sibling distances, provided that the differential reddening between the objects is non-negligible and the measurement uncertainties are low enough.In our sample, however, the pairs with the most significant differential reddening had also the largest uncertainty on A129, page 12 of 18   their distance and vice versa.Hence, it is currently not possible to estimate the scale of this effect on the basis of our four sibling pairs.

Summary
In this work, we investigated the consistency of the EPM distances of sibling Type II SNe, calculated based on the spectral fitting method introduced by Vogl et al. (2020).According to our analysis, the method not only yields precise results with an estimated individual distance uncertainty that is better than ∼10%, despite using literature data from a wide range of sources, but the resulting distances are consistent for SNe that are within the same host galaxy.The degree of consistency depends on how well the supernovae were observed, with the worst results occurring when the quality of the data barely allows for the EPM analysis, while the best consistency is achieved when both supernovae are similarly well observed and have a good constraint on the time of explosion.In the cases of M 61 and NGC 6946, highly consistent distances are derived with mismatches of less than 5%.
The high level of consistency between siblings also shows that if there are any other systematics present between the various SNe, they are likely subdominant compared to the effect of observation quality.Such systematics may include different reddening, CSM interaction, explosion energetics, and metallicities.Checking for the effect and scale of such systematics will require a larger set of siblings with data quality that matches that of the siblings of M 61 or NGC 6946.
Apart from checking for the internal consistency of the method, we obtained precise distances for the four investigated host galaxies, even though we used literature data with a highly varying level of calibration.Furthermore, the tailored-EPM provides absolute distances is physics-based and is affected by different systematics than the other existing distance estimation techniques; thus, it can be utilised completely independently of them.These properties and the presented analysis in this paper demonstrate the potential of the tailored-EPM to provide not only precise distances but a precise Hubble parameter as well.Such an analysis is currently being conducted on high-quality spectral time series obtained for SNe II in the Hubble Flow by the adH0cc collaboration (accurate determination of H0 with core-collapse supernovae8 ), which will provide important clues to the Hubble tension.
Figure C.1 shows the interpolated light curves and Table

Fig. 2 .
Fig. 2. Distribution of best-fitting reddening estimates for all the epochs of SN 2017eaw (blue bars) and the constructed KDE distribution (background shape).Dashed grey line at E(B − V) = 0.3 mag denotes the Galactic reddening towards NGC 6946 based on the Schlafly & Finkbeiner (2011) dust map, which sets the lower limit for the KDE.The green and grey curves show the Gaussian kernels for a single observation obtained by our setting and Silverman's rule, respectively.A normalisation procedure was applied to the KDE histogram to obtain a better comparison.For more details on this supernova, see Sect.4.4.

Fig. 3 .
Fig. 3. Exponential fit on the ROTSE light curve of SN 2008in (left) and the ZTF early light curves of SN 2020jfo (right).The blue shaded areas shows the obtained t 0 posteriors.The black curves show the mean fit, namely, the fit that results by taking the mean of the posterior distributions.In case of SN 2020jfo, the g-band light curve was further rescaled after the normalization to improve clarity on the figure.The shaded regions denote the 68% and 95% confidence intervals.

Fig. 13 .
Fig. 13.Exponential fit on the early light curves of SN 2004et (left) and SN 2017eaw (right), along with the determined time of explosion posteriors (blue shaded regions).

Fig. 16 .
Fig. 16.Distance posteriors obtained from the EPM analysis of the supernovae for the various host galaxies.The calculated Bayes factors are displayed in the top-left corners of the panels.
Fig. A.1.Comparison of the spectra of SN 2004et taken at the HCT and the Lick observatory.The discrepancy between the two spectral sets below the 7000 Å is clear.By comparing the spectra of SN 2004et taken at the HCT and Lick Observatory at similar epochs, we found that the former show a significant flux deficit in the wavelength range covered mainly by the Bessell V and R bands, as shown in the Fig. A.1.This trend could not be corrected by a linear flux correction based on the photometry.To fix this issue, we took the ratios of the uncorrected HCT and the re-calibrated Keck and Lick spectra taken at close epochs (which otherwise are too old for the emulatorbased fitting, see top panel of Fig. A.2), then averaged the obtained trend, assuming that it is the same for every HCT spectrum, and fitted it using Generalized Additive Models (GAMs, Hastie & Tibshirani 1990, by applying the pyGAM package from Python, Servén et al. 2018), as shown in the bottom panel of Fig. A.2.We then applied the fitted trend on all the HCT spectra by dividing them with this calibration curve (assuming this calibration trend remains the same regardless of epoch).With this empirical correction, we obtained the HCT spectral sequence shown in Fig. A.3, which we passed to the regular relative flux calibration routines described in Sect.3.1.

Fig
Fig. A.2. Differences between the HCT and Lick observatory spectra shown in A.1.Top panel: Ratios of selected HCT and Keck/Lick spectra that are close in time.Bottom panel: Average curve calculated from the four individual ones on the top panel, along with its GAM fit.

Table 2 .
Parameter range covered by the extended spectral emulator.

Table 3 .
Summary table of the EPM results.Notes.t 0 denotes the time of explosion as obtained by the EPM regression.The E(B − V) values listed here are those determined by the spectral fitting.The values in brackets denote the uncertainties of the given quantities.In the starred (*) cases we either had only one spectrum to fit or both spectra yielded the same reddening value, which did not allow us to give a reasonable empirical uncertainty.

Table 4 .
Bayes factors for the comparisons of the distance posteriors of the various SN pairs.