The Delay Times of Type Ia Supernova

TL;DR

This paper measures the delay time distribution of Type Ia supernovae in field galaxies using SDSS data, finding a power-law decline with time and comparing methods to ensure consistent results.

Contribution

It introduces an improved technique for measuring the delay time distribution, including normalization, and compares it with existing methods for the first time.

Findings

Power-law index s = -1.34 with uncertainties

Normalization of supernova rate A determined

Results consistent with cluster galaxy distributions

Abstract

The delay time distribution of Type Ia supernovae (the time-dependent rate of supernovae resulting from a burst of star formation) has been measured using different techniques and in different environments. Here, we study in detail the distribution for field galaxies, using the SDSS DR7 Stripe 82 supernova sample. We improve a technique we introduced earlier, which is based on galaxy color and luminosity, and is insensitive to details of the star formation history, to include the normalization. Assuming a power-law dependence of the supernova rate with time, ${\rm DTD}(t)=A(t/{\rm 1\,Gyr})^{s}$, we find a power-law index $s= -1.34 ^{+0.19} _{-0.17}$ and a normalization $\rm{log}\ A= -12.15 ^{+0.10} _{-0.13}\, {\rm dex}(M_\odot ^{-1}\, \rm{yr}^{-1})$, corresponding to a number of type Ia supernovae integrated over a Hubble time of $PE = 0.004^{+0.002} _{-0.001}\, M_\odot ^{-1}$. We also implement a method used by Maoz and collaborators, which is based on star formation history reconstruction, and find that this gives a consistent result for the slope, but a lower, marginally inconsistent normalization. With our normalization, the distribution for field galaxies is made consistent with that derived for cluster galaxies. Comparing the inferred distribution with predictions from different evolutionary scenarios for type Ia supernovae, we find that our results are intermediate between the various predictions and do not yet constraint the evolutionary path leading to SNe~Ia.

Equations

$$\textrm{DTD}(t)\equiv\left\{\begin{array}[]{ll}0&\textrm{if}\,\,\,t<t_{\rm{WD}}\\ A\,\left(\frac{t}{1\,\rm{Gyr}}\right)^{s_{1}}&\textrm{if}\,\,\,t_{\rm{WD}}\leq t\leq t_{\rm{c}}\\ B\,\left(\frac{t}{1\,\rm{Gyr}}\right)^{s_{2}}&\textrm{if}\,\,\,t\geq t_{\rm{c}}.\end{array}\right.$$

$$r(t)\equiv \text{DTD}(t)\ast \text{SFR}(t)={\int }_{t_{\mathrm{WD}}}^t\text{DTD}(t^{\prime })\text{SFR}(t-t^{\prime })d{t}^{\prime },$$

$$r_i\approx m_{i,1}{\psi }_1+m_{i,2}{\psi }_2+m_{i,3}{\psi }_3,$$

$$⟨{\psi }_j⟩\approx ⟨{\mathrm{DTD}}_{\mathrm{j}}⟩=\frac{1}{{\mathrm{t}}_{\mathrm{f},\mathrm{j}}-{\mathrm{t}}_{\mathrm{i},\mathrm{j}}}{\int }_{{\mathrm{t}}_{\mathrm{i},\mathrm{j}}}^{{\mathrm{t}}_{\mathrm{f},\mathrm{j}}}\mathrm{DTD}(t|A,{\mathrm{s}}_1,{\mathrm{s}}_2)\mathrm{dt},$$

$$r_i\approx \sum _{j=1,2,3}{m}_{i,j}\times {\psi }_j(\overline{t_j}|A,s_1,s_2).$$

$$PE\equiv {\int }_{t_{\mathrm{WD}}}^{13.7\mathrm{Gyr}}\mathrm{DTD}(\mathrm{t})dt.$$

$${\bar{t}}_1=\frac{t_{\mathrm{WD}}(\mathrm{Gyr})+0.42}{2}\mathrm{Gyr},{\bar{t}}_2=\frac{0.42+2.4}{2}=1.41\mathrm{Gyr}\mathrm{and}{\bar{t}}_3=\frac{2.4+14}{2}=8.2\mathrm{Gyr}.$$

$$r_{i,1}$$

$$r_{i,2}$$

$$r_{i,3}$$

$$r_i(A,s_1,s_2)=$$

$$s_i(A,s_1,s_2)=$$

$${\lambda }_i(A,s_1,s_2)=$$

$$t_i=\frac{269}{365.25(1+z_i)}[\mathrm{yr}],$$

$$\epsilon_{i}=\left\{\begin{array}[]{ll}0.72&\textrm{if}\,\,\,0\leq z_{i}\leq 0.175\\ -3.2\,z_{i}+1.28&\textrm{if}\,\,\,0.175\leq z_{i}\leq 0.4.\end{array}\right.$$

$$P_i(n_i|{\lambda }_i(A,s_1,s_2))\approx \frac{{\lambda }_i^{n_i}{e}^{-{\lambda }_i}}{n_i!}.$$

$$L(A,s_1,s_2)=\prod _i{P}_i(n_i|{\lambda }_i(A,s_1,s_2)),$$

$$\ln L(A,s_1,s_2)=\sum _i\ln \left(\frac{{\lambda }_i^{n_i}{e}^{-{\lambda }_i}}{n_i!}\right).$$

$$\ln L(A,s_1,s_2)=-\sum _{i\ne k}{\lambda }_i+\sum _k\ln \left({\lambda }_k{e}^{-{\lambda }_k}\right)=-\sum _{i\ne k}{\lambda }_i+\sum _k\ln \left({\lambda }_k\right)-\sum _k{\lambda }_k=-\sum _i{\lambda }_i+\sum _k\ln \left({\lambda }_k\right).$$

$$\ln L(A,s_1,s_2)=-N_{\exp }+\sum _k\ln \left({\lambda }_k\right).$$

$$\ln L(A_0,s_1,s_2)=-N_{\mathrm{obs}}+\sum _k\ln \left({\lambda }_k(A_0,s_1,s_2)\right).$$

$$\begin{array}{rl}\ln L(A,s_1,s_2)& =-\sum _{i}f{\lambda }_i(A_0,s_1,s_2)+\sum _k\ln [f{\lambda }_k(A_0,s_1,s_2)]\\ & =-f\sum _i{\lambda }_i(A_0,s_1,s_2)+\sum _k\ln {\lambda }_k(A_0,s_1,s_2)+\sum _k\ln (f).\\ & =-f{N}_{\mathrm{obs}}+[N_{\mathrm{obs}}+\ln L(A_0,s_1,s_2)]+N_{\mathrm{obs}}ln(f),\end{array}$$

$$\ln L(A,s_1,s_2)=N_{\mathrm{obs}}[1-f+ln(f)]+\ln L(A_0,s_1,s_2).$$

Full text

supernovae: general

††software: Astropy (Astropy Collaboration et al., 2013), KCORRECT (version 4.3 Blanton & Roweis 2007), FSPS (version 3.0 Conroy et al. 2009, 2010

1 Introduction

Carbon-oxygen white dwarfs (CO WDs) are exceptionally stable, long-lived stars. However, if they are heated or compressed by some external mechanism, or if mass transfer causes their mass to approach the Chandrasekhar (1931) limit, the result can be a thermonuclear runaway, consuming the entire star and leaving no compact remnant behind (for a review, see Livio & Mazzali 2018).

What triggers the explosion of an inert CO WD as a SN Ia remains a point of contention. It is currently accepted that a binary system is required, but neither the stellar evolutionary path (sometimes referred as the “SN Ia channel”), nor the explosion mechanism have been identified. Unlike core-collapse supernovae, which originate in massive and hence luminous stars, no direct direct imaging of the progenitor systems of SNe Ia has been obtained.

In this paper, we attempt to constrain the delay time distribution (DTD), which is defined as the SN Ia rate for a simple stellar population (SSP) per unit stellar mass and as a function of age. This quantity is relevant because proposed progenitor systems of SN Ia predict distinct DTDs (e.g. Greggio 2005; Toonen et al. 2014); therefore one may gain insight into the correct channel by constraining the parameters of this distribution.

For instance, in double-degenerate scenarios (DD; e.g. Tutukov & Yungelson 1981; Webbink 1984), the expected rates are determined primarily by the distribution of orbital separations after a common-envelope phase (Ivanova et al., 2013) and by the corresponding timescales for orbital decay due to gravitational radiation, leading to a generic prediction a power–law like DTD with a slope of $-1$. On the other hand, in most single-degenerate scenarios (SD; e.g. Whelan & Iben 1973; Nomoto et al. 2000), the predicted rates are sensitive to a number of parameters, including the mass transfer efficiency, the distribution of binary mass ratios, etc. (e.g. Greggio 2005). Still, in many of those scenarios, supernovae are produced only when the white-dwarf companions are relatively massive, which implies DTDs with relative low rates at late times.

It should be noted that many SNe Ia have peculiar signatures, leading to events being classified into subclasses (e.g. Branch et al. 2009) which in turn may have distinct channels. For our main analysis, we are interested in subluminous, normal and overluminous events, which form the bulk of SNe Ia and which may result from a common progenitor system (e.g. Nugent et al. 1995; Heringer et al. 2017b, but see also Wang et al. 2013, Childress et al. 2015 and Polin et al. 2018).

The true form of the DTD may be complicated, but based on generic model predictions (e.g. Greggio 2005; Ruiter et al. 2011), we assume a broken power-law,

[TABLE]

This equation has three free parameters: a scale factor $A$ (in units of $[{\rm dex}(M_{\odot}^{-1}\,\rm{yr}^{-1})]$) and two slopes, $s_{1}$ and $s_{2}$, at young and old ages, respectively. The transition between slopes occurs at a “cutoff” time $t_{\rm{c}}$, and the rate is continuous in value at $t_{\rm{c}}$, such that $B=A\,(t_{\rm{c}}/1\,\rm{Gyr})^{s1-s2}$. The DTD is null prior to the time required for the first WDs to form ($t_{\rm{WD}}$, which we also refer to as the “onset” time). Following Heringer et al. (2017a)[hereafter H17], we typically adopt $t_{\rm{WD}}=100$ Myr and $t_{\rm{c}}=1$ Gyr; these values are based on constraints on the lifetime of the WD or its companion.

Below, we summarize a few of the techniques that have been used to probe the DTD:

• •

Cluster evolution: Clusters of old galaxies at different redshifts are surveyed for SNe Ia, thus allowing one to measure the SN Ia rate at multiple lookback times (e.g. Maoz et al. 2010). Under the assumption that the galaxies can be described as SSPs, one can recover rates as a function of cluster age. This method assumes that the galaxies in the sample formed at nearly the same redshift ($z\sim 3$–$4$, e.g. Andreon et al. 2016), and therefore only the late time DTD is probed ($t\gtrsim 3.2$ Gyr, Friedmann & Maoz 2018). This method also assumes that the observed supernovae did not originate from residual young populations (Schawinski, 2009).

• •

Remnants: Maoz & Badenes (2010) used the local star formation rate (SFR) near the location of SN Ia remnants in the Large and Small Magellanic Cloud to infer rates. The advantage of this technique is that it makes use of resolved stellar populations. However, the data has to be binned in time (owing to the small number of events); in addition, it is difficult to recover many SNe Ia from old populations with this method.

• •

Volumetric rates: One can measure volumetric rates up to a fixed redshift and derive a DTD by assuming a cosmic SFH (e.g. Graur et al. 2011; Perrett et al. 2012; Frohmaier et al. 2019).

• •

SFH reconstruction (SFHR): Given a large sample of galaxies that are surveyed over a few years, this method attempts to recover the SFH of each galaxy via spectral fitting. Once the distribution of masses as a function of (binned) ages is obtained, one can statistically infer the slope of the DTD based on which galaxies have hosted SNe Ia and which have not (e.g. Brandt et al. 2010; Maoz et al. 2011). This method is discussed in detail in §3.

• •

Color-luminosity (CL): H17 noted that, when expressed in terms of color and luminosity, specific rates are nearly independent of SFH. This allows one to constrain the DTD parameters from a large sample of galaxies without reconstructing their SFHs. This method has the advantage of not requiring binned data, but has poor sensitivity to rates at early times. This method is also described in detail in §3.

Some recent results suggest that the DTD is well-constrained, with a power-law slope close to $s=-1$; this is often interpreted as supportive of the DD channel (Totani et al., 2008; Maoz et al., 2010; Graur et al., 2011; Maoz & Mannucci, 2012; Sand et al., 2012; Graur & Maoz, 2013; Graur et al., 2014; Rodney et al., 2014). However, other works have indicated steeper slopes; in particular, rates assessed via the cluster evolution method suggest both a steeper slope and a higher production efficiency111Number of SN Ia per unit mass expected during a Hubble time; see §4.3. (e.g. Maoz & Graur 2017; Friedmann & Maoz 2018). Also, when applied to a similar sample of galaxies, the CL method indicated a steeper slope than the SFHR method.

This paper has three parts. First, we compare the SFHR and CL methods; D. Maoz has kindly provided us with the data used in Maoz et al. (2012)[hereafter M12], thus allowing a direct comparison. Second, we use the CL method to provide an up-to-date characterization of the DTD, including its normalization. Third, we investigate whether the production efficiency of SN Ia and the slope of the DTD need to be different in field and cluster galaxies.

This paper is organized as follows: in §2 we describe the samples, all drawn from the SDSS supernova survey. §3 introduces both the CL and the SFHR methods, and in §4 we present the results from both of these approaches. We evaluate some of the systematic uncertainties in the CL method in §5. Our conclusions are presented in §6, and we summarize our work in §7.

2 Data

In order to constrain the DTD, it is necessary to construct a well-defined sample of galaxies for which multiple SNe Ia are observed. As with M12 and H17, we use data collected by the SDSS-II Supernova Survey (Frieman et al., 2008), which covered a $\sim$300 deg2 region of Stripe 82 ($-60^{\circ}\lesssim\alpha\lesssim 60^{\circ}$, $-1\fdg 25\lesssim\delta\lesssim 1\fdg 25$) over nine months during 2005-2007. Specifically, we use the final compilation of the observed SN Ia and their matched host galaxies of Sako et al. (2018)[hereafter S18].

Starting with the sample of galaxies, a variety of data cuts (e.g. in limiting redshift and/or magnitude) are used to optimize both the completeness and purity of the sample. Additional cuts may be adopted to remove spurious or undesired objects, such as stars or QSOs.

Also important is the selection of the SN Ia sample. This is a non-trivial task because there is a range of reliability for reported events: Some have poorly sampled lightcurves, while others have no or low signal-to-noise (S/N) photometry. Ideally, one wants to maximize the number of SNe Ia, while making sure that the contamination by other events, such as core collapse SNe, is minimized (e.g. Dilday et al. 2010).

The S18 SNe Ia come in three flavours: (i) those which are confirmed via spectroscopy and whose hosts also have spectra [SNIa], (ii) those which are typed photometrically (e.g. Niemack et al. 2009) and whose hosts have spectra [zSNIa], and (iii) those which are typed photometrically and whose hosts do not have spectra [pSNIa]. As in M12 and H17, we only include galaxies that have been observed spectroscopically.

Furthermore, identifying the host galaxy can be complicated. M12 defines the host as that galaxy whose center is closest in angular separation to the SN Ia. Alternately, Sullivan et al. (2006) (see also Gao & Pritchet 2013; Sako et al. 2018) measure an isophotal radius for all candidate hosts, and find the galaxy with the smallest ratio of angular separation to isophotal radius. In both cases, the redshifts of the SN Ia and its host must be consistent.

Since we will be comparing our results to M12 and H17, who have also used SNe Ia from the SDSS-II Supernova Survey, we next discuss selection criteria and data cuts in these works; in addition, we describe a “standard sample” which will be used throughout this paper.

2.1 M12

M12’s sample of galaxies is derived from that of Brandt et al. (2010), which covers most of Stripe 82, with $-51^{\circ}<\alpha<59^{\circ}$ and $-1\fdg 25\leq\delta\leq 1\fdg 25$, and contains $\sim$83 000 galaxies with spectroscopic data. Of these, about 77 000 remained after excluding active galactic nuclei and QSOs. M12 further trimmed this sample to $\sim$66 400 galaxies by imposing three conditions: (i) A lower limit in velocity to remove possible misidentified stars; (ii) An upper limit of $\alpha<57^{\circ}$ to match the sky region adopted by Dilday et al. (2010)[hereafter D10], from which M12 retrieved part of their SN Ia sample. (iii) A redshift limit of $z<0.4$ based on the D10 efficiency estimates.

M12’s sample of 312 SN Ia is retrieved from D10 (spectroscopic SNe) and Sako et al. (2011) (photometric SNe), and includes objects with at least one photometric measurement with S/N $>5$ in the $g$, $r$ and $i$ bands, photometric observations near and after maximum, and light curves well-fitted by a standard template. For host identification, M12 select host candidates with consistent redshifts and projected physical separations less than 30 kpc, and, in the case of multiple candidates, select the one with the smallest projected angular separation. They find 61 matches.

In the redshift range of interest for this work ($0.01\leq z\leq 0.2$), there are 47 728 galaxies and 141 SNIa222D10 list 140 in their Table 1, but there is one additional without an IAU designation., of which 51 have hosts matched by M12.

2.2 H17

In H17, the sample of galaxies contained all objects with spectra for which $-51^{\circ}\lesssim\alpha\lesssim 59^{\circ}$, $-1\fdg 258<\delta<1\fdg 258$, $0.01<z<0.2$ and $14<r<17.77$ (after correction for Milky Way extinction). The redshift criterion ensured reliable SN Ia identification (Dilday et al., 2010; Sako et al., 2018), while the host magnitude cut ensured a nearly complete magnitude-limited sample of galaxies with spectroscopy (Strauss et al., 2002). H17 imposed further cuts of $g<22.2$ (after extinction correction) and $g_{\rm{err}},\,r_{\rm{err}}<0.05$ to limit the number of spurious objects, leaving 20 707 galaxies.

For the SNe Ia, H17 retrieved 53 hosts from Gao & Pritchet (2013), which were found using the Sullivan et al. (2006) technique of minimizing the ratio of angular separation to isophotal size. An upper limit of 3.8 to this relative distance was enforced, below which contamination was estimated to be small ($<8\%$). In addition, 3 events observed during the survey engineering time in 2004 were included (we will not use these in the present work, since the survey cadence and duration in 2004 were not adequate for a study of rates; Sako et al. 2018). All these events are classified as SNIa in the S18 terminology.

Because the method H17 employed to estimate rates relied on colors, another cut was made by imposing $-0.4\leq\Delta(g-r)\leq 0.08$, where $\Delta(g-r)$ is the galaxy color relative to the color of the red sequence; this reduced the number of galaxies and hosts to 17 539 and 52, respectively.

2.3 Comparison between M12 and H17 samples

Under the assumption of a single power-law DTD, M12 and H17 derived slopes that were not in agreement. It is therefore important to understand how their data samples differed.

M12’s sample is much larger than H17’s: 66 400 vs 17 539 galaxies. This is mostly because of the different redshift ($0.4$ vs $0.2$) and magnitude limits (no cut vs $14<r<17.77$): with the redshift and photometric cuts of H17 and, for completeness, the slightly more restrictive sky region of M12, the M12 and H17 samples would contain 19 329 and 17 253 galaxies, respectively, with 15 381 galaxies in common.

Another significant difference is the data source. M12’s objects were retrieved from the SpecPhotoAll Table database, while H17’s were retrieved from the SpecPhoto View database, both for the DR7 data release333Queries in this database were performed with https://skyserver.sdss.org/CasJobs/default.aspx.. While the implications of each choice are not obvious, it is likely that the Table database includes additional imaging in the Stripe 82 region besides that of the Main Galaxy Sample (Strauss et al. 2002, see Brandt et al. 2010 for more details.).

Yet another difference is the sample of SNe Ia hosts in each work: in the redshifts of interest, M12 and H17 identified 51 and 52 hosts, respectively. While these numbers seem consistent, there are only 20 SNe Ia in common.

The majority of SNIa that are present in H17 but not in M12 are also not present in the SNIa sample of D10 (which M12 is based on; see above). We are not sure why this is the case. For example, the D10 list misses even the spectroscopically confirmed SNe Ia SN2006fd and SN2006er, both of which had 15 observed epochs with S/N $>$5. Regarding the SNIa that are present in M12 but not present in H17, these are all, except 5, explained by the photometric cuts enforced of the galaxy sample in H17 and would have otherwise been included.

2.4 This work

We now define a “standard sample”, based on the galaxies used in M12, which we will use to derive our most likely DTD parameters and to study systematic uncertainties. We summarize the selection criteria for this dataset below.

$-51^{\circ}<\alpha<57^{\circ}$ and $-1\fdg 25<\delta<1\fdg 25$. 2. 2.

$0.01<z<0.2$. 3. 3.

$g_{\rm{err}},\,r_{\rm{err}}<0.1$ 4. 4.

Exclude galaxies with duplicate SDSS DR7 IDs. 5. 5.

Exclude SNe Ia flagged as peculiar in S18.444Classified as sn00cx, sn02ci or sn02cx, following the notes given in Table 4 of S18.

For consistency, the choice of sky region follows that of M12 and we relax the photometric limits imposed in H17; However, given the considerations in §2.3, we adopt the redshift cuts of H17, while still enforcing an upper limit to photometric errors. The list of both spectroscopically and photometrically typed SNe Ia and their respective hosts is taken from S18555The hosts in S18 are indicated with a SDSS DR8 ID, which, when possible, we match to a SDSS DR7 by querying under the DR8 context, using the PhotoPrimaryDR7 Table.. In the redshift range of interest, S18’s list contains a total of 215 SNIa and 70 zSNIa. There are only 3 (5) SNIa in the lists of M12 (H17) that are not present in the S18 compilation.

In total, this sample contains 43 895 galaxies and 94 SNe Ia, out of which 2 SNe Ia are of peculiar type and therefore not included, leaving 76 SNIa and 16 zSNIa. A color-magnitude diagram of this sample is shown in Fig. 1.

Before continuing with the above sample, we note that we have tested different combinations of the galaxy samples from M12 and H17 and the hosts from M12, H17 and S18. As discussed in Appendix D, when the same data cuts are applied, we find DTD slopes that are in agreement within the derived uncertainties, though the normalizations differ.

3 Methods

The SN Ia rate per unit of mass is, by definition, calculated as

[TABLE]

where $t_{\rm{WD}}$ is the time at which the first white dwarfs appeared (typically $t_{\rm{WD}}=100\,\rm{Myr}$, unless we are comparing with M12’s results, in which case we adopt their value of $t_{\rm{WD}}=40\,\rm{Myr}$), and the form of the DTD is given by Eq. 1, with a cutoff time $t_{\rm{c}}=1\,\rm{Gyr}$, unless otherwise stated.

From Eq. 2, one sees that for a SSP, for which the SFR approaches a delta function, the expected rate will be equivalent to the DTD itself. For the general case, however, the situation is more challenging: the SN rate of each galaxy depends on its SFH.

In this paper, we adopt and compare two methods of analysis: one based on galaxy colours with respect to the red sequence (the CL method), and one based on inferred galaxy masses (the SFHR method), following the prescriptions given in H17 and M12, respectively. Below we summarize each method, including possible adjustments, and discuss their respective strengths and shortcomings.

3.1 The Color-Luminosity (CL) Method

The CL method of H17 uses the fact that the specific supernova rate per unit luminosity ($sSNR_{L}$) for a galaxy of a given color is quite strongly dependent on the assumed DTD, but nearly independent of the SFH. In this method, we first use a stellar population code (FSPS; Conroy et al. 2009, 2010) to derive luminosity and color vs. age for some assumed SFH666The SFH is normalized such that 1 M⊙ is formed over $\sim 14\ \rm{Gyr}$., and then compute rates as a function of age via Eq. 2. These results are combined to give a normalized rate per unit luminosity ($sSNR_{L}$) as a function of color.

The FSPS stellar population is characterized by properties such as the metallicity, initial mass function, etc. These are similar to the assumptions made when reconstructing SFHs and are discussed below (§ 5) as part of the systematic uncertainty budget. To mitigate the impact of some of these unknown parameters, colors are calculated with respect to the red sequence, which is computed as the color of an SSP with age 10 Gyr. As discussed in H17, we choose to use the $g-r$ color as the locus of the red sequence is clear in these filters and the uncertainties are typically smaller than in $u-r$.

Fig. 2 shows predicted $\Delta(g-r)$–$sSNR_{L}$ relations. Two general forms for the SFH are tested: an exponential case, where $\rm{SFR}(t)\propto e^{-t/\tau}$ (left panel) and a delayed exponential case, where $\rm{SFR}(t)\propto t^{-1}\,e^{-t/\tau}$ (right panel), with timescales ranging from $\tau=1$ to $10$ Gyr. One sees that the specific rates are nearly independent of the SFH, but are strongly dependent on the DTD. H17 showed this to be true also for more complicated SFHs. The CL method is most sensitive to the slope at later times, $s_{2}$; to constrain the slope at early ages ($s_{1}$), one requires very blue galaxies, with $\Delta(g-r)\lesssim-0.8$.

For simplicity we ignore (as did H17) the effects of marginalizing our rate predictions over the uncertainty in $\Delta(g-r)$ and $r$-band magnitude. This is justified where the power-law slope of the DTD is shallow, because the supernova rates are nearly independent of $\Delta(g-r)$ (see Fig. 2) and the $r$-band error is typically small. On the other hand, for steeper power-law slopes, one sees larger scatter around the $\Delta(g-r)$–$sSNR_{L}$ relation, which likely dominates the uncertainty of these models.

To infer likelihoods for the DTD parameters from the observations, we follow a procedure similar to that of H17, described in detail in Appendix B.

3.2 Improvements to the H17 CL Method

• •

We define our model to be the median of the rates predicted for different SFH’s (instead of interpolating the predicted rate per unit luminosity at $10$  Gyr – the blue squares in Fig. 2). For the few galaxies that are redder than the color range in our models, we simply adopt the rate at the reddest colour, i.e. $sSNR_{L}(\Delta>\Delta_{\rm{max}})=sSNR_{L}(\Delta_{\rm{max}})$. The final $\Delta(g-r)$–$sSNR_{L}$ relations are shown as green curves in Fig. 2. With this change, we can use the full range of colours, making it unnecessary to impose cuts on $\Delta(g-r)$.

• •

For a more straightforward comparison with the SFHR method (see §3.3), our models assume a Kroupa IMF (Kroupa, 2007), rather than the Chabrier IMF (Chabrier, 2003) used in H17. (This change makes very little difference.)

• •

Originally, the CL method was used to probe only the shape of the DTD. Here we extend our models to fit the normalization parameter $A$ in Eq. 1.

• •

We improve our treatment of $K$-corrections. We use the KCORRECT package (Blanton & Roweis, 2007), which requires the input magnitudes to be converted to “maggies”, a flux unit on a linear scale. In H17, we had made this conversion as if the SDSS photometry were given in magnitudes, when in reality they are given in “luptitudes” (Lupton et al., 1999); this has been fixed. In addition, following standard routines of the KCORRECT package, we now also perform photometry corrections to the AB system777These are small: $0.012$ ($0.010$) in the $g$ ($r$)-band. and add small errors in quadrature in each band.

• •

We now take into account effective visibility times, as discussed in Appendix B (Eqs. B3 and B4), and we adopt $M_{r,\odot}=4.65$888Consistent with the FSPS value; see https://python-fsps.readthedocs.io/en/latest/filters/ (e.g. Willmer 2018), rather than $M_{r,\odot}=5$.

• •

We have also changed the cosmological parameters to make them consistent with those used to compute galaxy masses in M12 (from WMAP5; Komatsu et al. 2009). The new (old) parameters are $H_{0}=70.5$ ($67$) km s*-1* Mpc*-1*, $\Omega_{\Lambda}=0.726$ ($0.68$) and $\Omega_{m}=0.274$ ($0.32$). We have ensured that these parameters are consistently used in all our calculations.

Several of the modifications described above are necessary for this method to estimate the scale factor $A$. By comparing stellar masses with those in an independent study (Chang et al., 2015), we argue in Appendix C that such a measurement is viable.

We have redone the analysis in H17 after incorporating the improvements above, finding that the H17 conclusions remain unchanged. The routines used to produce the $\Delta(g-r)$–$sSNR_{L}$ relationships are publicly available999https://github.com/Heringer-Epson/ssnral..

3.3 The SFH Reconstruction (SFHR) Method

This method measures the DTD by estimating masses and ages for each galaxy in the dataset. This procedure is similar to that adopted by Brandt et al. (2010) and by Maoz et al. (2011), but we focus on the work of M12.

Masses in each of a set of age bins (defined in Eq. A1) are estimated using the VESPA code, by matching a galaxy’s spectrum with templates for stellar populations with a Kroupa (2007) IMF and either an exponential or a dual-burst SFH (Tojeiro et al., 2007). M12 uses the masses made available by Tojeiro et al. (2009).

M12 combined the VESPA masses in three age bins: 0–0.42 Gyr, 0.42–2.4 Gyr and $>\!2.4$ Gyr. They also scaled all masses by a factor of 0.55 to account for light outside of SDSS fibers, etc., verifying that their final masses are consistent with other, independent estimates.

For each galaxy $i$ in their sample, M12 assign an expected supernova rate summed over the three age bins,

[TABLE]

where the numeric indices denote the three time bins. The rates $\psi$ are estimated using the same likelihood calculation as for the CL method, summarized in appendix B, except that in Eqs. B6–B12, the likelihood is kept as a function of $\psi_{1-3}$. In order to determine the slope of the DTD, M12 fitted a power-law function to the most likely binned rates $\psi$.

Following M12’s analysis101010Eq. 7 of M12 contains a typo: the factor $t_{i}^{2}$ should be $(\epsilon t_{i})^{2}$. we derived maximum likelihood rates that agreed with their Table 2 to within $1\%$. Fitting a DTD slope to these rates yielded slightly different results: $s=-1.23\pm 0.19$, rather than $s=-1.12\pm 0.08$. In any case, we show below (in §3.5) that one can assess the DTD most likely parameters directly from the Bayesian analysis, thus avoiding the additional step of fitting a power-law to rates $\psi$.

3.4 Comments on the SFHR Method

While the rate recovery method is reliable, one needs to make important assumptions in order to derive DTD parameters from those:

i) The rate in each bin was assumed to be representative of a single age, which was taken to be the linear mean age of the bin (Fig. 1 of M12).

ii) The observed rates were approximated as a binned mass multiplied by the intrinsic DTD rate (see Eq. 3). It is difficult to quantify how much this approximation affects the DTD slope, given that it bypasses the convolution in Eq. 2, and given that only three age bins were used.

iii) VESPA masses in the youngest bin are overestimated, because they include stars in the age range [math]–$70$  Myr which cannot be distinguished from later ages ($70$–$420$  Myr).

iv) The SFHR method uses only 3 age bins, and hence depends strongly on the assumed onset time of SN Ia, $t_{\rm{WD}}$, which is poorly known. Note that CL method is also sensitive to this variable, and its effects are explored in §5.

To illustrate the point that fitting a power-law to binned rates may not be an optimal way to recover DTD parameters, and that the the fitted power-law depends on the assumed onset time, we calculate mock mean rates as

[TABLE]

from a known DTD with an assumed slope $s=-1$ and a scale factor $A=10^{-12}$ ${\rm dex}(M_{\odot}^{-1}\,\rm{yr}^{-1})$; $t_{i,j}$ and $t_{f,j}$ are the initial and final ages of the $j$-th bin. These rates are then fitted by a power-law, as shown in the top panel of Fig. 3, for three choices of $t_{\rm{WD}}$ (color coded). The bottom panel shows the recovered slope and scale factor.

We find that the recovered slope is sensitive to the choice of $t_{\rm{WD}}$; the scale factor is over-predicted in all cases. We also checked that better agreement would be achieved with a larger number of bins. This exercise is not entirely equivalent to M12’s procedure, since here the mean DTD rates have been taken directly from the input DTD, without computing binned masses (which would add to the uncertainties). Therefore we will not use this fitting procedure.

3.5 Adjustments to the SFHR Method

To avoid some of the problems associated with fitting a power law to recovered rates, we compute likelihoods which depend directly on the DTD parameters, by expressing rates as $r_{i}(t_{i}|A,s_{1},s_{2})$ rather than $r_{i}(t|\psi_{i,1},\psi_{i,2},\psi_{i,3})$. It is still necessary to use the approximation that the expected rates in each galaxy depend on binned quantities. Furthermore, we still need to assign an age to each bin; for consistency with M12, we will approximate this as the linear mean in each bin, $\bar{t_{j}}$. In summary, we replace Eq. 3 with

[TABLE]

Note that this formulation allows us to fit a broken power law ($s_{1}\neq{}s_{2}$), but the rates in the intermediate age bin, for which $0.42<t<2.4$ Gyr, become awkward if the cutoff time falls within it, as is the case for our default choice of $t_{\rm{c}}=1\,$Gyr. For the comparison below (§4.1), we simply try cut-off times at the start and end of the intermediate bin.

The bin age limits and mean ages and an expanded form of Eq. 5 are given in appendix A; Likelihoods are computed following the description given in appendix B.

3.6 Method Comparison on Mock Data

We used FSPS to simulate the color and luminosity of galaxies according to a set of random input parameters, which consists of an SFH timescale, a formed stellar mass, and a galaxy age. With a known DTD slope and scale factor, we used Eq. 2 to compute the true rate for each galaxy, based on which a random sample of host galaxies was drawn. We computed the expected rates for each galaxy and confidence contours for the DTD parameters using both the CL and SFHR methods, finding good agreement with the true values, noting that some of the approximations in the SFHR method (see §3.5) cause the rates to be underestimated by $\sim 10\%$ and the derived scale factor to be overestimated at a similar level. This suggests that discrepancies for the best DTD parameters when these methods are applied to real data are likely due to possible systematic errors in either or both techniques, such as the mass estimation in the SFHR method or by having colors computed relative to the red-sequence in the CL method.

4 Results

4.1 Method Comparison

As a validation exercise, we first compare the likelihoods derived for the modified SFHR and CL methods. For this we use our “standard” sample, while adopting $t_{\rm{WD}}=40$  Myr, as in M12.

We note that since we are using a slightly different sample, the fitted DTD parameters will not be identical to the M12 and H17 results. We also caution that the only parameter that can be calculated and compared using the original methods of M12 and H17 is the continuous DTD slope $s=s_{1}=s_{2}$. The calculation of the scale factor $A$ and the $s_{1}$–$s_{2}$ likelihoods111111As shown in Appendix B, one does not need to marginalize the $s_{1}$–$s_{2}$ likelihood over the scale factor because the integration over $A$ does not depend on the choice of slopes. arise from adjustments discussed in §3. For the latter calculation, we try cutoff ages at the start and end of the intermediate bin (see §3.5), while the standard case for which $t_{\rm{c}}=1$ Gyr is discussed below.

The left panel of Fig. 4 shows the $A$-$s$ parameter space, calculated via the CL method (green) and the SFHR method (orange). The most likely parameters are retrieved directly from the posterior distribution using Eq. 4. If a continuous slope is assumed, then $s=-1.25^{+0.16}_{-0.15}$ ($s=-1.24^{+0.17}_{-0.16}$) when applying the CL (SFHR) method; these slopes are in excellent agreement. The scale factors however differ at the $>2\sigma$ level: $\rm{log}\ A=-12.24^{+0.09}_{-0.11}$ ${\rm dex}(M_{\odot}^{-1}\,\rm{yr}^{-1})$ for the CL method and $\rm{log}\ A=-12.48^{+0.07}_{-0.09}$ ${\rm dex}(M_{\odot}^{-1}\,\rm{yr}^{-1})$ for the SFHR method. These results indicate a discrepancy by a factor of $\sim 2$ for the most likely DTD normalization; we discuss possible sources for this difference in §6.

The $s_{1}$–$s_{2}$ parameter space is shown in the middle and right panels of Fig. 4, according to the choice of cutoff time. One can see that neither method is capable of constraining well the early slope if $t_{\rm{c}}=0.42$ Gyr, because of the small fraction of blue galaxies. At the $68\%$ confidence level, we find $-1.46<s_{2}<-0.45$ ($-1.60<s_{2}<-0.72$) for the CL (SFHR) method. If $t_{\rm{c}}=2.4$ Gyr, then, for both methods, the early slope becomes better constrained, while a larger range of $s_{2}$ is acceptable. We find $-2.23<s_{1}<-0.94$ and $s_{2}>-1.54$ ($-1.90<s_{1}<-0.93$ and $s_{2}>-1.96$) for the CL (SFHR) method.

Overall, we conclude that the two methods give consistent information on the shape of the DTD, but not on its normalization.

4.2 Results From the Default Sample

We use our “standard” sample (defined in §2.4) to characterize the most likely DTD parameters using the CL method. The confidence contours are shown in Fig. 5 differ somewhat from those in Fig. 4 because we use our preferred $t_{\rm{WD}}=100$ Myr.

Assuming a continuous slope ($s=s_{1}=s_{2}$; left panel), we find $\rm{log}\ A=-12.15^{+0.10}_{-0.13}$ ${\rm dex}(M_{\odot}^{-1}\,\rm{yr}^{-1})$ and $s=-1.34^{+0.19}_{-0.17}$. If instead the slopes at early and late times are allowed to vary, then $s_{1}$ is poorly constrained due to the small number of SNe Ia in blue galaxies. As a consequence, $s_{2}$ will also span a large range and, at the $68\%$ confidence level, we find $s_{1}<-0.6$ and $-1.54<s_{2}<-0.30$. Note that the scale factors (dotted lines) are not the same across this parameter space, but vary such that the predicted number of SNe Ia matches the observed one.

4.3 Production Efficiency

The production efficiency is the number of SN Ia per unit mass expected from a burst of star formation over a Hubble time (e.g. Friedmann & Maoz 2018). It is defined as

[TABLE]

This quantity is more physical than either $A$ or $s$, is less sensitive to correlations between $s$ and $A$, and is better for comparing DTDs with different shapes.

By applying CL (SFHR) method to our default sample under the assumption of a continuous slope $s$, we find $PE=0.004^{+0.002}_{-0.001}\,M_{\odot}^{-1}$ ($PE=0.002^{+0.001}_{-0.001}\,M_{\odot}^{-1}$), noting these values do not change appreciably if $t_{\rm{WD}}=40\,\rm{Myr}$ or $t_{\rm{WD}}=100\,\rm{Myr}$. We defer a discussion of these values to §6.

5 Model uncertainties

In this section we examine how different model assumptions might impact the likelihoods for our standard sample. We use the CL method and focus on the $A$-$s$ parameter space, as its likelihood is much better constrained. The “fiducial” model parameters include: a Kroupa IMF (Kroupa, 2001), an exponential SFH, solar metallicity ($Z=0.019$), no blue horizontal branch stars ($f_{BHB}=0$), the BASEL spectral library (Lejeune et al., 1997, 1998; Westera et al., 2002), the PADOVA isochrone library (Marigo & Girardi, 2007; Marigo et al., 2008), an evolutionary factor $Q=1.6$ (Wyder et al., 2007), and a white dwarf onset time $t_{\rm{WD}}\ =100$  Myr.

We perform several tests, changing the standard model parameters above one at a time. The changes for these tests include: a Chabrier (Chabrier, 2003) or Salpeter (Salpeter, 1955) IMF, a delayed-exponential SFH, a metallicity of $Z=0.015$ or $Z=0.03$, $f_{BHB}=0.2$, the MILES spectral library (Sánchez-Blázquez et al., 2006), the MIST isochrone library (Dotter, 2016), $Q=0$, and $t_{\rm{WD}}\ =40$  Myr.

For reference, in Fig. 6 we show the confidence contours obtained for the fiducial model parameters. This confidence contour is similar in shape and size for all of our tests, and so the figure shows only the changes in the most likely value of the fit for each test. Furthermore, for clarity we show only those tests which result in a significant change in slope or scale factor ($\Delta s\geq 0.03$ or $\Delta\rm{log}\ A\geq 0.05$). In general we find from the tests that the systematic uncertainies are $\Delta s\sim 0.1$ $\Delta(\rm{log}\ A)\sim 0.2$. In more detail:

Relative to a Chabrier or Kroupa IMF, the Salpeter IMF has more low-mass stars (bottom-heavy distribution). Other works have indicated that galaxy masses derived assuming a Salpeter IMF tend to be roughly 1.4 times larger than assuming a Chabrier or Kroupa IMF (da Cunha et al., 2008; Maoz et al., 2012). This explains the lower scale factor (orange segment) seen in Fig. 6.

Changing the metallicity has little effect on the shape of the $\Delta(g-r)$–$sSNR_{L}$ relation shown in Fig. 2, but can lead to a small offset $\Delta(\rm{log}\ sSNR_{L})\sim 0.1$ in the rate at $\Delta(g-r)\sim 0$. For instance, lower metallicities will shift the confidence contours upwards, such that these lower rates are compensated for by slightly shallower slopes (see the -1/-1 and -1.5/-1.5 cases in Fig. 2). Also, since the cumulative rates become smaller, a slightly higher scale factor will be preferred. This exercise suggests that the CL method is sensitive to the precise $sSNR_{L}$ for colors redder than the RS; the assumption of a constant specific rate in that regime may be an oversimplification.

The effect of adopting a shorter $t_{\rm{WD}}$ is that $sSNR_{L}$ rates become higher for bluer galaxies, but are unchanged for redder objects. As a result, the cumulative rates will be higher, thus leading to a smaller scale factor $A$. Moreover, a shallower slope is preferred because it would boost the relative rates expected for red galaxies, thus compensating for the change in rate for the blue galaxies due to a shorter $t_{\rm{WD}}$.

5.1 Dust

Our stellar population models are calculated using the FSPS package (Conroy et al., 2009, 2010), for which different choices of dust parameters are allowed. As in H17, we adopt the default mode for dust absorption (dust_type=0), which follows the two-component model of Charlot & Fall (2000). This model assumes that all stars are affected by an uniform dust screen and, in addition, young stars also suffer from a power-law attenuation curve. This dust treatment is similar to that adopted in VESPA (see §2.2.2. of Tojeiro et al. 2009). We performed additional tests with different choices of attenuation curves, finding that the most likely slope and scale factor did not change appreciably.

6 Discussion

We have demonstrated that both the CL and the SFHR methods yield approximately the same power-law slope, $s\simeq-1.35$ for the DTD if no break is assumed. This value is inconsistent with that reported by both M12 ($s\sim-1.1$) and H17 ($s\sim-1.5$), mainly because of the difference between adopted samples (see §2.3). Under the assumption of $t_{\rm{WD}}=40$ Myr, if we relax the photometric and redshift limits so that the CL (SFHR) method is applied to the entire sample of M12, we obtain $s=-1.02^{+0.15}_{-0.13}$ ($s=-1.18^{+0.13}_{-0.13}$).

The scale factors derived from both methods are marginally inconsistent, at the $>2\sigma$ level, with the SFHR method suggesting a lower normalization, by a factor of $\sim$2. Both methods are applied to the same sample and thus differences in the data cannot explain this discrepancy. In Appendix C we find that, for a subset of red-sequence galaxies, the mass scaling based on FSPS is consistent with independent source. M12 also finds the VESPA masses to be consistent with another independent source. In § 3.6, we have verified that both methods are able to accurately recover the DTD parameters, which suggests that the scale factor discrepancy must arise from systematic uncertainties in either or both methods. For instance, VESPA masses are based on magnitudes computed by integrating the SDSS spectra over multiple bands, while the CL method uses photometric magnitudes solely from the $r$-band. Also, due to the nature of the CL method, we can only verify that the masses from galaxies in the red-sequence are consistent, but it is possible that this agreement does not hold for bluer galaxies.

The slope derived by M12 is in good agreement with those estimated in volumetric surveys (Perrett et al., 2012; Frohmaier et al., 2019), which yielded $s\sim-1$. However, we show that when the SFHR method is applied to a more conservative sample (i.e. $z<0.2$), a steeper slope is recovered. Moreover the production efficiency are not in agreement: Perrett et al. (2012) obtained ($4.4\pm 0.2$)–($5.2\pm 0.2$) $10^{-4}\,M_{\odot}^{-1}$, while we obtain value roughly a factor of 10 (4) larger using the CL (SFHR) method. It is intriguing to us that the cluster and volumetric measurement of the DTD could yield such discrepant production efficiencies, while the SFHR and CL methods seemingly recover intermediate values.

6.1 DTD Dependency on Environment

We find that our results are in marginal agreement with recent inferences of both the slope and the scale factor of the DTD in cluster galaxies. In particular, Friedmann & Maoz (2018) obtained $s=-1.45^{+0.51}_{-0.38}$ and $\mathrm{log}\,A=-12.00^{+0.34}_{-0.52}{\rm dex}(M_{\odot}^{-1}\,\rm{yr}^{-1})$, with a production efficiency of $0.009^{+0.029}_{-0.007}\,M_{\odot}^{-1}$ (see their Table 7121212. For fairness of comparison, we use the result that assumes no knowledge of the cluster iron content, and use the values derived using a single redshift bin; the values derived using two bins are similar.), while for our standard sample of field galaxies, the CL method yields, for the same assumed $t_{\rm{WD}}=40$ Myr, $s=-1.25^{+0.16}_{-0.15}$ , $\rm{log}\ A=-12.24^{+0.09}_{-0.11}$ ${\rm dex}(M_{\odot}^{-1}\,\rm{yr}^{-1})$  and $PE=0.004^{+0.002}_{-0.001}\,M_{\odot}^{-1}$. In contrast, previous production efficiencies derived for field galaxies were much lower; for instance, M12 obtained $PE=0.00130\pm 0.00015\,M_{\odot}^{-1}$ (reflecting the lower scale factor $A$ discussed above) and Maoz & Graur (2017) derived $PE=0.0013\pm 0.0001\,M_{\odot}^{-1}$.

Our higher production efficiency may influence conclusions about abundances. For instance, Maoz & Graur (2017) have proposed that the spatial distribution of iron and $\alpha$-elements in the Milky Way (MW) may be explained by a simple model: the halo and thick disk consist of an old population, formed in a burst at $z\sim 3$, with the steeper DTD slope and higher SN Ia production efficiency of cluster galaxies, while the thin disk would be explained by a more extended SFH, with the DTD derived for field galaxies. It would be interesting to see how a different production efficiency changes these results.

6.2 DTD Scenarios

Next we consider which DTD scenario is favored using our slope and scale factor. A value of $s\sim-1.35$, assuming no DTD break, is intermediate between the generic predictions of the double-degenerate and the double detonation (Ruiter et al., 2011) scenarios. This slope may also be possible for the SD channel with an appropriate distribution of binary mass ratios (see Fig. 7 of Greggio 2005 and the findings of Moe & Di Stefano 2015), or if there exists a mechanism that delays the explosion. Therefore it is not possible to constrain the DTD channel from our analysis.

Furthermore, the DTD does not discriminate between scenarios if a break in the DTD is allowed: given the limitations of the models described in §3, only a part of the $s_{1}$–$s_{2}$ parameter space can be constrained by either the CL or the SFHR method (see Fig. 4). Both methods seem to allow for a shallower slope at young ages, and a steeper slope at later times, although it is never steeper than $s_{2}=-2$ at the $68\%$ confidence level. Using our default sample and assuming $t_{\rm{c}}=1$ Gyr, we find $-1.54<s_{2}<-0.30$.

Finally, it is possible that more than one channel is realized in nature; in this case, the effective DTD need not have a simple form. For instance, Wang et al. (2013) have suggested the possibility that high velocity (HV) SNe Ia have a distinct progenitor different from that of normal velocity (NV) events; it would therefore be interesting to repeat the analysis of this paper separately for NV and HV events.

7 Summary

We have studied the CL and SFHR methods for probing the DTD of SN Ia. The CL technique, first introduced by H17, avoids the use of stellar masses and ages, two quantities that are notoriously difficult to measure. Furthermore, this method is insensitive to the SFH ($\it cf.$ the cluster method), performs the full Eq. 2 convolution to compute rates, and uses unbinned data ($\it cf.$ the SFHR method).

In order for the CL method to fit the normalization of the DTD, we improved the photometric corrections necessary to compute absolute magnitudes. These corrections were indirectly checked against an independent source, exhibiting excellent agreement. We showed that the SFHR and the CL methods lead to a consistent DTD slope, while the scale factors disagree at the $95\%$ confidence level (using the same dataset).

By applying the CL method to our default dataset of field galaxies, and under the assumption of a continuous power-law slope, we obtained $\rm{log}\ A=-12.15^{+0.10}_{-0.13}$ ${\rm dex}(M_{\odot}^{-1}\,\rm{yr}^{-1})$ and $s=-1.34^{+0.19}_{-0.17}$. Unlike previous inferences for field galaxies, we find values that are in reasonable agreement with the production efficiency of Friedmann & Maoz 2018 for galaxy clusters. A study of systematic uncertainties indicates that the slope $s$ and scale factor $\rm{log}\ A$ can be determined to within $\sim 0.1$ and $\sim 0.2\,$dex, respectively. Based on these results, we do not favor any particular DTD scenario.

There are several ways in which the measured DTD and its interpretation could be improved. A much larger sample of SNe could be obtained by relaxing the requirement that host galaxies have spectroscopy; this would improve the statistical errors in the broken power-law slope measurements, and provide an improved constraint on proposed DTD scenarios. A larger sample would furthermore allow a test of whether SN events with different photometric and spectroscopic properties have different DTD properties. A specific example of this is the prediction that NV and HV events form two distinct populations of SN Ia (Wang et al. 2013); with a large enough sample of NV and HV events, one could test whether these events are described by distinct DTD parameters and hence evolutionary scenarios.

We thank Dani Maoz for his insightful comments and for providing us with the VESPA data used in his original work. We also thank the organizers of the “Observational Signatures of Type Ia Supernova Progenitors III” workshop, which led us to improve the methods to measure the delay times of SNe Ia. ), Python-FSPS (version 2017.07.05 Foreman-Mackey et al. 2014).

Appendix A VESPA Rates as a Function of DTD Parameters

For completeness, here we expand Eq. 5 to a form where it depends only on the DTD parameters. The mean ages are set by the domain of the age bins in M12, except that we leave the unknown early age of the first bin as a free parameter; the mean ages are then

[TABLE]

Eq. 5 therefore becomes

[TABLE]

Appendix B Determining Most Likely Parameters

Here we show in detail how our likelihoods are computed. This is similar to the implementation of Maoz & Badenes (2010), M12, Gao & Pritchet (2013) and H17.

We first assign a rate to each galaxy (denoted by the subscript $i$) in the data set. This is accomplished by using the SFHR or the CL methods:

[TABLE]

where $j$ runs through each age bin ($j={1,2,3}$) in the SFHR analysis; $\mathcal{V}$ denotes the application of the SFHR relation from Eq. 3, and $\mathcal{S}$ the application of the $\Delta(g-r)$–$sSNR_{L}$ relation, as established in Fig. 2.

Next, these rates are converted to unitless absolute rates through the following relations:

[TABLE]

where $z_{i}$ is the galaxy redshift and $t_{i}$ corresponds to the approximate observing time window for each galaxy. As in M12, we adopt

[TABLE]

while $\epsilon_{i}$ is a detection efficiency estimated by Dilday et al. (2010) and parameterized by M12 as

[TABLE]

Note that M12 also tried a slightly different detection efficiency function, which would start declining past $z=0.184$ ($z=0.166$) for star forming (passive) galaxies. We do not adopt these extra corrections because (i) the values we use for comparison (from their Table 2) do not use this modified function, (ii) these corrections are likely small for a sample trimmed at $z=0.2$, and (iii) information regarding whether a galaxy is star forming or not is not always available.

To derive likelihoods, one can treat the SN Ia rate problem as a Poisson experiment in which the expected number of events during the survey time is very small. Thus, following the approach of Cash (1979), the probability of a galaxy hosting $n_{i}$ events given an expected rate $\lambda_{i}$ is

[TABLE]

The corresponding likelihood $L$ of a given DTD model will then be the simple multiplication over all galaxies of the probabilities in Eq. B5:

[TABLE]

which can be simplified in logarithmic space to

[TABLE]

Because all galaxies in the sample have either $n_{i}=0$ (non hosts) or $n_{i}=1$ (hosts; numbered with the subscript $k$), Eq. B7 simplifies to:

[TABLE]

Given the identity $N_{\rm{exp}}\equiv\sum_{i}\lambda_{i}$, Eq. B8 becomes

[TABLE]

Next, we show that the computational cost of surveying the parameter space of ($A$,$s_{1}$,$s_{2}$) can be reduced by decoupling the scale factor $A$ and treating this parameter analytically. First, assume that the likelihood $L$ is computed for a scale factor $A_{0}(s_{1},s_{2})$, which ensures that the predicted number of events matches the observations, $N_{\rm{exp}}=N_{\rm{obs}}$. Therefore, from Eq. B9,

[TABLE]

Because the specific rates in Eq. B1 are linear with respect to $A$, one can define $A\equiv f\times A_{0}$ so that the predicted absolute rates in Eq. B2 can be written as $\lambda_{i}(A,s_{1},s_{2})=f\times\lambda_{i}(A_{0},s_{1},s_{2})$. Eq. B8 then becomes

[TABLE]

which simplifies to

[TABLE]

Eq. B12 is well behaved and, for a given pair $s_{1},s_{2}$, reaches a maximum when $f=1$ (i.e. $A=A_{0}$). This demonstrates the expected relation that the maximum likelihood occurs when $N_{\rm{exp}}=N_{\rm{obs}}$ and that marginalizing over the scale factor parameter would only correspond to a constant factor that does not depend on $s_{1}$ or $s_{2}$.

Appendix C Mass Comparison

Here we investigate whether the CL method is able to produce reliable estimates for the scale factor $A$. For this purpose we need to ensure that our photometric corrections and absolute magnitudes are realistic. We do so in an indirect way, via mass comparisons with VESPA masses, and with masses obtained in an independent study (Chang et al., 2015). Masses can be reported as “present-day” ($M_{\ast}$) or “formed” stellar masses ($M_{\rm{T}}$). Taking a SSP as an example, its present-day stellar mass will decrease with time, but its formed mass remains unchanged.

The CL method does not directly assess galaxy masses to constrain the DTD parameters because the SFH itself is not fixed. However, given the computed absolute magnitude of a galaxy, a mass range can be estimated if an SFH is assumed. For this exercise, we select a subset of galaxies that belong to the red sequence by imposing $-0.1<\Delta(g-r)<0.1$ and by requiring VESPA masses in the young bins to be negligible ($m_{i,1}=m_{i,2}=0$). We further assume that these objects may be described by an exponential SFH with a short timescale of $\tau=1$  Gyr. We assign an age of 10  Gyr to these galaxies and estimate their masses via the mass-to-light ratio predicted by FSPS. A mass uncertainty is computed by assigning ages of 2.4 and 14  Gyr, corresponding to the age range of the oldest mass bin in M12.

In Fig. 7 we show the comparison between stellar masses from the CL method and Chang et al. (2015) (left panel), and the comparison between formed masses from CL method and VESPA (right panel). For consistency, $sSNR_{L}$ models adopted a Chabrier (2003) IMF to compute present-day stellar masses and a Kroupa (2001) IMF to compute formed masses. Notwithstanding small biases, the agreement between our masses and those of Chang et al. (2015) is excellent, and more than sufficient to validate our photometric corrections.

Appendix D Data Set Comparison

In this appendix we compare the confidence contours calculated via the CL and SFHR methods for different datasets. This is shown in Fig. 8: the top panels use supernovae classified as SNIa, whereas the lower panels accept both SNIa and zSNIa classifications. The left panels use hosts as defined in the original papers, whereas the right panels adopt the hosts of S18. In order to apply the CL method, we impose redshift and photometric error limits to all datasets (see §2.4). We note that the contours change little if $z_{\rm{max}}=0.4$ had been chosen instead of $z_{\rm{max}}=0.2$.

The slopes $s$ in Fig. 8 are consistent, and therefore we focus on the scale factor. First, when both methods are applied to the same sample (red and green contours), the SFHR method leads to smaller scale factors (as also noted in §4). One possible explanation is that the VESPA masses tend to be higher than the estimated FSPS masses for the most massive galaxies in the sample. This trend is hinted at in the right panel of Fig. 7, and we note that the effects of binning masses in the SFHR method may suggest an even lower scale factor (see Fig. 3). Understanding this discrepancy is beyond the scope of this work, but, as argued in appendix C, we believe the scale factor derived from the CL method is reliable.

Second, we have applied the same method (CL) to control galaxies from M12’s and H17’s datasets (green and blue contours, respectively). If the hosts are retrieved from the original works, then the scale factor assessed for M12’s dataset is much lower (top left panel). This makes sense, given that both datasets contain approximately the same number of hosts (in the chosen redshift range), but M12’s sample contains almost twice as many galaxies as H17’s sample. This discrepancy is likely because there are some hosts missing in M12’s sample and is resolved when using hosts from S18’s list.