Issue 
A&A
Volume 682, February 2024



Article Number  A121  
Number of page(s)  21  
Section  Extragalactic astronomy  
DOI  https://doi.org/10.1051/00046361/202245497  
Published online  12 February 2024 
Intrinsic and extinction colour components in SNe Ia and the determination of R_{V}
^{1}
Université de Lyon, Université Lyon 1, CNRS/IN2P3, IP2I, 4 rue Enrico Fermi, 69622 Villeurbanne, France
email: g.smadja@ipnl.in2p3.fr
^{2}
Centre de Physique des Particules de Marseille, AixMarseille Université, CNRS/IN2P3, 163 avenue de LuminyCase 902, 13288 Marseille Cedex 09, France
^{3}
MaxPlanckInstitut fûr Astrophysik, KarlSchwartzschildStr. 1, 85748 Garching, Germany
^{4}
Princeton University, Department of Astrophysics, 4 Ivy Lane, Princeton, NJ 08544, USA
^{5}
Tsinghua Center for Astrophysics, MongManWai building, Tsinghua University, Beijing 100084, PR China
Received:
18
November
2022
Accepted:
25
July
2023
Context. The colour fluctuations of type Ia supernovae (SNe Ia) include intrinsic and extrinsic components, which both contribute to the observed variability. Previous works proposed a statistical separation of these two contributions, but the individual intrinsic colour contributions of each SN Ia were not extracted. In addition, a large uncertainty remains on the value of the parameter R_{V}, which characterises the dust extinction formula.
Aims. Leveraging the known parameterisation of the extinction formula for dust in our Galaxy, and applying it to the host galaxy of SNe Ia, we propose a new method of separation –valid for each SN– using the correlations between colour fluctuations. This also allows us to derive a wellconstrained value of the extinction parameter R_{V} with different, possibly smaller systematic errors. We also define a threedimensional space of intrinsic colour fluctuations.
Methods. The key ingredients in this attempt at separating the intrinsic and extinction colour components for each SN –and subsequently measuring R_{V}– are the assumption of a linearized dependence of magnitude on the extinction component of colour, a onedimensional extraintrinsic colour space (in addition to Ca II H&Kλ3945 and Si IIλ4131 contributions) over four independent colours, and the absence of correlation between the intrinsic and extrinsic variabilities.
Results. We show that a consistent solution is found under the previous assumptions, but the observed systematic trends point to a (small) inadequacy of the extinction formula. Once corrected, all systematic extinction effects can be cancelled by choosing a single scaling of the extinction colour component as well as an appropriate value of R_{V} = 2.181 ± 0.117. The observed colours are described within an accuracy of 0.025 mag. The resulting magnitude variability is 0.13 over all UBVRI bandpasses, and this fluctuation is shown to be independent of the bandpass to within 0.02 mag.
Key words: instrumentation: spectrographs / supernovae: general / dust / extinction
© The Authors 2024
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.
1. Introduction
Type Ia supernovae (SNe Ia) have been used for a long time as standardisable candles in the extraction of cosmological parameters (Perlmutter et al. 1999; Riess et al. 1998), and the largest correction in the standardisation scheme is related to the extinction by the host galaxy. However, a significant uncertainty remains on this reddening correction, which affects the potential of SNe Ia for the accurate determination of cosmological parameters. Here, we try to answer some remaining questions related to the use of the semiempirical extinction formulae derived for our galaxy by Fitzpatrick (1999). A widely employed standardisation method using SNe Ia is provided by SALT2 (Guy et al. 2007; Betoule et al. 2014). This latter incorporates purely empirical correlation between colour and magnitude with the general form m_{B} = M_{0} + αx_{1} + βc, where x_{1} characterises the shape of the light curve, and c ∼ (m_{B} − m_{V}) the colour of the supernova (SN). The parameters α and β are tuned so as to improve the standardisation of the SNe Ia in a given sample. This description is simple, but explicitly abstains from separating the intrinsic and reddening contributions to the observed colour. The reddening is caused by a variable mix of gaseous contributions (H_{2}, H, CH, etc.) and dust grains of variable sizes. Such a complex mixture would not be expected a priori to lead to a simple extinction law, but detailed work by several authors extending over 50 years (Rieke & Lebofsky 1985; Cardelli et al. 1989, hereafter CCM; Fitzpatrick 1999, hereafter FM99) led to universal formulae –depending on a single parameter R_{V}– related to the cross sections of the diffusion centres (molecules or dust grains) for photons. R_{V} is operationally defined as the ratio of the extinction in the V bandpass A_{V} to the reddening indicator E′∼E(B − V), either measured or adjusted by a fit to the observations.
A detailed study by Schlafly et al. (2016, 2017) confirmed a directional as well as distance variability of R_{V} in our galaxy (from 3 to 3.75), and provided confirmation of the results of FM99. All attempts to determine the parameter R_{V} from SNe assume that the extinction formulae derived for our galaxy apply directly to all host galaxies, leaving the possibility of different values of R_{V}, as the averaging which occurs when considering a sample of host galaxies differs from that performed when deriving the extinction formulae from stars of the Milky Way. The same universal extinction formula is assumed to be valid for the host galaxies of the SNe Ia, but previous investigations actually relied on a combination of photometric measurements with slit spectroscopy, which lacked the spectrophotometric information provided by the Nearby Supernova Factory (SNfactory, Aldering et al. 2020) collaboration. The contribution of the intrinsic variability was not directly monitored, and the applicability of the extinction formula to other galaxies could not be checked accurately. All the groups restrict their extinction analysis to redshifts larger than 0.01 in order to avoid a significant contribution of peculiar velocities.
2. Previous investigations
There have been many attempts to separate the intrinsic and extinction colour components of SNe, and we describe a few of them here. In many instances, an a priori distribution of the extinction is assumed; for instance an exponential distribution with a cutoff. This introduces a bias in the determination of the extinction of each SN, as the reddening involves an average over many different geometrical configurations of the SN with respect to its host. The exponential behaviour is a simple but unrealistic assumption. In the present work, an extinction colour and an intrinsic colour are extracted for each SN. The intrinsic colour contribution to E(B − V) is frequently derived from the lightcurve shape, whether it is the SALT2 x_{1} or the stretch parameter. In our case, an extra intrinsic component is introduced (for each SN).
Lira et al. (1998) and later Phillips et al. (1999) noticed that the colour evolution of all SNe Ia in the 30–90 days post Vmaximum is universal, which allowed them –by selecting SNe in E or S0 galaxies (without dust)– to derive a universal intrinsic SN Ia colour E(B − V)_{0} at any reference date in the interval from 30 to 90 days; the magnitude error quoted is 0.05. The galaxy reddening E(B − V)_{tail} at late epochs is then found by subtracting the intrinsic colour from the observed colour. After selecting a sample of SNe with low extinction (E(B − V)_{tail} < 0.06), a timedependent correction provided by the light curve allows them to evaluate the intrinsic value E(B − V)_{max, 0} at maximum B luminosity.
Another photometric technique was used by Wang et al. (2003), who analysed the colour–magnitude diagrams using the data from Hamuy et al. (1996) and Riess et al. (1998). The B magnitude varies linearly as a function of colour in the 10–30 day period post B maximum. The scatter of the residuals from the straight line is typically 0.05. It is claimed that extinction does not affect the shape of these diagrams (an approximation, as the spectra of the SNe evolve). The magnitude of the SN Ia at colour B − V = 0.6 is used as a reference value for the standardisation. To obtain a ‘reference’ magnitude at maximum, the linear behaviour is extended to B − V = 0, although some SNe Ia show a ‘bump’ departing from linearity. In the absence of a such a bump, the comparison with the observed allows a value of the reddening E(B − V) to be derived, which has a dispersion of about 0.1 mag with respect to the determination of Phillips et al. (1999) for the same SN Ia.
Jha et al. (2007) also relied on photometric bandpasses and K corrections, but took into account the colour evolution along the light curve of the SNe using the MLCS2k2 software (Multicolor Light Curve Shapes). A preliminary fit extracts the intrinsic Gaussian fluctuation (with σ = 0.049, and an average value of zero) and the exponential reddening distribution at date +35 days (with a colour decay constant of τ_{E} = 0.138). In each UBVRI bandpass, the same reddening distribution is used in the modelling of the light curve. Each of the five light curves is described by five parameters: date of peak B luminosity, distance modulus, time evolution, a unique reddening scale (A_{V}) (with the B − V shape found in the previous step), and R_{V}. A prior R_{V} = 3.1 with σ = 0.4 is assumed for the CCM law (the parameter R_{V} is needed to convert A_{V} into an extinction in different filters). These latter authors do not claim to describe the intrinsic distribution at maximum luminosity, nor to measure R_{V}.
The multicolour light curves of a sample of 80 SNe measured by different groups are analysed by Nobili & Goobar (2008). The light curve shape is characterised by its stretch, that is, the time dilatation factor s(sn) bringing the mean lightcurve width to the SN lightcurve width. The stretch is similar to the x_{1} variable mentioned earlier for SALT2. Given two bandpasses, X and Y, and stretch s: X − Y = b_{XY}(t)+a_{XY}(s − 1)+c_{XY}E(B − V). The data are Kcorrected for the changes of restframe bandpasses with redshifts, and E(B − V) is the average over all epochs of the B − V colour curve of each SN. R_{V} is first extracted from the measured coefficients c_{XY} found (e.g. c_{BV} = R_{B} − R_{V}). The best fit is found (with the CCM law) for R_{V} = 1.01 ± 0.25. When the excess dispersion of the residuals along the light curve after correction for the extinction is assigned to the intrinsic colour contribution, the same analysis finds R_{V} = 1.75 ± 0.27 using CCM. This difference strongly suggests that the intrinsic colour components should be taken into account.
The ‘twin’ SNe SN2014J and (unextincted) SN2011fe were compared in the optical and nearinfrared (NIR) range by Amanullah et al. (2014) from data collected by the HST (UV bands), at Mauna Kea (NIR), and NOT (optical). The comparison of the 12 light curves was performed using the extinction formula of FM99. Amanullah et al. (2014) find E(B − V) = 1.37 ± 0.03, R_{V} = 1.4 ± 0.1.
The Carnegie Supernova Project (Burns et al. 2014) emphasises the use of a ‘pseudostretch’ s_{BV}, which describes the time dependence of the B − V colour instead of the usual ‘stretch’, which relates directly to the light curve. This pseudostretch allows a prediction of the maximum B − V colour. The analysis proceeds to analyse the eight observed colours B − m_{i} in bandpasses (u, g, r, i, V, Y, J, H) of each SN as a sum c_{i}(sn) = P(s_{BV})+ΔA_{i}(E(B − V),R_{V}), where P(s_{BV}) is a secondorder polynomial, and ΔA_{i} is the predicted reddening for colour c_{i} from the extinction formula. From a set of 75 SNe, Burns et al. (2014) find R_{V} = 2.14 ± 0.16. With one intrinsic colour and one (measured) host extinction per SN, and without any assumption on the distribution of the extinction, their method and their result are relatively similar to ours. In Sect. 9, we suggest using the SALT2 variable x_{1} in the way s_{BV} is used by these authors.
Amanullah et al. (2015) include photometric measurements of seven SNe Ia between 0.2 and 2 microns from date −10 to +50. A mean spectral template is used at each date to evaluate the extinction correction in all filters with a standard CCM extinction law where E(B − V) and R_{V} are parameters to be adjusted to the observed colour. Amanullah et al. (2015) assume an intrinsic U − V dispersion of 0.1 mag (compatible with our own spectrophotometric measurements), and a larger dispersion of 0.3 mag for wideband filters at lower wavelengths. The different intrinsic colours have a mean value of zero and are supposed to be uncorrelated at any given date, but fully correlated at different phases. The light curves are then fitted in all filters as in Amanullah et al. (2014) to extract R_{V} and E(B − V) for each SN. The values of R_{V} for the different SNe vary from 1.4 ± 0.1 to 3.8 ± 1.5.
Mandel et al. (2017) analysed a sample of 250 nearby SNe Ia (0.01 < z < 0.10). They adapted the SALT2 formalism (Betoule et al. 2014) at maximum B light to include an additional intrinsic contribution of the intrinsic colour c^{int} to the B magnitude as well as an explicit reddening R_{B}E_{s}. The intrinsic B band absolute magnitude is given by with a contribution of the intrinsic colour c^{int} to the B magnitude in addition to the reddening R_{B}E_{s}. The intrinsic lightcurve shape parameter x, similar to SALT2 x_{1}, is assumed to have a normal distribution x^{int} ≈ 𝒩(−0.40, σ = 1.2), the intrinsic colour is with , and the extra (grey) magnitude dispersion ϵ ≈ 𝒩(0.,σ = 0.1) is introduced. The ‘grey’ fluctuation is the dispersion in the magnitude of SNe that have identical or quasiidentical spectral shapes after reddening correction; it is bound to be achromatic by definition, and its origin is not understood. The intrinsic parameters of SN s are thus . The host galaxy reddening E_{s} is positive, with an exponential dependence on colour exp(−E_{s}/τ), as in Jha et al. (2007), with a magnitude scale defined by R_{B}. The model is adjusted to the outcome of a SALT2 fit (m_{B}, c, x_{1}); the Hubble residuals are 0.16 mag. The 11 parameters of the model are obtained by a hierarchical Bayesian method. The observed values, which include reddening and distance modulus, are . The main results of the fit in Mandel et al. (2017) are: , σ^{int} = 0.100 (intrinsic magnitude dispersion), (intrinsic colour dispersion), R_{B} = 3.73 ± 0.31, τ_{E} = 0.069 (extinction decay rate). The dependence of the magnitude with intrinsic colour is observed with β^{int} = 2.25 ± 0.25. The dispersion of the intrinsic B − V colour found by these latter authors is four times larger than our result.
Thorp et al. (2021) used photometric data to investigate the potential dependence of the reddening as a function of host mass. Their model includes the Gaussian grey magnitude fluctuation as well as a spectral variability function (largely spectral lines) with an amplitude scale as a parameter. The R_{V} values of individual SNe are drawn from a truncated Gaussian distribution (R_{V} > 0.5). The Hubble residuals are improved from 0.16 (SALT2) to 0.12. Thorp et al. (2021) finds no significant difference in R_{V} between the large galactic mass and lowmass samples. The full sample averages to R_{V} = 2.70 ± 0.25.
As in Jha et al. (2007), Brout & Scolnic (2021) describe the observed colour as a sum of intrinsic and extinction colours, with a Gaussian distribution of the intrinsic colour, and an exponential shape of the dust. Intrinsic and extinction colour components are added into the observed colour. Their model also includes a Gaussian R_{V} distribution. When splitting their sample into lowmass and highmass galactic hosts, Brout & Scolnic (2021) find two different values of R_{V}, with R_{V} = 2.0 ± 0.25 (low mass) and R_{V} = 3.0 ± 0.4 (high mass).
A recent study by Wojtak et al. (2023) also investigates the impact of the host mass on the standardisation of SNe Ia within a Bayesian approach using multifilter light curves of SNe Ia. These authors introduce a shape parameter for the extinction distribution rather than an exponential with a cutoff and find strong evidence for two populations of SNe. The two populations, with probabilities 0.38 and 0.62, differ in their intrinsic colour and their reddening distribution. Wojtak et al. (2023) introduce an (unaccounted) intrinsic magnitude dispersion of 0.11 mag. The authors also define two intrinsic colour dispersions for the two populations: and . Both these values are larger than our own determination, which includes more spectroscopic information.
The only spectrophotometric data available are provided by the SNfactory with the SNIFS instrument. A first attempt by Chotard et al. (2011, hereafter C11), analysed the magnitudes of 76 SNe Ia at B maximal luminosity in the five synthetic tophat filters UBVRI. The magnitudes were corrected for the Si II and Ca II H&K equivalent widths, and the extinction correction was then derived from the CCM law with scale A_{V}. A value of R_{V} = 2.8 ± 0.30 was obtained. The extinction scale A_{V} was a free parameter adjusted to the data for each SN in this latter study. The magnitude uncertainties were estimated from the measurement errors (dominated by the flux calibration uncertainty) without including the (unaccounted) grey fluctuation of 0.11 mag. Colours are better measured than magnitudes as some systematic uncertainties cancel out, such as the absolute flux calibration uncertainty of 0.03 mag, and (largely) the error on the redshift (measurement and peculiar velocity). The present work is an extension of this latter investigation, but with an increased sample and several improvements.
Huang et al. (2017) uses the twin SNe SN2012cu (with a significant reddening) and SN2011fe (without extinction) measured by the SNfactory to derive the parameter R_{V} by comparing their magnitudes at different epochs. Although the spectra are very similar after extinction correction (the two SNe are ‘twins’), some spectral bins near absorption lines differ. These bins are averaged by a Gaussian convolution, and deweighted so that the χ^{2} in these regions is unity (a ‘floor’ error of 0.03 is added); the residual fluctuations of the deweighted residuals are of the order of 0.01 mag, and the value of R_{V} is R_{V} = 2.952 ± 0.081 for SN2012cu. Léget et al. (2020) and Boone et al. (2021) used the SNfactory spectrophotometric data to show that the largest spectral features occur with a remarkable intrinsic variability in particular around Ca II H&K and Si II.
While Léget et al. (2020) singles out the correlation of spectral lines with magnitudes, Boone et al. (2021) introduced three intrinsic variables ξ_{i}, which provide an accurate description of the full spectra for 173 SNe Ia (after selection). The formula of FM99 is used to correct for reddening. As in Huang et al. (2017), the parameter R_{V} is found by deweighting the spectral ranges with large intrinsic variability by a factor 1/σ(λ) at each wavelength λ, with the result R_{V} = 2.40 ± 0.16 (the method is called Read Between The Lines). The amplitude of the residual intrinsic variability of magnitudes away from spectral lines after reddening correction is generally of the order of 0.08 mag, but it is as low as 0.02 mag between 6600 and 7200 Å.
A large uncertainty remains today on the (mean) value of the parameter R_{V}, which plays a substantial role in the standardisation of SNe. Part of the discrepancy arises from underestimated systematic measurement errors on the path from detector to magnitude: the absolute calibration, the bandpass wavelength dependence, the redshift accuracy, the peculiar velocities, and the corrections for galactic reddening (host and ours). Another source of variability lies in different definitions of R_{V} as a parameter of an extinction formula, or as A_{V}/E(B − V), where A_{V} and B − V depend on the spectral distribution of the stars used (HD stars and SNe) and on the restframe bandpasses selected. Peculiar velocities are part of the redshift accuracy and influence the rest frame bandpasses as well as the absolute magnitude.
The spectrophotometric quality of the SNfactory data (Aldering et al. 2020) helps alleviate some of these uncertainties, and the difference in R_{V} between Huang et al. (2017) and Boone et al. (2021) –who average over 173 SNe– as well as the results of Amanullah et al. (2015) –who use a somewhat less homogeneous sample– may well point to an actual physical distribution of the value of R_{V}. Here, to quantify the smooth spectral evolution of the reddening, we favour the use of colours rather than magnitudes, as the latter exhibit an ‘unaccounted’ grey fluctuation of 0.08–0.10 mag from one SN to another, and the evaluation of the ratio of the correlated to uncorrelated colour variabilities, including the intrinsic contribution, is not straightforward.
We assume in the present work that the reddening of SNe by the host galaxy can be factored out and consider the observed colour fluctuations (from one SN to another) as a sum of intrinsic and extinction contributions, as done previously by Jha et al. (2007), Mandel et al. (2017) or Brout & Scolnic (2021). As Si IIλ6355 Å and Ca II H&Kλ3945 Å are dominant components of the spectral intrinsic variability, we first consider their contribution –as in C11–, but substituting Si IIλ4131 Å for Si IIλ6355 Å as the dependence of the B magnitude on the equivalent width of Si IIλ6355 Å is quadratic. We shall assume in addition that there is a single dominant source of intrinsic variability on top of the Ca II H&K and Si II absorption lines, so that the intrinsic variation of the different ‘colours’ is linked by three ‘intrinsic’ couplings. For each SN, we can then extract the extinction and intrinsic colour components using the extinction formula as leverage, and we revisit the determination of R_{V} using the correlations of these colour variations.
There are three advantages to using the particular data sample chosen here: the highquality spectrophotometric measurement, with wellmeasured spectral lines and subpercent errors on colours, and the sizeable homogeneous sample of 165 SNe, which allows us to use small colour fluctuations with respect to a mean template rather than the absolute values of these colours.
3. Sample of supernovæ
The data used were gathered by the SNfactory collaboration using their SuperNova Integral Field Spectrograph (SNIFS, Lantz et al. 2004), an automated instrument optimised for the observation of point sources on a structured background, with a spectral resolution of 0.25–0.30 nm and a good efficiency from 330 to 900 nm. SNIFS consists of a multifilter photometric channel used to monitor the transmission in nonphotometric nights and provide an image for guiding, a lenslet integral field spectrograph covering a field of view of with a grid of 15 × 15 spaxels, and an internal calibration unit (continuum and arc lamps). A more complete description of SNIFS, its operation, and the data processing can be found in Aldering et al. (2006); updated in Scalzo et al. (2010). The measured wavelength range in the SN restframe extends from 3300 to 8400 Å. It is divided into five logarithmically distributed tophat filters UBVRI, with central wavelengths at 363.9, 438.6, 528.7, 637.4, and 768.4 nm. The light curves are reconstructed in these five synthetic restframe bandpasses, and the present study is limited to the spectrum closest in time (within a window of ±2.5 days) to the maximal B luminosity, as found from a SALT2 fit (Betoule et al. 2014) to the light curves. The redshift range covered by our measurements extends from z = 0.02 to z = 0.11. Although there are only four independent colours, we consider all pairs of filters in the present analysis, which leads to ten possible filter combinations. The initial sample consists of 172 SNe Ia, which were selected as those yielding an acceptable SALT2 fit (Betoule et al. 2014), allowing us to define the date of maximum light: at least 5 nights of observations, nMAD^{1} of residuals < 0.12 mag, and a satisfactory phase coverage: at least four epochs from −10 to +35 days from maximum light, with one epoch between −10 and +7 days, and one between +7 and +20 days. In addition, we required in this work one spectrum within 2.5 days of maximum light. We eliminated seven SNe Ia where the subtraction of the host galaxy signal in our version of the data processing left a brightness gradient in either the B channel of the spectrograph (up to a wavelength of 5000 Å) or the R channel (above 5000 Å) over the 225 microlenses larger than 0.05 mag. The analysis is performed on the remaining 165 SNe Ia from the SNfactory. This sample is almost the same as that described in Aldering et al. (2020). As explained in the following section, the SALT2 magnitudes for our five synthetic bandpasses are not used; only the shape of the light curve near maximum light is used, as provided by the SALT2 model after the determination of x_{0} and x_{1} from the data.
Magnitudes at maximum light and errors. We require the magnitudes produced by the SALT2 fit in each bandpass to be consistent with the data near peak. From each restframe spectrum within a ±2.5 day window around the Bband maximum light, the distancecorrected magnitudes (brought back to an arbitrary common redshift of 0.05) are derived by integrating the photon count in each of the five UBVRI tophat bandpasses described in Table 1. The index T in Table 1 is used to distinguish the bandpasses from Bessell filters, and is dropped later in the text. The magnitude in each filter is rescaled to an interpolation of the actual observations as described in C11. The SALT2 fit extracts the date of the B maximum, the values of x_{0} (absolute magnitude scale) and x_{1} –which relate to the light curve–, and the colour c, which is very close to B − V at maximum B light. These parameters are then used to integrate the SALT2 spectral templates over our bandpasses, and the corresponding SALT2 light curves are reconstructed at the observation dates, for each SN, according to the values of , x_{1} and c parameters found in the lightcurve fit.
SNfactory tophat filter bandpasses (restframe).
In each bandpass, the shape of the SALT light curve near maximal B luminosity is then used as an interpolating curve in each filter and scaled as described in Eqs. (1) and (2) to fit our synthetic photometry from spectra within 5 days of maximum luminosity, as in C11. The same integration over the bandpass is applied to the SALT2 spectral template so that the corresponding SALT2 light curves are reconstructed at the observation dates, for each SN, according to the values of x_{0}, x_{1}, and c parameters found in the light curve fit. The SALT magnitudes reconstructed from the light curve in each bandpass F are then shifted at all dates p according to Eq. (2) so as to match the observed distancecorrected magnitude m_{F, p}. The shift ϵ_{F} is found by minimising
The values of ϵ_{F} are typically about 0.01 mag within the 2.5 day window around maximum light t_{max}, with a spread of 0.02 to 0.05 mag around the SALT2 light curve for UBVR, and somewhat larger for I. The rescaled magnitudes m_{F} at maximum B luminosity are then obtained by shifting the SALT maximum by this (averaged) ϵ_{F}:
Our magnitudes are then largely independent of the SALT colour law and of the SALT lightcurve parameters. Given the measurement errors of the magnitudes for the different filters at each date, and the covariance of the measurement errors between the different bandpasses (mostly a grey absolute calibration factor and a redshift uncertainty from peculiar velocities), we can then estimate the covariance matrix of the magnitudes in the different filters at maximum B luminosity. The rescaled SALT2 magnitudes at maximum are corrected for the extinction in our galaxy.
4. Methodology
The data analysis presented here starts from the processing described in Aldering et al. (2020). It is followed by an evaluation of the magnitudes at maximum luminosity in each bandpass, and a measurement of the equivalent widths of the significant spectral features as performed in C11 and Chotard (2011). The first step of the analysis, as in the earlier work of C11, is to subtract in Sect. 4.1 the intrinsic colour variability correlated to Si IIλ4131 Å and Ca II H&Kλ3945 Å signals using their equivalent widths ew^{Si} and ew^{Ca}. What is referred to as ‘colour’ in our study –except in this first step– is actually the ‘colour difference’ between the measured colour and the average value of the same colour over our sample of 165 SNe after subtraction of the contribution of these two lines. It differs from the standard colour excess, measured with respect to a sample of SNe without extinction by a constant, and can be negative. This constant offset (different for each filter) has no impact on any of our results. Each colour difference e_{i} is the sum of two components, an intrinsic part I_{i} and an extinction part X_{i}. In our method, X_{i} is also an extinction variation with respect to its average value and may assume negative values. For each filter, it also differs from the full extinction by an (irrelevant) constant. We then express the ten colour differences observed (after correction for ew^{Si} and ew^{Ca}) in terms of the following relations (without any summation of the indices):
and for each colour i:
The coefficients γ_{ij} and δ_{ij} are common to all SNe, labelled by index n; ϵ_{i} is the residual associated with the description of the colour difference e_{i} within the present model. While γ_{ij} relates intrinsic colours and does not depend formally on R_{V}, its determination will be weakly dependent on the value of R_{V} assumed. The previous relations do not involve convolutions with measurement errors, as these ‘colour’ errors are small enough with respect to the physical quantities I_{i} and X_{i} to be ignored at this stage.
Our notations are summarised in Table 2. The colours are numbered from zero to 9 in the sequence (U − B, U − V, U − R, U − I, B − V, B − R, B − I, V − R, V − I, R − I), and the differences in their measurement accuracies is ignored. The coefficients δ_{ij}(R_{v}) relating the extinction colour components X_{i} are derived from the extinction formula as described in Appendix A, while the 45γ_{ij}, which relate the intrinsic colour variations I_{i}, can be expressed as a function of three of them, selected (arbitrarily) to be γ_{10}, γ_{20}, γ_{30}. Thus, the proposed model is characterised at this stage by four parameters: R_{V}, and the three intrinsic coefficients γ_{10}, γ_{20}, γ_{30}. However, below, we find it necessary to introduce four corrections to the extinction formula in the bandpasses BVRI, which results in a total of eight parameters. For presentation purposes, an extra parameter (not required by our analysis, where it is unity) will be introduced for the extinction scale. For each SN, two coordinates are derived from the data using the previous parameters: X_{0}(n) and I_{0}(n). We call them coordinates rather than parameters as they are computed with a fixed algorithm from the measured spectra, without any tuning, once the parameters involved in our modelling of the full sample have been found. The extinction analysis proceeds in eight stages, with the successive application of different fits. To assess the validity of the model used, we have chosen to control each successive stage and apply judgement, rather than combine the different χ^{2}, which occur along the path into a single χ^{2} or a likelihood:
Definitions.
1. We subtract the contributions of the silicon and calcium equivalent widths ew^{Si} and ew^{Ca} to the colour fluctuations to ensure they are small in Sect. 4.1, allowing a linear approximation.
2. The extinctions in filter F derived from the extinction formula 11 are rescaled to A_{F} (Sect. 4.2.2) to ensure mathematical consistency. The colour coefficients δ_{ij} of Eq. (3) are invariant in this rescaling.
3. For a given set of γ_{ij}, the intrinsic (I_{i}) and extinction (X_{i}) colour components are found by a least square fit to the observed colours (Sect. 5.2).
4. The values of the intrinsic colour correlation coefficients γ_{ij} are found by iteration, and these are tuned in Sect. 5.3 by minimising the scatter of the output ratios of the 10 × 9/2 = 45 pairs of colourratios found from the colour fit in Eq. (29).
5. The extinction formula (FM99) is modified so as to suppress any correlation between the residuals ϵ_{i} derived from the fit to the extinction colour. This is achieved in Sect. 5.4 by minimising . The consistency rescaling is then applied to the corrected coefficients.
6. Iterations are performed over steps 4 and 5 on the values of the filter extinction corrections and the intrinsic coefficients γ_{ij}, until neither nor the scatter can be improved by a change of 10^{−3} of the extinction corrections or the coefficients (γ_{10}, γ_{20}, γ_{30}).
7. The previous steps do not depend on an overall scale factor applied to the (filtercolour) extinction coefficients δ_{Fi}. For each value of R_{V}, this scale factor is unity within our model, but it is introduced in almost all extinction corrections under the name of A_{V}, where the extinction at wavelength λ is A_{V}ϕ(λ), with the function ϕ standing for the extinction formula. To compare with models that do not enforce the normalisation of the extinction, we introduce in Sect. 7 a (illegitimate) scale parameter s multiplying the extinction magnitude A_{F}, so that the extinction in bandpass i will now be s δ_{i4} X_{4}, where colour 4 ≡ B − V (with our tophat filters).
8. It is then observed in Sect. 7 that there is a unique value of R_{V} for which the extinctioncorrected magnitudes in the other filters UBRI do not depend on the extinction colour. This provides our determination of R_{V}.
This derivation of R_{V} relies on the suppression of the (wellmeasured) colour dependencies on R_{V}, rather than on the fit of magnitudes involving a χ^{2} (subject to the ‘grey’ fluctuation of 0.08–0.12 mag and other correlated effects, such as modelling and calibration). There is no attempt at any stage in the previous sequence to minimise the colour residuals or the residuals of magnitudes in the Hubble diagram.
4.1. Calcium and silicon contributions
Nugent et al. (1995) showed that the impact of Ca II H&K and Si II absorption lines can be characterised by their equivalent widths and that they correlate strongly with magnitudes, colours, and lightcurve shapes. This section addresses a mathematical issue: how to suppress the part of the colour variation that depends on ew^{Si} and ew^{Ca} in the presence of a correlation between these two widths. C11, and other earlier works referenced therein, observed that there is a strong correlation between the equivalent widths ew^{Si} of Si IIλ4131 Å, ew^{Ca} of Ca II H&Kλ3945 Å, and the different colours. The equivalent widths ew^{Si} and ew^{Ca} are measured as in Bronder et al. (2008), using an algorithm described in C11, and in more detail in Chotard (2011): following Nugent et al. (1995), a spectral range centred on the spectral line is defined from λ_{1} to λ_{2}, and a reference level f_{c}(λ) is defined by a straight line connecting the two spectral values associated to λ_{1} and λ_{2}. The actual spectrum observed is s(λ). The equivalent width is
The errors on the equivalent widths are derived using a simulation that takes into account the photon statistics as well as the algorithm used to select the boundaries in the integration over the line width. As these two Ca and Si spectral features have a strong colour contribution over the whole spectral range, it is advantageous to treat them separately from the extraintrinsic colour and take advantage of the detailed spectral information. The analysis described in the following sections refers to the intrinsic colour component remaining after subtraction of the effect of Si IIλ4131 Å and Ca II H&Kλ3945 Å on the ten different colours^{2}. The small (formal) variations of the ‘observed’ restframe colour around its mean value are assumed to depend linearly on the variations (δew^{Si},δew^{Ca}) of the equivalent widths (with other contributions treated as noise).
The subtraction of the ew^{Si} and ew^{Ca} contribution to the colour variability is performed sequentially, taking into account their correlation. Let the restframe colour be a function of the two equivalent widths with and . The small fluctuations (from SN to SN) are defined as δ ew = ew − ⟨ew⟩ (where ⟨ ⋅ ⟩ denotes the averaging over all SNe). The associated change in ‘colour difference’ is given by the partial derivatives:
As we are only interested in the derivatives, the subtraction of the averages is unnecessary. The coefficients σ_{i} and κ_{i} are not directly observed from the data as a consequence of the correlation between the two equivalent widths, which needs to be removed, and can be described by two linear relations:
The coefficients a_{Si/Ca} and a_{Ca/Si} are not inverse of each other as a consequence of measurement errors and uncorrelated physical fluctuation. They are found directly from the observed equivalent widths by two different linear fits to the observed distribution in the (ew^{Si},ew^{Ca}) plane. Two χ^{2} are minimised:
The dependence of ew^{Si} on ew^{Ca} gives a_{Si/Ca} = 0.0379 ± 0.0771, and the dependence of ew^{Ca} on ew^{Si} gives a_{Ca/Si} = 0.513 ± 0.286. The observed values are not weighted by the measurement errors, which do not include the intrinsic variability of the equivalent widths. The correlation a_{Si/Ca} of ew^{Si} as a function of ew^{Ca} is compatible with zero, and ignoring it would not have a significant impact on the corrections. We nevertheless take it into account to ensure that when the procedure described below is applied, there is no residual dependence whatsoever of colours on ew^{Si} or ew^{Ca}. We stress that we do not claim that either a_{Si/Ca} or a_{Ca/Si} represents the physical correlation between the equivalent widths, which is described by a single number . Extracting would require a full understanding of the convolution of measurement errors and variability, which we do not actually need.
The large difference between the two coefficients reflects the large difference in the equivalent widths of the two lines. The observed dependencies of colours on the variations of Ca and Si equivalent widths are the total derivatives of the restframe colours with respect to ew^{Si} and ew^{Ca}, and , which can be expressed in terms of the partial derivatives and :
The total derivatives K^{Si} and K^{Ca} are obtained from σ_{i} and κ_{i} in Eq. (9), once the coefficients a_{Ca/Si} and a_{Si/Ca} are known. Inversely, the coefficients σ_{i} and κ_{i} can then be derived in Eq. (10) from the observed dependencies of the colours on the equivalent widths, as seen in Figs. 1 and 2:
Fig. 1. Correlations between the first five colours , , , , , (mag) and ew^{Si} (Å). The slopes of the straight line fits (red lines) give the value of K^{Si} directly. The correlations of the colours can be derived from the first four. 
Fig. 2. Correlations between the same colours as Fig. 1, corrected for the ew^{Si} correlation (mag), as a function of ew^{Ca} (Å). The slopes of the straight line fits (red lines) give the value of K′^{Ca}. 
The denominator is zero when a_{Ca/Si} = 1/a_{Si/Ca}, but as explained previously, this would only occur if there were no fluctuations in the correlation between the two equivalent widths, and the issue addressed in this section would vanish.
We now use Eq. (10) to correct for the dependence of colours on the equivalent widths in the data. If the ew^{Si}corrected colour is first applied, the observed dependence with ew^{Ca} will be . The values of K^{Si} and K′^{Ca} are directly obtained from the correlation between the colours and the equivalent widths observed in Figs. 1 and 2. The ew^{Si} and ew^{Ca} contributions to each colour are then subtracted to yield . In all the following steps of the analysis, we use the colour differences e_{i} = c_{i} − ⟨c_{i}⟩, where ⟨ ⋅ ⟩ represents –as above– the mean over all the SNe Ia of our sample^{3}. To preserve the algebraic relations between the ten colours, the ew^{Si} and ew^{Ca} corrections are performed on the first four colours, U − B, U − V, U − R, and U − I, while the others are obtained by combining them. The coefficients K^{Si}, K′^{Ca}, σ, and κ are summarised in Table 3 for the first four colours; the other colours can then be obtained by subtraction from the appropriate pair.
Colour correction from ew^{Si} and ew^{Ca}.
As shown in Sect. 9, we found that the correction of colours for spectral features was necessary. The parameters σ_{i} and κ_{i} could be considered as additional parameters of the model implemented, though they are not tuned but measured; we show in Sect. 9 that the CaSi colour correction is closely related to the SALT variable x_{1}.
4.2. Evaluation of extinction coefficients δ
4.2.1. From extinction formula to coefficients δ
This section describes the contribution of the reddening to the colour variations. Our reddening indicator is the extinction part X_{4} (defined in Eqs. (3) and (4)) of e_{4} = E(B − V)−⟨E(B − V)⟩ after the ew^{Ca} and ew^{Si} corrections have been applied; its mean ⟨X_{i}⟩ is zero, and it can have a positive or negative sign, but its variations are the same as the standard one (up to the intrinsic contributions of ew^{Si}, ew^{Ca}, and I_{i}, which we do want to exclude from the extinction correction). The SNfactory tophat filters are defined in Table 1, and the extinction formula is (at first) assumed to be perfectly known, as given by FM99, once a value of R_{V} is chosen, which we assume until Sect. 7.
The reddening corrections and of Eq. (12) for the bandpass F are then evaluated for each value of R_{V} in a grid, and the analysis described in this section is repeated over the grid. The difficulty that we address is that the extinction indicator E_{n}(B − V) (with or without subtraction of the mean) is always tied to the choice of a wide bandpass, tophat (our case), or Bessell. The reddening corrections at each wavelength depend on wavelength and R_{V}, but also on the choice of wide band filter, and on the spectral shapes of the selected stars through the wideband integration. The restframe photoncount spectrum s_{n}(λ) for SN n is
With the colour difference definition E_{n}(B − V) = (B − V)_{n} − ⟨B − V⟩, the spectrum s_{0}(λ) at null colour difference is the mean SN Ia photoncount spectrum; the function ϕ(λ) is provided by Cardelli et al. (1989), O’Donnell (1994), Fitzpatrick & Massa (2005). This extinction indicator E_{n}(B − V), or ‘colour difference’, differs by a constant term from the standard E(B − V), which is referenced to the zero extinction value rather than the mean; the extinction indicator E_{n}(B − V) can be positive or negative. This constant shift (bandpassdependent) affects colours, but not their variations, with which we deal exclusively. Equation (11) is a slight extension of the standard extinction formula of FM99 on two counts: the extinction formula relies on a specific set of blue stars, not SNe, and the bandpasses are defined by Bessell filters, rather than tophat filters as described in Table 1. This issue is addressed in Sect. 4.2.2. Hereafter, we use the extinction component of E_{n}(B − V), that is, X_{4}(n) = X_{B − V}(n), as a substitute for the colour difference E_{n} in Eq. (11). It is shown in Appendix A that, in the linear approximation, the extinction correction in filter F is given for a SN with photoncount spectrum s(λ) by:
The coefficient is the exact derivative of the reddening correction (n) with respect to the extinction colour component X_{4} at X_{4} = 0. We show in Table A.1 that the relative error resulting from this approximation remains smaller than 0.01 over the full range of values of X_{4}. The reddening correction coefficient should actually be weighted by the spectrum of each SN according to Eq. (12), but as we only consider a linear reddening correction, the SN dependence would bring a negligible secondorder correction. The term X_{B − V}(n) = X_{4}(n) in Eq. (12) is the extinction content of the colour difference e_{4} = E_{n}(B − V) for SN n. As noted previously, the extinction colour X_{i}(n) can be positive or negative, but should be a better reddening indicator than E_{n}(B − V) as the intrinsic colour variation has been removed.
For each bandpass F = UBVRI, the extinction coefficient is obtained from Eq. (12). The associated colour differences are e_{i} = c_{i} − ⟨c_{i}⟩. The extinction coefficients (j ≠ 4) are then derived from the integrated extinction formula. The value of for another colour i can be derived from ; for example (with e_{3} ≡ U − I) is given by:
The colour–colour coefficients relate the (extinction) colour excesses of all colours; they are obtained from the two bandpasses defining each colour, from the analogue of Eq. (13); for example, for e_{2} ≡ U − R and e_{6} ≡ B − I, one has X_{2} = δ_{24}X_{4} = δ_{26}X_{6} with
The extinction correction for all the other colours is similar. The colour indicators used in this study is actually the extinction component X_{B − V} of B − V computed in Sect. 5.2 rather than the measured value e_{B − V}. For instance, the correction to e_{5} ≡ B − R using X_{6} ≡ B − I is
Using for instance the extinction part X_{4} of colour 4 ≡ B − V or the extinction part X_{6} of colour 6 ≡ B − I:
(We recall that has been brought to a common reference redshift).
4.2.2. Rescaling of extinction coefficients: from δ′ to δ
The extinction formulae, such as those of Cardelli et al. (1989) or FM99, have not been tuned for SN spectra, nor to our set of tophat bandpasses. In addition, the extinction law was established for a sample of mostly HD stars, where the extinction of each spectral bin is weighted by a spectrum that differs from the SNe Ia spectrum. It should not be assumed that the extinction formulae provided by Cardelli et al. (1989), O’Donnell (1994), Fitzpatrick & Massa (2005) –with different filters (Bessell) and different stars (HD)– apply directly to SNe Ia, given the presence of the wide band term E_{n}(B − V) or X_{4}(n) in Eq. (11). The minimal correction proposed here is a rescaling of the formula so as to ensure consistency. If we apply Eq. (12) to e_{4} = E_{n}(B − V), we find (or ). We must expect
If the extinction in the V bandpass is used –which is a frequent occurrence– a different consistency condition should be implemented: to ensure that A_{V} is the actual extinction in this bandpass. Equation (18) is sometimes written (e.g. Mandel et al. 2017; Jha et al. 2007) as R_{B} − R_{V} = 1. However, we find (for R_{V} = 2.25) with the previously defined tophat bandpasses , while with O’Donnell (1994) extinction parameters, . Most of this difference is due to the Bessell bandpasses assumed in the extinction formula, while we are using tophat filters in the present study. Indeed, when using Bessell filter weighting (Bessell & Murphy 2012), , which is closer to but still incompatible with unity. We show in Eq. (A.5) of the Appendix that the expression of the coefficient in Eq. (12) is the exact derivative of the reddening correction with respect to the extinction colour at X_{4} = X_{B − V} = 0, and in Table A.1 that the linear approximation is accurate to within 5 10^{−3} mag over the full range of extinctions. Another contribution to this mismatch may be the different spectra of stars used in FM99 and in the present work (SNe Ia). We assume that an acceptable correction ensuring the ‘consistency condition’ in Eq. (18) is an overall rescaling of the extinction formula at all wavelengths, so as to force R_{B4} − R_{V4} = 1. All the coefficients are rescaled:
The extinction coefficient of each filter is then increased by the factor . No such global rescaling was applied in C11 where an extinction parameter A_{V} was fit independently for each SN. As mentioned above, Eq. (11) is only an approximation. Even if one accepts that different molecular and grain compositions can be summarised by a single coefficient R_{V} (related to their cross section for light), its value could differ in the halo and in the disc (if there is one), and the observations of Schlafly et al. (2016) in our own galaxy suggest this possibility. The δ_{ij} coefficients of Eq. (4) verify that δ_{ij} = 1/δ_{ji}. The five values of the extinction coefficients δ_{F4}(R_{V}), which provide the extinction in filters for F = UBVRI from X_{4} = X_{B − V}, are given in Table 4. The extinction coefficients δ_{F}(R_{V}) are subject to a first correction (Dex_{F}), which cancels the correlations of the colour residuals with the extinction components as discussed in Sect. 5.4; they are then rescaled for algebraic consistency as described in Eq. (19). Finally, for arbitrary values of R_{V}, the suppression of the correlation of extinctioncorrected magnitudes with the extinction colour in the V bandpass requires an extra overall rescaling s, which is a function of the extinction parameter R_{V}, as explained in Sect. 7. Such a scaling destroys the consistency condition in Eq. (18), but is introduced in practice by almost all authors (Boone et al. 2021; Chotard et al. 2011; Amanullah et al. 2015; Saunders et al. 2018; Thorp et al. 2021) who substitute where A_{V} is arbitrary and ϕ(λ, R_{V}) stands for one of the reddening formula. The analysis of the extinction is performed over a grid of nine values of R_{V}: 1.95, 2.05, 2.15, 2.20, 2.25, 2.35, 2.50, 2.60, and 3.10. We show below that for our value of R_{V}, the scaling factor s is indeed compatible with unity, which confirms the consistency of the model, as well as the value of R_{V} obtained.
Extinction coefficients from FM99 and UBVRI corrections to the extinction formula for the case R_{V} = 2.25.
5. Extraction of intrinsic and extinction colour components
5.1. Intrinsic colour coefficients γ
As mentioned in Eq. (3), the colour difference e_{i}(n) of SN n (after correction for ew^{Si} and ew^{Ca}) is the sum of two components e_{i}(n) = I_{i}(n)+X_{i}(n)+ϵ_{i}(n), and the intrinsic component I_{i} of colour i is assumed to belong to a onedimensional space of intrinsic colour fluctuations, as suggested by Léget et al. (2020) and Boone et al. (2021). These are then interrelated by a set of coefficients γ_{ij} common to all supernovæ n as introduced in Eq. (4). As shown below, small corrections to the extinction formula are required. An initial value for the coefficients γ_{ij} can be derived directly from the requirement that intrinsic and extrinsic colours be uncorrelated, though this value can differ significantly from the result of the χ^{2} minimisation discussed later. For any pair of leading (i) and auxiliary (j) colour, two equations link X_{i}, X_{j} and I_{i}, I_{j}:
The initial values of γ_{i0} can be obtained from the previous equations, assuming the absence of correlations between the I_{i} and the X_{i}, and neglecting the contribution of the residuals ϵ_{i} and ϵ_{ij},
Here, δ_{ij} is the ratio of the extinction components X_{i}/X_{j}, as defined in Table 2 –and is not the Kronecker δ–, and the index n of the SN is dropped.
Combining Eqs. (23) (×δ_{10}) and (24), one gets . Similarly, ; dividing the two previous relations one obtains an initial value of γ_{10} (with large errors). More generally, an initial value for any pair of colours would be:
The initial values of the three coefficients γ_{10}, γ_{20}, and γ_{30} are respectively 0.671, 0.759, and 0.401. The statistical errors on these initial values are relatively large, and an improved determination is obtained from the iterative fits described in the following section. (We recall that all colour differences have zero average, and that intrinsic and extinction components are assumed to be uncorrelated.)
Given the relation δ_{ij} = 1/δ_{ji}, we find that (as expected) γ_{ij} = 1/γ_{ji}. The coefficients γ_{ij} can be expressed according to Eq. (4) in terms of three of them, which have been selected as γ_{10} = I_{(U − V)}/I_{(U − B)}, γ_{20} = I_{(U − R)}/I_{(U − I)}, and γ_{30} = I_{(U − I)}/I_{(U − B)}. The remaining 42 coefficients can be obtained from the following relations, which result from the algebraic relations between colours:
5.2. χ^{2} for intrinsic and extinction colour components
The model described by Eq. (20) should reproduce the restframe colours. We consider the colour difference for SN n: e_{i}(n) = c_{i}(n)−⟨c_{i}⟩. The two components (I_{i}(n),X_{i}(n)) can be found by requiring that all colours e_{j}(n) are properly described by Eqs. (20) and (21), and the coordinates I_{j}(n) = γ_{ji}I_{i}(n), X_{j}(n) = δ_{ji}X_{i}(n). Colour components X_{2}(n) and I_{2}(n) are, for example, obtained by minimising –for each SN n– the in Eq. (29) as a function of the ‘coordinates’ I_{2}(n) and X_{2}(n). The reddening correction A_{V}(n) in bandpass VA_{V}(n) = R_{V4}X_{4}(n), which varies from SN to SN, plays no role in the minimisation: the extinction correction is A_{F4}(n) = R_{F4}X_{4}(n) for filter F and SN n (=A_{Fi}(n) = R_{Fi}X_{i}(n) using another colour):
as a function of I_{2} and X_{2}, under the constraint
Equation (30) imposes the decorrelation of the extraintrinsic colour from extinction. The index n is omitted in Eq. (27). For each colour i and each SN, the residual ϵ_{i}(n) defined by Eq. (27) should not be confused with the magnitude shift at maximum light defined in Eq. (2). We are aware that Eq. (30) could be criticised under certain circumstances, such as in the presence of circumstellar matter with a composition differing from the average host galaxy. Such matter, as pointed out by Borkowski et al. (2009), Ferretti et al. (2017), might invalidate this assumption, but if the circumstellar dust is created by a nondegenerate companion, as suggested in one case by Nagao et al. (2017), it is not necessarily correlated to the SN Ia properties. While not rejecting the possibility, we wanted to explore the accuracy that could be reached within the assumption of decorrelation.
Each pair (I_{i}(n),X_{i}(n)) of components allows the prediction of all the other colours of the SN. When defining χ^{2}(n) in Eq. (29), each colour i enters with equal weight, which implicitly assumes that the errors on all colours are similar. The solution of Eq. (29) in terms of (X_{2}(n),Y_{2}(n)) is found by requiring :
with
The colours excesses e_{1}(n) and e_{2}(n) are measured, meaning that Y_{0} and Y_{1} are known. The intrinsic and extinction components I_{2}(n), X_{2}(n) are obtained by inverting Eq. (31). The Lagrange multiplier λ arises from the assumption of uncorrelated I_{i}, X_{i}, and is found from Eq. (30); it depends on the colour, but is the same for all SNe. We introduce the set Aux of auxiliary colours used to extract the intrinsic and extinction components of each colour. In the optimisation of the intrinsic couplings γ_{ij} below, the full set of ten colours cannot be used in view of the appearance of algebraic constraints, as explained in the following section. Table 5 presents the actual list of auxiliary colours used for each colour. There is an arbitrariness in the choice of the participating colours in Aux; we selected the one that yielded the smallest residuals among five trials.
For each SN n, this minimisation yields ten intrinsic colours I_{i}(n), and ten extinction colours X_{i}(n), once the intrinsic colour coefficients (γ_{01}, γ_{02}, γ_{03}) are known. Equation (29) does not take into account the correlation between colours, but as the solution gave acceptable residuals for all colours, we kept this simplified expression for χ^{2}. We verified that including the obvious algebraic correlations between the different colours does not change the solution, but it does add significant complications. As the mean of e_{i} is zero (by definition), the mean values of I_{i} and X_{i} obtained by linear combinations of e_{i} will also be zero.
5.3. Determination of intrinsic components and intrinsic couplings
The colour differences which we now use were subtracted from the contributions of Ca and Si, but leave room for an extraintrinsic component I_{i} of each colour, with all of them dependent on three parameters (γ_{10}, γ_{20}, γ_{30}) as required by the assumption of a single extraintrinsic contribution. Up to this stage, the γ_{ij} are obtained from their initial values in Eq. (25). The solution of Eq. (31) is only acceptable if the ratio of the solutions of Eq. (31) I_{i}(n) to I_{j}(n) for two colours i and j is equal to γ_{ij} for all SNe and all pairs of colours. Similarly, the ratio of X_{i}(n) to X_{j}(n) should be equal to δ_{ij}.
This is achieved by a sequence of iterations. At each iteration, and for each colour, the input coefficients are the parameters introduced in Eqs. (4) and (28), while the linear coefficients of the straight line fits in the plane (I_{i}, I_{j}) as in Fig. 3 are the output value of the intrinsic couplings. The output values of (γ_{10}, γ_{20}, γ_{30}) are reinjected into Eq. (29) until Fig. 3 is satisfactory. The scatter is seen to range from 0.002 mag to 0.01 mag depending on the colours, and the measured values of the ratios are within 0.01 of the corresponding set of , as obtained from (γ_{10}, γ_{20}, γ_{30}) using Eq. (26). The δ_{ij} coefficients are left unchanged in this sequence, as their values are fixed by the extinction formula (at this stage).
Fig. 3. Intrinsic and extinction ratios as defined in Eq. (31) (i = 2, 6, 9) compared to the input and calculated δ_{0i} for R_{V} = 2.25. The compatibility is at the level of 0.001 for the extinction, and the relative accuracy of 0.05 for the (smaller) intrinsic colour (worst case). The red lines are linear fits to the data, and γ^{out}, δ^{out} are the values of the linear coefficients. 
We then turn to Fig. 4, which compares the full set of input values in Eq. (31) in the secondtolast iteration to the output values, that is, of I_{i}/I_{j} from the linear fits in Fig. 3. The final values of the three coefficients (, , ) are found by minimising the scatter (RMS) of the output coefficients with respect to the input coefficients (horizontal scale). This last step brings in the constraint of consistency with Eq. (26) for the whole set of γ_{ij}, and modifies the result of the previous sequence of iterations; however, as seen in Fig. 3, the outcome is still satisfactory with respect to the convergence of the solution for each γ_{ij}
Fig. 4. Relation between input and output (⟨I_{i}/I_{j}⟩) values for the coefficients γ_{ij}. The overall agreement for all γ_{ij} (excluding those involving colours 5 ≡ (B − R) and 7 ≡ (V − R)) is within 0.01. 
The relation between and is shown in Fig. 4 at the last iteration. The linear coefficient is expected to be close to unity if convergence has been achieved, while we find 1.00069, with an RMS of 0.0104. This dispersion arises from a combination of measurement errors, SiCa correction errors, and modelling errors. Two precautions were implemented to extract the three intrinsic couplings from Eq. (31) and Fig. 4: we eliminated the auxiliary colours 4 ≡ B − V and 7 ≡ V − R that have very small intrinsic components (cf. Table 9) from the 10 × 9/2 = 45 γ_{ij} couplings, 28 remain after the elimination of colours 5 and 7. The 28 output values of the parameters at the last iteration are compared in Fig. 4 to the input values derived from (γ_{10}, γ_{20}, γ_{30}) using Eq. (26). The equality of input and output values was not forced into the optimisation process. The ratios δ_{ij} = X_{i}/X_{j} are also seen to have a small scatter of 0.002 mag (in Fig. 3), with averaged output values close to the input coefficients, supporting the model Eqs. (27)–(29).
Nevertheless, the RMS in Fig. 4 fails as an indicator when the same set of colours (auxiliary and main) is used to find the intrinsic content, as the algorithm defined by Eq. (31) yields . For each colour, we eliminated at least one extra auxiliary colour to avoid such duplicated sets. The choice of the auxiliary colours has an impact on the values found for the three intrinsic couplings. We selected the one that gave the smallest mean residual scatter (RMS) in Fig. 4. The scatter of in Fig. 4 is 0.0103, and the ratio is within 0.001 of unity. After the sequence of iterations described in Sect. 5.4, the optimal values in Fig. 4 for R_{V} = 2.20 are given in Table 6.
Coupling constants of the intrinsic colours.
The errors are derived at this stage from the spread of the results in four different choices of auxiliary colours and do not include the statistical error. The comparison of our intrinsic colours to previous results is postponed to Sect. 9, as the contribution of the equivalent widths ew^{Ca} and ew^{Si} must be reintroduced.
5.4. Modifications to the extinction formula
For all values of R_{V}, the colour residuals ϵ_{i} described in the previous section (after Eq. (29)) show a strong correlation with the extinction colour component when the FM99 extinction formula is used. This (unexpected) correlation is apparent in Fig. 5, and there is no value of R_{V} that suppresses this effect. As our subsequent determination of R_{V} relies on the absence of any correlation between the extinctioncorrected magnitudes and the extinction, it is crucial to suppress the correlation observed in Fig. 5. This suppression depends on the coefficients γ_{ij} and δ_{ij} in Eq. (32) but in no way on the scale of the reddening correction (nor on the rescaling applied in Eq. (19)).
Fig. 5. Correlation of the colour residuals (mag) for R_{V} = 2.25 with the extinction colour for colour differences e_{0}, e_{2}, e_{5}, e_{6}, e_{8}, e_{9} (no extinction formula correction). The strong correlations observed would forbid a safe evaluation of R_{V}. 
We can cancel the effects observed in Fig. 5 by introducing corrections Dex_{F} to the extinction coefficients R_{F4} of the five F = UBVRI bandpasses in Eq. (12). The minimised quantity is
The extinction corrections Dex_{F} in the five filters are applied to the values of from Eq. (12): and fitted for each value of R_{V} so as to minimise , and cancel the coefficients sl_{i} in Fig. 5. The result is given in Table 7. As there is an overall scale degeneracy of the extinctions in the colour analysis, we arbitrarily set the U bandpass correction to be 0.0675. At each value R_{V} of the grid, the corrections are evaluated. As these corrections (partially) compensate for the distortions introduced by the difference between the ‘optimal R_{V}’ and the value used in the grid, they tend to blow up when R_{V} is close to the limits of the range considered in this work, namely as R_{V} → 1.95 or R_{V} → 3.10. The offsets are applied in the second line of Table 4. The impact of these corrections on the correlation between the residuals and the extinction is shown in Fig. 6, and for all colours in Table 8. The large value of sl_{9} observed (e_{9} ≡ R − I) in Fig. 5 reflects the significant corrections to the extinction formula in bandpasses R and I, with opposite signs. The errors given in these figures are only indicative, as somewhat arbitrary ‘floor’ errors of 0.006 mag are imposed on the determination of the extinction colour.
Fig. 6. Correlation of residuals of colour differences (mag) (e_{0}, e_{2}, e_{5}, e_{6}, e_{8}, e_{9}) with extinction after correcting the extinction formula with R_{V} = 2.25. There are only four independent colours, meaning that the four extinction corrections allow us to cancel the correlations of all ten colours. 
UBVRI corrections (mag) to the extinction formula.
Residual (mag) dependence on the extinction colour R_{V} = 2.25 (after correction of the extinction formula).
The corrections to the extinction formula are extremely unlikely to be caused by the limitations of our model: The linear approximation involved is numerically correct to better than 0.005 mag in the whole range of e_{4} = E(B − V) considered. There could be a second intrinsic colour involved; that is, rather than onedimensional (beyond Ca and Si), as in Eq. (4), the intrinsic colour space could be twodimensional, with a second set of γ_{ij} coefficients, but intrinsic colours are uncorrelated (or are at best weakly correlated) to extinction, and this new intrinsic colour would have to be larger than the one already introduced to make up for the 0.10 mag discrepancies observed in Fig. 5. There is no room in the residuals for such a large extraintrinsic component. On the other hand, the need for an adaptation of the extinction formula to SNe and Tophat filters is shown in Sect. 4.2.2, where the presence of a rescaling factor is first introduced.
The numerical value of the corrections are dependent on the specific bandpasses involved, on the spectral template, and also possibly on the dust properties of our galaxy, which may differ from the averaged dust extinction. Figure 6 shows that, after correction, the remaining correlations of the residuals ϵ_{i} with extinction are negligible. This can be achieved regardless of the value of R_{V}. The largest of these effects amounts to a residual linear dependence of 0.006 for B − V. The coefficients of all other colours are smaller than 0.0036. The statistical errors on the offsets away from the extinction formula for our sample of SNe will be estimated from a simulation in Sect. 8.4.
It seems unnatural to constrain the γ_{ij} coefficients with , although they are (weakly) sensitive to the value of the extinction corrections. Our choice is instead to optimise the dispersion σ(γ) displayed in Fig. 4, as stated in Sect. 5.3. At each value of the parameters (Dex_{F} or γ_{i0}), the I_{i} and X_{i} are obtained from Eq. (31), and the γ_{i0} are retuned in order to minimise σ(γ) in Fig. 4. This sequence is repeated four or five times to reach the final values of the seven unknown parameters (four Dex_{F} and three γ_{i0}). Within each cycle (extinction corrections Dex_{F} or γ_{i0} fixed), convergence is reached when or σ(γ) increase whenever any one of the seven parameters is changed by 10^{−3}. As the implemented model is still imperfect, the minima of and σ(γ) are not perfectly compatible; the differences between the two minima are nevertheless approximately ten times smaller than the quoted errors.
6. Quality of the colour reconstruction
The quality of the reconstruction of the colour differences e_{i} from its intrinsic I_{i} and extinction X_{i} components is shown in Table 9 and in Fig. 7, where the distribution of the difference ϵ_{i} = c_{i} − (I_{i} + X_{i}) is given for six selected colours. The average colour residual over the ten colours is 0.0176. As the colour measurement errors do not exceed 0.006, this figure (together with Table 9) suggests that the dominant contributions arise from modelling error and from the subtraction of the ew^{Si} and ew^{Ca} contributions. Given the accuracy reached, the phase of the spectrum might also contribute.
Fig. 7. Residuals ϵ_{i} (mag) after the reconstruction of six colours (c_{0}, c_{2}, c_{5}, c_{6}, c_{8}, c_{9}). The full list of residuals is given in Table 8. 
Range of the intrinsic and extinction colours for observation, simulation (seed2), and residuals (mag).
7. Determination of R_{V}
7.1. Reddening correction of magnitudes
In Sect. 5.4, we ensure that the colour residuals are independent of the extinction component for each selected value of R_{V}. However, after correction for reddening, the magnitudes and the colours themselves may still depend on the extinction –even though the dependences in data and model must be the same after the corrections of Sect. 5.4. For each choice of R_{V}, the magnitude m_{F} found from the data at a reference redshift of 0.05 is corrected for its dependence on ew^{Si} and ew^{Ca} –as carried out in Sect. 4.1 for colours– in order to obtain . The dependence of on the intrinsic component I_{i} with the (measured) coefficient z_{F1} is then taken into account, and we allow for a single overall scaling factor s(R_{V}), which mimics the arbitrary A_{V} used in most analyses (it should be unity in our framework) to define the extinctioncorrected magnitude of each SN in the bandpass F using colour i:
Once the extinctions have been scaled according to Eq. (19), as explained in Sect. 4.2.2, the magnitude should not depend on X_{i}. The extinction scale s(R_{V}) in Eq. (34) is chosen so as to cancel the dependence of the corrected data magnitude in the V (tophat) bandpass, and we expect s(R_{V}) = 1 for the physical value of R_{V}. It is seen in Fig. 8 that this is true when R_{V} ≈ 2.25. The extinction in filter F is now A_{F} = s(R_{V})×R_{F3} × X_{3}. However, if we turn to the other filters, we see in Fig. 9 that even after the extra scaling factor s(λ) has been introduced in Eq. (34) to ensure , the derivatives of the other bandpasses are finite when R_{V} lies away from the 2.15 − 2.3 range. The dependence of as a function of R_{V} is shown for the bandpasses UBRI in Fig. 9. We find the extinction scale s to be s(R_{V} = 2.20) = 1.0269 and s(R_{V} = 2.25) = 1.0014 (close to our ‘effective’ value of R_{V}). The errors in Figs. 8–11 are strongly correlated over the range of wavelengths, and for different values of R_{V} as the same events are used. The error is evaluated from the dispersion of the scale factor in the ten samples of the simulation. From sample to sample, the scale parameter moves up and down in Fig. 8, but the smooth behaviour is preserved as a function of R_{V} as the same events are involved for a given sample.
Fig. 8. Extinction scale factor s as a function of R_{V} for dm_{V}/dX_{3} = 0. For each value of R_{V}, there is a value of the extinction scale factor s that cancels the dependence of the V magnitude on extinction. 
Fig. 9. for the UBVRI bandpasses as a function of R_{V}. Although the dependence of the V bandpasses as a function of extinction has been cancelled for each value of R_{V}, the extinction in other bandpasses is not adequately corrected. 
Fig. 10. d(U−I)^{corr}/dX^{3} as a function of R_{V}. The difference between U and I magnitudes (colour 3) is the most sensitive to the choice of R_{V}. The error bar is the dispersion found in the simulation. 
Fig. 11. UBVRI magnitudes corrected for I_{3}, X_{3} with R_{V} = 3.100, s = 0.73344. The derivative is compatible with zero in each bandpass, but the variation from U to I is quite significant, taking correlations into account. The red line is weighted by errors and the blue straight line is unweighted. 
Each bandpass provides a value of the parameter R_{V} that cancels the derivative . For UBRI, these values are respectively 2.283, 2.183, 2.137, and 2.240. We average the U and I measurements –which have a stronger R_{V} dependence– and obtain (with the error from the simulation):
The scale factor s for this value of R_{V} is s = 1.001 ± 0.045, which is additional independent confirmation of the result for R_{V}. We used magnitudes rather than colours to derive the previous result, but it can be reframed as a colour measurement: Fig. 9 shows that the largest mismatch in the value of is obtained by comparing the U and I bandpasses. We take advantage of this observation in Fig. 10 by imposing the condition that the derivative of the U − I colour with respect to the extinction colour X_{3} be zero after application of the extinction correction. The value of R_{V} obtained in Fig. 10 is now R_{V} = 2.265, which is almost the same as in Eq. (35), but the error is slightly smaller, as expected; the contributions from flux calibration error (0.03 mag), the redshift error, and the ‘grey’ fluctuation are suppressed. With this determination, the corresponding scale factor for the ‘extinction scale’ s = 1.000. If we turn to the usual ‘operational’ definition of R_{V}, as A_{V}/E(B − V), it is seen in Table 4 that R_{V} = A_{V}/X_{4} = R_{V4} = 2.7. However, this numerical value is specific to SNe and our tophat filters.
The corrected magnitudes in the UBVRI bandpasses are shown in Fig. 12 for the value R_{V} = 2.25. As expected, no residual dependence of the (corrected) magnitude on the extinction is observed in any filter. The RMS (which was not minimised), has a value of 0.13 mag in all filters, which is smaller than the scatter of 0.15 mag typically seen in the SALT2 analyses of our sample (Saunders et al. 2018), but larger than the result from refined analyses designed to minimise this fluctuation, as in Boone et al. (2021). The blue straight line in Fig. 12 is a weighted fit and the red line is unweighted; the two slopes are almost identical. The magnitude fluctuation is ‘grey’, that is, it is almost identical over all bandpasses to a remarkable accuracy. For comparison, we show the magnitude correlations obtained for R_{V} = 3.1 in Fig. 11. The scale factor must now be set to s = 0.857 to ensure that the derivative . The error on the derivative of d(m_{U}−m_{I})/dX_{3} in Fig. 10 is 0.04, as evaluated from the spread of the simulation samples, meaning that the significance of the rejection of R_{V} = 3.1 (and higher) is actually 3.5 standard deviations.
Fig. 12. UBVRI magnitudes corrected for I_{3}, X_{3} with R_{V} = 2.25, s = 1.0024. As for the previous figure, the red straight line is weighted by errors. 
7.2. Deriving R_{V} from the grey fluctuation
We mention in the previous section that the fluctuation of the magnitudes around the common ‘grey’ fluctuation of each SN Ia i, , was remarkably small. This observation can be turned around and used to derive R_{V}: in each bandpass F, we compute the RMS σ_{F} of the dispersion of the corrected magnitudes around m_{grey}. We show how this dispersion varies with R_{V} in Fig. 13: a minimum is reached for R_{V} = 2.277, which confirms the previous result. The dispersion in each UBVRI bandpass with respect to the average offset over all filters is remarkably small, that is 0.0114, 0.0197, 0.0145, 0.0148, and 0.0060, respectively, which is everywhere smaller than 0.02. This result quantifies the level at which the remaining magnitude fluctuation is grey, and also confirms the validity of the modelling implemented (up to this ‘unaccounted’ grey fluctuation).
Fig. 13. RMS of the offsets (with respect to m_{grey}) of the corrected UBVRI magnitudes as a function of R_{V}. 
8. Simulation
The previous results can be summarised by ten numbers: the five extinction corrections to the extinction formula in the five bandpasses used for our analysis, the three intrinsic colour ratios, the extinction scale s, and the extinction parameter R_{V}. These are the ‘effective’ values found within our algorithms, but not necessarily the ‘true’ values. The biases and the errors are evaluated below based on a simulation that mimics the observations. We use the values of the parameters tuned to describe the data to generate the events, and then proceed with the same analysis. The differences from the reconstructed values allow us to find the methodological biases, and the scatter among the different generated samples helps us find the statistical error.
8.1. Generation of colours with their error
We generate nine samples of 165 SNe Ia UBVRI magnitudes. The generated intrinsic and extinction components of colour 1 ≡ U − V follow the observed distribution in the data after smoothing by a sliding local averaging (to avoid counting the statistical fluctuations twice): in each histogram bin, the data distribution of I_{1} and X_{1} was replaced by the average of three adjacent bins. Whenever the bin statistics was below four events, the average of five bins was used. The other intrinsic colours obey , and the extinction colours follow from . The coefficients δ_{ij} used in the generation are the ones found in the data analysis; that is, they include the correction to the FM99 formula described in Table 7. The simulation shows that the analysis does recover the input corrections, with a small bias. As the measurement errors are included in the observed data distribution, they should have been removed from the generation in the simulation. We instead introduced a single simulation scaling factor of the distribution, tuned so that the reconstructed intrinsic components should match the observations once the analysis is applied; it depends on the simulated sample, but its difference from unity never exceeds 0.01. This effect is expected to be smaller for the extinction colours X_{i}, which have a larger range, and the observed distribution of the data has been used.
8.2. Colour noise
The noise is derived from the residuals observed in the data, and it combines measurement and modelling errors. The residuals ϵ_{i} measured from the data in four colours (i.e. 2 ≡ U − R, 5 ≡ B − R, 7 ≡ V − R, 9 ≡ R − I) are used to derive a 4 × 4 covariance matrix ⟨ϵ_{i}ϵ_{j}⟩. The square roots of the dominant eigenvalues are σ_{1} = 0.026525, σ_{2} = 0.017947; the other eigenvalues are negligible (∼10^{−4} and ∼3 10^{−5}). The Gaussian noise of the two eigenvectors is then projected onto the four bandpasses to obtain a simulated ‘colour noise’ for the four colours. We find that this noise has to be multiplied by 1.004 in the simulation to reproduce the observed residuals. This effect may arise from the contribution of modelling errors to the data, while the model is exact in the simulation.
8.3. Generation of magnitudes
The R magnitude is used as the reference. A correction reproducing its (small) correlation with the intrinsic colour is added. The magnitudes in the UBVI bandpasses are obtained by adding the U − R, B − R, V − R, and I − R colours, as found from their intrinsic and extinction components, and including the ‘colour noise’. Finally, a Gaussian grey magnitude fluctuation with σ(grey) = 0.12 is added, though its value is irrelevant for all the results of the present study. The generated magnitudes of event n in bandpass R and U are therefore described as
where u_{1}(n) and u_{2}(n) are centred random Gaussian variables with a standard deviation of unity, n_{21} and n_{22} are the components of the corresponding eigenvectors on the U − R colour, and z_{R1} is the observed correlation between the R magnitude and the intrinsic colour 1. The generated magnitudes are processed by the same algorithms as the data (after ew^{Ca} and ew^{Si} corrections), providing reconstructed values for I_{i} and X_{i}, reconstructed extinction corrections, reconstructed intrinsic couplings, and a reconstructed R_{V}.
8.4. Simulated and observed distributions
As we want to evaluate the errors and biases of the data analysis, with the help of the simulation, both the distributions and the errors must be similar in the observations and the simulation. The (I_{i}, X_{i}) obtained through the colour analysis (as in the data), and the generated values are compared for three colours in Fig. 14 for the simulated sample 1 and in Table 9. There is a very small scatter between the generated and reconstructed colour components, and their ratio is close to unity. This scatter varies from 0.002 to 0.006 for I_{i}, and from 0.001 to 0.006 for X_{i}. The systematic trend in the ratio of measured and generated colour is partially statistical and fluctuates from one sample to another in a range of 0.02, which justifies using the observed distributions as an input for the simulation. The scatter in Fig. 14 is surprisingly smaller than the RMS of the observed residuals (in average 0.0176 mag). This is a consequence of the averaging over different colours performed in Eq. (31). The simulated accuracy is tuned to reproduce the observed residuals by the introduction of two ad hoc factors in the noise and the rescaling of the distribution of I_{1} (of colour U − V) used in the generation. The distribution of the residuals is shown in Fig. 15 and Table 9.
Fig. 14. Comparison of intrinsic and extinction components (mag) in generation and reconstruction for colours c_{0}, c_{5}, and c_{9}. The small differences found justify using the (smoothed) observed distributions of I_{1} and X_{1} as input. 
For the same three colours (c_{0}, c_{5}, c_{9}), the intrinsic and extinction colour distributions for the data and one of the simulated samples are shown in Fig. 16. The results for all colours are summarised in Table 9. The residuals in the data and the simulation are compared in Fig. 15 and Table 9. The intrinsic and extinction colours in the simulation and observation agree, and the trends of the errors are also well reproduced by the simulation.
Fig. 15. Residuals in observations and simulations for colours c_{0}, c_{2}, c_{5}, c_{6}, c_{8}, and c_{9} (R_{V} = 2.25). The RMS of the residuals of all colours are given in Table 9. The blue shade is the simulated distribution, the orange shade the data. 
Fig. 16. Distributions of intrinsic and extinction components of c_{0}, c_{6}, and c_{9} in observation and simulation for R_{V} = 2.25 (mag). The extinction is asymmetric, with an extended tail to high values. The intrinsic part is symmetric. The RMS for all colours is given in Table 9. The blue shade is the simulated distribution, the orange shade the data. 
8.5. Errors on the extinction formula corrections
We rely on the nine simulated samples to control the bias in the measurements arising from the implemented algorithms, as well as the errors on the extracted parameters: extinction corrections Dex_{F}, the value of R_{V}, and the extinction scale s_{X}. The analysis was performed on each of the samples, on the same grid of values of R_{V}. The mean value of the reconstructed extinction corrections and their dispersion is given in Table 10.
Extinction correction for R_{V} = 2.25.
When comparing to Table 4, it is seen that the corrections to the extinction are statistically significant, rising from 1 to 10 percent over the visible spectrum once the bias is taken into account. The algorithm used in deriving the corrections is able to extract their value, with some biases. Accepting the biases found in the simulation, our final values for the corrections to the extinction coefficients are given in the last column of Table 10 for R_{V} = 2.25. The error on the bias correction is σ(Dex_{F})/3 as nine samples are generated.
8.6. Errors on the intrinsic coupling coefficients γ
For each simulation sample and each value of R_{V} in the grid, the three intrinsic coupling coefficients γ_{10}, γ_{20}, and γ_{30}, are evaluated as in the data analysis, that is by minimising the scatter in the plot comparing the input values to the results. The standard deviation of the three coefficients from one sample to the next is respectively 0.0075, 0.0061, and 0.0041. The values averaged over all samples allow us to estimate the bias arising from the analysis; it is found to be of the same order as the statistical error. The largest contribution to the error is by far the one discussed in Sect. 5.1, which arises from different choices of auxiliary colours, and is given in Table 6. The statistical error will be added into the final uncertainty.
8.7. Errors on the extinction parameter R_{V}
The measurement of R_{V} is also performed on the simulation. The value in the generation was R_{V} = 2.20, and the average value found is R_{V} = 2.284 ± 0.112. We correct the result obtained in the data analysis, of 2.265, to R_{V} = 2.265 − 0.084 = 2.181 ± 0.117 (the error is increased to account for the uncertainty on the bias correction). As the values found in the four filters BVRI are almost fully correlated, we use the error on R_{V} derived from the U − I colour. We add the impact of an estimated systematic error in the measurement of colours and other corrections (earth atmosphere, galactic extinction, host galaxy subtraction, silicon and calcium contribution). As we only use differences between mean and measured colours, many instrumental systematic errors cancel out. To evaluate the impact on R_{V}, all colours were modified by a shift proportional to their value and reaching 0.005 mag over their range of variation. The outcome is a change in R_{V} by 0.05. The final result for the parameter R_{V} (including the error on the bias) is then
9. No calcium or silicon line corrections
We tried to reproduce our analysis without performing the initial correction for the ew^{Si} and ew^{Ca} spectral lines. The intrinsic colours are expected to be much larger than in the previous analysis. They reach 0.10 in colours (U − B, U − V, U − R). The mean colour residual is then 0.0402 (instead of 0.0176 with the colour corrections from ew^{Ca} and ew^{Si}). The square roots of the dominant eigenvalues of the covariance matrix from the residuals are 0.0724 and 0.0296, meaning that it would be necessary to add two spectral coordinates for each SN in order to reduce these residuals, instead of tolerating these contributions as ‘noise’ as is done in the present work. The two eigenvectors are (0.5842, 0.5974, 0.4373, −0.3324) and (0.5972, −0.7405, 0.2908, 0.1013). In the absence of ew^{Si} and ew^{Ca} corrections, the model we used fails to describe the observations, as we are missing key spectral information. Instead, the approach we take is to add back the colour contributions of Si IIλ4131 and Ca II H&Kλ3945 –which are subtracted in Sect. 4.1– to obtain a ‘full’ intrinsic component of colour variability; this is shown in Fig. 17.
Fig. 17. Distributions of six ‘full’ intrinsic colours including Si II and Ca II H&K contributions. 
The range of the full intrinsic colour component is given in Table 11 for all colours. The small size of the intrinsic part of e_{4} ≡ B − V justifies using it as a simplified indicator of extinction for SNe Ia.
Range of the full intrinsic colour component and correlations between the full calcium–silicon colour correction DCaSi and the SALT x_{1} variable (R_{V} = 2.25).
The standard choice of colour 4 ≡ B − V as an extinction indicator is fortunate, as the intrinsic content happens to be small. However, a large part of the calcium and silicon correction can be recovered, with information extracted from the light curve. The SALT2 stretch variable x_{1}, as described in Betoule et al. (2014), is strongly correlated to the combined calcium–silicon corrections DCaSi(i) = σ_{i}ew^{Si} + κ_{i}ew^{Ca}, as shown in Fig. 18 and in Table 11 for all colours. The dispersion of the full intrinsic colour can be compared to the results of previous determinations: (Jha et al. 2007), ; , (Nobili & Goobar 2008), (Mandel et al. 2017), , (2populations/decline rate, Wojtak et al. 2023). The values of found in this work for the two populations of the H_{α} subsample (galaxies with high and low H_{α} signal, as discussed in item 10 of Sect. 10) are: , .
Fig. 18. Correlations between the combined calcium–silicon correction DCaSi and the SALT2 variable x_{1} for four colours. 
The smallest scatter of DCaSi is found for colour 4 ≡ B − V, with σ = 0.00349. This provides an estimate of the contribution of measurement errors to colours, all colours being processed in the same way. Extra physical variabilities worsen the evaluation of the DCaSi contribution of other colours.
10. Conclusions
The focus of this work is the extraction for each SN of our sample of extinction and intrinsic colour components, using the extinction formula for dust in our galaxy as leverage. This colour analysis is well adapted to the determination of the reddening, but does not make use of the full power of the spectral information for the intrinsic part. The linearity that we find in the correlations between different colours is nevertheless an interesting feature by itself.

The model presented in this work assumes the factorisation of the extinction by the host galaxy, without any correlation with the intrinsic properties of the SN. We acknowledge that the presence of circumstellar dust might invalidate this assumption.

The modelling of colour fluctuations requires the inclusion of four corrections to the extinction formula in the BVRI bandpasses, three intrinsic coefficients correlating the different colours, and one value for R_{V}. An additional scale parameter for the extinction should be set to unity. The colour accuracy achieved varies from 0.011 to 0.028 mag depending on the colour. All colours discussed in this work are ‘restframe’ colours. Such corrections could arise from different average dust contents in host galaxies and ours. One of the assumptions of the model (factorisation of the extinction, linearity of the correlations of intrinsic colours) could be wrong, though plausible, but we are surprised by the large difference found in the corrections to the R and I bandpasses Dex_{R} and Dex_{I} shown in Table 4 and observed directly in Fig. 5. The simulation confirms that the algorithm we use is indeed sensitive to these small corrections. The errors arising from the linear approximation cannot account for the size of the effects found, and additional intrinsic contributions are extremely unlikely to induce a correlation with extinction as strong as the one observed in Fig. 5.

A single extraintrinsic contribution beyond Si IIλ4131 Å and Ca II H&Kλ3945 Å is needed to account for the colour correlations within the accuracy described immediately above We can then extract intrinsic and extinction colour components for each SN.

We show that the remaining magnitude fluctuation of 0.13 mag is independent of the five UBVRI bandpasses to an accuracy of better than 0.02 mag, reaching 0.007 mag for the I bandpass. There is no attempt to minimise this magnitude dispersion directly in our procedure.

When the extinction formula of FM99 is used as leverage to extract this extra intrinsic (onedimensional) component, and the associated extinction component, the intrinsic coupling coefficients between the different colours are consistent over the full sample of SNe Ia at the level of 0.01–0.05. The values of the three independent coefficients are (for R_{V} = 2.25)
The errors quoted include the systematic errors from the colour choice, and the statistical error derived from the simulation. The extra ‘intrinsic’ content of U − I is half the size of U − B or U − V, with a large uncertainty.

After correcting for the analysis bias found in the simulation, the optimal mean value of the FM99 parameter R_{V} over our sample is found to be
This value is obtained within our specific treatment of the extinction scale, which is not universally adopted. The usual ‘operational’ definition R_{V} = A_{V}/E(B − V) = R_{V4} then leads to R_{V4} = 2.7 according to Table 4, which includes the effect of our different filter definition.

The associated scaling factor for the extinction component is measured to be 1.001 ± 0.044, and, as desired for consistency, is compatible with unity.

A previous SNfactory analysis in C11 found R_{V} = 2.8 ± 0.30, with a different method and a substantially different data set. A few improvements are brought here, namely the rescaling of the extinction coefficients, ensuring consistency with the B and V SN spectra; the presence of an ‘extraintrinsic’ colour component; and we avoid an arbitrary extinction parameter for each SN (A_{V}). The presence of the grey fluctuation was ignored by C11, and the extraction of the extinction colour component was not explicit. The R_{V} value that we obtain is dependent on the method used, with an emphasis on colour correlations rather than a minimisation of Hubble residuals. The smaller number of parameters involved and the accuracy of the colour description is an argument in favour of the present method.

Within the accuracy of the present observations, we see no hint of a variation of the parameter R_{V} from one SN to another. As we find corrections to the extinction formula to be necessary, for all values of R_{V} in the range 1.95–2.60, the meaning of this parameter is slightly blurred. Simulations have been performed with samples including a Gaussian distribution of the value of R_{V} around the value of 2.20 and with an RMS of 0.5. The results do not show any significant change, and we would thus be insensitive to such a fluctuation.

Many instances have been given in the past of a correlation between the star formation rate (Hα signal, Rigault et al. 2013) or the host galaxy mass (Hamuy et al. 2000; Kelly et al. 2010) and the absolute magnitudes of SNe, as well as their intrinsic properties, suggesting two populations depending on their age. We checked the impact on our analysis of including a sample of 85 SNe from the SNfactory collaboration common to the Rigault et al. (2013) sample. A striking difference is seen in the amount of extinction of the two groups. The young SNe (high Hα signal of the host galaxy) have significantly more dust, with ⟨X_{3}⟩ = 0.0124 ± 0.036, whilst the older ones (low galactic Hα signal) have less extinction and ⟨X_{3}⟩= − 0.0577 ± 0.025. Our corrections slightly reduce the discrepancy in absolute magnitude of the two samples from 0.095 to 0.087 mag. The statistics of the two groups (47 and 37 SNe) are too small to evaluate a difference in the associated R_{V}.
The colour model used here achieves an impressive level of accuracy in the description of the ’grey’ magnitude fluctuation, which may help interpret the effect. The accuracy reached on the (smaller) intrinsic colour might be improved by using additional spectral information, in particular Si IIλ6355 Å, which is better measured. This will require dealing with the additional complexity arising from nonlinear effects.
To recover the full intrinsic colours, the suppressed intrinsic colour fluctuations associated to Si and Ca contributions should be added back, as in Sect. 9.
The simulation in Sect. 8, which ensures the consistency of the analysis, starts from the Si and Cacorrected colours.
Acknowledgments
We thank the Nearby Supernova Factory collaboration for providing access to its data. The contribution of N. Chotard was crucial in the processing of the data, the estimate of the magnitudes, and the evaluation of the equivalent widths used in this analysis. We are grateful to the technical and scientific staff of the University of Hawaii 2.2 m telescope, to the Palomar Observatory, and to the High Performance Research and Educational Network (HPWREN) for their assistance in obtaining these data. We also thank the people of Hawaii for the access to Mauna Kea. We thank D. Birchall for observing assistance. This work was supported by the Director, Office of Science, Office of High Energy Physics of the US Department of Energy under contract DEAC0205CH11231. This work was supported in France by CNRS/IN2P3, CNRS/INSU, and PNC. Support in Germany was provided by the DFG through TRR33 “The Dark Universe”, and in China from Tsinghua University grant 985 and NSFC grant 11173017. Some results were obtained using resources and support from the National Energy Research Scientific Computing Center, supported by the Director, Office of Science, Office of Advanced Scientific Computing Research, of the US Department of Energy under contract DEAC0205CH11231. HPWREN is funded by National Science Foundation Grant ANI0087344, and the University of California, San Diego. We thank G. Aldering and A. Kim for their remarks, which have helped clarify the present text. We thank particularly the referee for his careful reading and numerous useful questions and comments. The authors take responsibility for the remaining inadequacies.
References
 Aldering, G., Antilogus, P., Bailey, S., et al. 2006, ApJ, 650, 510 [NASA ADS] [CrossRef] [Google Scholar]
 Aldering, G., Antilogus, P., Aragon, C., et al. 2020, Res. Notes Am. Astron. Soc., 4, 63 [Google Scholar]
 Amanullah, R., Goobar, A., Johansson, J., et al. 2014, ApJ, 788, L21 [NASA ADS] [CrossRef] [Google Scholar]
 Amanullah, R., Johansson, J., Goobar, A., et al. 2015, MNRAS, 453, 3300 [Google Scholar]
 Bessell, M., & Murphy, S. 2012, PASP, 124, 140 [NASA ADS] [CrossRef] [Google Scholar]
 Betoule, M., Kessler, R., Guy, J., et al. 2014, A&A, 568, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Boone, K., Aldering, G., Antilogus, P., et al. 2021, ApJ, 912, 70 [NASA ADS] [CrossRef] [Google Scholar]
 Borkowski, K. J., Blondin, J. M., & Reynolds, S. P. 2009, ApJ, 699, L64 [NASA ADS] [CrossRef] [Google Scholar]
 Bronder, T. J., Hook, I. M., Astier, P., et al. 2008, A&A, 477, 717 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Brout, D., & Scolnic, D. 2021, ApJ, 909, 26 [NASA ADS] [CrossRef] [Google Scholar]
 Burns, C. R., Stritzinger, M., Phillips, M. M., et al. 2014, ApJ, 789, 32 [Google Scholar]
 Cardelli, J. A., Clayton, G. C., & Mathis, J. S. 1989, ApJ, 345, 245 [Google Scholar]
 Chotard, N. 2011, PhD thesis, Université Claude Bernard– Lyon I, Institut de Physique Nucléaire de Lyon, France [Google Scholar]
 Chotard, N., Gangler, E., Aldering, G., et al. 2011, A&A, 529, L4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Ferretti, R., Amanullah, R., Bulla, M., et al. 2017, ApJ, 851, L43 [NASA ADS] [CrossRef] [Google Scholar]
 Fitzpatrick, E. L. 1999, PASP, 111, 63 [Google Scholar]
 Fitzpatrick, E. L., & Massa, D. 2005, AJ, 130, 1127 [NASA ADS] [CrossRef] [Google Scholar]
 Guy, J., Astier, P., Baumont, S., et al. 2007, A&A, 466, 11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hamuy, M., Phillips, M. M., Suntzeff, N. B., et al. 1996, AJ, 112, 2391 [Google Scholar]
 Hamuy, M., Trager, S. C., Pinto, P. A., et al. 2000, AJ, 120, 1479 [Google Scholar]
 Huang, X., Raha, Z., Aldering, G., et al. 2017, ApJ, 836, 157 [NASA ADS] [CrossRef] [Google Scholar]
 Jha, S., Riess, A. G., & Kirshner, R. P. 2007, ApJ, 659, 122 [NASA ADS] [CrossRef] [Google Scholar]
 Kelly, P. L., Hicken, M., Burke, D. L., Mandel, K. S., & Kirshner, R. P. 2010, ApJ, 715, 743 [Google Scholar]
 Lantz, B., Aldering, G., Antilogus, P., et al. 2004, in Optical Design and Engineering, Proc. SPIE, 5249, 146 [NASA ADS] [CrossRef] [Google Scholar]
 Léget, P.F., Gangler, E., Mondon, F., et al. 2020, A&A, 636, A46 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Lira, P., Suntzeff, N. B., Phillips, M. M., et al. 1998, AJ, 115, 234 [NASA ADS] [CrossRef] [Google Scholar]
 Mandel, K. S., Scolnic, D. M., Shariff, H., Foley, R. J., & Kirshner, R. P. 2017, ApJ, 842, 93 [NASA ADS] [CrossRef] [Google Scholar]
 Nagao, T., Maeda, K., & Yamanaka, M. 2017, ApJ, 835, 143 [NASA ADS] [CrossRef] [Google Scholar]
 Nobili, S., & Goobar, A. 2008, A&A, 487, 19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Nugent, P., Phillips, M., Baron, E., Branch, D., & Hauschildt, P. 1995, ApJ, 455, L147 [NASA ADS] [CrossRef] [Google Scholar]
 O’Donnell, J. E. 1994, ApJ, 422, 158 [Google Scholar]
 Perlmutter, S., Aldering, G., Goldhaber, G., et al. 1999, ApJ, 517, 565 [Google Scholar]
 Phillips, M. M., Lira, P., Suntzeff, N. B., et al. 1999, AJ, 118, 1766 [Google Scholar]
 Rieke, G. H., & Lebofsky, M. J. 1985, ApJ, 288, 618 [Google Scholar]
 Riess, A. G., Filippenko, A. V., Challis, P., et al. 1998, AJ, 116, 1009 [Google Scholar]
 Rigault, M., Copin, Y., Aldering, G., et al. 2013, A&A, 560, 66 [Google Scholar]
 Saunders, C., Aldering, G., Antilogus, P., et al. 2018, ApJ, 869, 167 [NASA ADS] [CrossRef] [Google Scholar]
 Scalzo, R. A., Aldering, G., Antilogus, P., et al. 2010, ApJ, 713, 1073 [NASA ADS] [CrossRef] [Google Scholar]
 Schlafly, E. F., Meisner, A. M., Stutz, A. M., et al. 2016, ApJ, 821, 78 [NASA ADS] [CrossRef] [Google Scholar]
 Schlafly, E. F., Peek, J. E. G., Finkbeiner, D. P., & Green, G. M. 2017, ApJ, 838, 36 [NASA ADS] [CrossRef] [Google Scholar]
 Thorp, S., Mandel, K. S., Jones, D. O., Ward, S. M., & Narayan, G. 2021, MNRAS, 508, 4310 [NASA ADS] [CrossRef] [Google Scholar]
 Wang, L., Goldhaber, G., Aldering, G., & Perlmutter, S. 2003, ApJ, 590, 944 [NASA ADS] [CrossRef] [Google Scholar]
 Wojtak, R., Hjorth, J., & Osman Hjortlund, J. 2023, MNRAS, 525, 5187 [NASA ADS] [CrossRef] [Google Scholar]
Appendix A: Linearisation of the extinction formula
The observed photoncount distribution was given in Eq. 11 where the function ϕ(λ) is provided by Cardelli et al. (1989), O’Donnell (1994), Fitzpatrick (1999), Fitzpatrick & Massa (2005). To simplify the computations, we assume the extinction E(B − V) is small enough to allow a linear approximation under the exponentiation.
We define A(λ, R_{V}) = ϕ(λ, R_{V}) E(B − V), the observed spectrum (photon counts) is then
When the unextincted SN spectrum s_{0}(λ) is integrated over a spectral band, such as B, with transmission T(λ), the rest frame flux is:
In our synthetic photometry, unittransmission tophat filters are used, and so we omit T(λ). Using the previous linear approximation:
The logarithm can be expanded to give
In this approximation, the extinction is simply weighted by the spectrum over the filter bandpass. The maximal value of the extinction component X_{4} is 0.4, and the error arising from the linear approximation in the evaluation of the extinctions coefficient R_{F4} in the different (tophat) filters is given in Table A.1. The magnitude correction is obtained by multiplying the coefficients R_{F4} by X_{4}, and so the error on this correction never exceeds 0.005 mag, and cannot account for the effects observed in Fig. 5.
Error on the bandpass extinction coefficients arising from the linear approximation for X_{4} = 0.4.
All Tables
Extinction coefficients from FM99 and UBVRI corrections to the extinction formula for the case R_{V} = 2.25.
Residual (mag) dependence on the extinction colour R_{V} = 2.25 (after correction of the extinction formula).
Range of the intrinsic and extinction colours for observation, simulation (seed2), and residuals (mag).
Range of the full intrinsic colour component and correlations between the full calcium–silicon colour correction DCaSi and the SALT x_{1} variable (R_{V} = 2.25).
Error on the bandpass extinction coefficients arising from the linear approximation for X_{4} = 0.4.
All Figures
Fig. 1. Correlations between the first five colours , , , , , (mag) and ew^{Si} (Å). The slopes of the straight line fits (red lines) give the value of K^{Si} directly. The correlations of the colours can be derived from the first four. 

In the text 
Fig. 2. Correlations between the same colours as Fig. 1, corrected for the ew^{Si} correlation (mag), as a function of ew^{Ca} (Å). The slopes of the straight line fits (red lines) give the value of K′^{Ca}. 

In the text 
Fig. 3. Intrinsic and extinction ratios as defined in Eq. (31) (i = 2, 6, 9) compared to the input and calculated δ_{0i} for R_{V} = 2.25. The compatibility is at the level of 0.001 for the extinction, and the relative accuracy of 0.05 for the (smaller) intrinsic colour (worst case). The red lines are linear fits to the data, and γ^{out}, δ^{out} are the values of the linear coefficients. 

In the text 
Fig. 4. Relation between input and output (⟨I_{i}/I_{j}⟩) values for the coefficients γ_{ij}. The overall agreement for all γ_{ij} (excluding those involving colours 5 ≡ (B − R) and 7 ≡ (V − R)) is within 0.01. 

In the text 
Fig. 5. Correlation of the colour residuals (mag) for R_{V} = 2.25 with the extinction colour for colour differences e_{0}, e_{2}, e_{5}, e_{6}, e_{8}, e_{9} (no extinction formula correction). The strong correlations observed would forbid a safe evaluation of R_{V}. 

In the text 
Fig. 6. Correlation of residuals of colour differences (mag) (e_{0}, e_{2}, e_{5}, e_{6}, e_{8}, e_{9}) with extinction after correcting the extinction formula with R_{V} = 2.25. There are only four independent colours, meaning that the four extinction corrections allow us to cancel the correlations of all ten colours. 

In the text 
Fig. 7. Residuals ϵ_{i} (mag) after the reconstruction of six colours (c_{0}, c_{2}, c_{5}, c_{6}, c_{8}, c_{9}). The full list of residuals is given in Table 8. 

In the text 
Fig. 8. Extinction scale factor s as a function of R_{V} for dm_{V}/dX_{3} = 0. For each value of R_{V}, there is a value of the extinction scale factor s that cancels the dependence of the V magnitude on extinction. 

In the text 
Fig. 9. for the UBVRI bandpasses as a function of R_{V}. Although the dependence of the V bandpasses as a function of extinction has been cancelled for each value of R_{V}, the extinction in other bandpasses is not adequately corrected. 

In the text 
Fig. 10. d(U−I)^{corr}/dX^{3} as a function of R_{V}. The difference between U and I magnitudes (colour 3) is the most sensitive to the choice of R_{V}. The error bar is the dispersion found in the simulation. 

In the text 
Fig. 11. UBVRI magnitudes corrected for I_{3}, X_{3} with R_{V} = 3.100, s = 0.73344. The derivative is compatible with zero in each bandpass, but the variation from U to I is quite significant, taking correlations into account. The red line is weighted by errors and the blue straight line is unweighted. 

In the text 
Fig. 12. UBVRI magnitudes corrected for I_{3}, X_{3} with R_{V} = 2.25, s = 1.0024. As for the previous figure, the red straight line is weighted by errors. 

In the text 
Fig. 13. RMS of the offsets (with respect to m_{grey}) of the corrected UBVRI magnitudes as a function of R_{V}. 

In the text 
Fig. 14. Comparison of intrinsic and extinction components (mag) in generation and reconstruction for colours c_{0}, c_{5}, and c_{9}. The small differences found justify using the (smoothed) observed distributions of I_{1} and X_{1} as input. 

In the text 
Fig. 15. Residuals in observations and simulations for colours c_{0}, c_{2}, c_{5}, c_{6}, c_{8}, and c_{9} (R_{V} = 2.25). The RMS of the residuals of all colours are given in Table 9. The blue shade is the simulated distribution, the orange shade the data. 

In the text 
Fig. 16. Distributions of intrinsic and extinction components of c_{0}, c_{6}, and c_{9} in observation and simulation for R_{V} = 2.25 (mag). The extinction is asymmetric, with an extended tail to high values. The intrinsic part is symmetric. The RMS for all colours is given in Table 9. The blue shade is the simulated distribution, the orange shade the data. 

In the text 
Fig. 17. Distributions of six ‘full’ intrinsic colours including Si II and Ca II H&K contributions. 

In the text 
Fig. 18. Correlations between the combined calcium–silicon correction DCaSi and the SALT2 variable x_{1} for four colours. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.