Think Global, Act Local: The Influence of Environment Age and Host Mass on Type Ia Supernova Light Curves

The photometry of the environment around the site of the supernovae are taken from the "Scene Modeling" method described in Holtzman et al. (2008). The method incorporates all the images taken during the survey to build a photometric model at the location of the transient. The resulting photometry at the site is a convolution of the galaxy on the scale of the typical seeing. By applying this fixed angular aperture, we always get the most compact local environment possible. Finally, in order to keep the environment truly "local," a redshift cutoff was imposed, z < 0.2. With the average seeing for SDSS being 14, the maximum extent of a galaxy's local environment was 3 kpc in radius. At higher redshifts, the angular resolution encompasses a majority of typical galaxy, and there is little difference between local and global properties. The typical size of the projected aperture defining the environment at the supernova location was 1.5 kpc.

Since core collapse supernovae are less luminous than SN Ia, their contamination percentage increases in the low-redshift portion of any data set. Anticipating this higher percentage of core collapse supernova (CC) interlopers, we added further Hubble residual cut to minimize the CC contamination. In addition, the statistics used in this work are robust against the ≳3.9% CC contamination estimated in C13. A detailed study of CC contamination affecting SN Ia standardization with host-galaxy properties should be done because the ratio of SN Ia and CC is highly dependent on host properties.

This results in a final data set of 103 SN Ia. A partial list of the final data set used in the local environment analysis is visible in Table 1.

In addition to the photometry of the local environment, we also study the correlation between the supernova characteristics and the host properties as a whole (hereafter the global properties). From the global host properties, we can compare our results directly to the work presented in Gupta et al. (2011), hereafter G11, to check if our age estimator is more sensitive to younger stellar populations. Second, we can contrast the local and global properties of hosts to check if there is more information contained in the local environment rather than the more easily studied global properties.

The final analysis of G11 included 206 SN Ia and hosts. We looked at the 76 objects that passed our redshift and other quality cuts to use in our method validation. A sample of the data set used for this validation is visible in Table 2.

For each of the 103 hosts where we had local environment photometry, we gathered the SDSS-DR12 model magnitudes for the estimate of the global properties. A sample of these data set can be seen in Table 3.

The SDSS model magnitudes are unreliable for galaxies with a large angular extent. For the nearby galaxies with distances calibrated with Cepheid variable stars, we used aperture photometry to obtain both the local and global magnitudes. Images of the large galaxies were downloaded from the SDSS-DR12 and individual apertures were designed to capture 90% of the combined light in all the filters. After masking out stars projected on the galaxy, the aperture was then applied to the image of each filter. The magnitude was then calibrated using nearby stars. The photometry for these galaxies can be seen in Table 4.

We use the C13 supernova sample to provide light-curve properties and Hubble residual information. We use the Malmquist bias corrected distances derived from the best-fit cosmology (H₀ = 73.8 km s⁻¹ Mpc⁻¹, Ω_M = 0.24, and Ω_Λ = 0.76). When needed, we use these values for our assumed cosmology. For the maximum redshift in our sample, the Malmquist bias correction is ∼0.01 mag. However, C13 noted that the stretch correction coefficient, α they found for their full sample was significantly larger than typical and larger than the α derived from their spectroscopically classified sub-sample. After our cuts, we found a significant correlation between the supernova stretch parameters, x₁, and the C13 Hubble residuals, which would likely result in spurious correlations with our host-galaxy analysis. We corrected the C13 Hubble residuals using their spectroscopically derived α value of 0.16 and no longer detected a significant correlation between x₁ and our sample's Hubble residuals. The resulting Hubble residuals can also be found in Table 1.

For the nearby SN Ia used in the Cepheid calibration of SN Ia peak luminosity, we obtained light curves from the SNANA database² and fit them using SALT2.4 implemented from the sncosmo.³ The model also corrected for Milky Way dust extinction from the dust maps of Schlegel et al. (1998) and Schlafly & Finkbeiner (2011) via sfdmap.⁴

Many of the FSPS parameters were set at their default values, but a number of key settings were adjusted to produce the desired model space.

To control the metallicity of our stellar population, we set zcontinuous = 2. This setting convolves the SSPs (simple stellar populations) with a metallicity distribution function. The metallicity distribution is defined as

$$\begin{eqnarray}&&{\left({{Ze}}^{-Z}\right)}^{{\mathtt{pmetals}}}\end{eqnarray} \tag{ 1 }$$

with

$$\begin{eqnarray}&&Z\equiv \displaystyle \frac{z}{{z}_{\odot }{10}^{\mathrm{log}(z/{z}_{\odot })}},\end{eqnarray} \tag{ 2 }$$

where z is the metallicity in linear units and z_⊙ = 0.019. This metallicity distribution is governed by two more FSPS parameters: pmetals and logzsol (i.e., $\mathrm{log}(z/{z}_{\odot })$). We left pmetals at its default value of 2, and during the fitting process logzsol was allowed to vary but was marginalized over when the age probability distribution was determined.

The next set of parameters govern the treatment of dust. We used the default power-law dust model as explained in Conroy et al. (2009), based off of Charlot & Fall (2000). The attenuation curve of a star, as a function of stellar age, is defined as

$$\begin{eqnarray}{\tau }_{\lambda }(t)=\left\{\begin{array}{ll}{\tau }_{1}{\left(\lambda /5500\mathring{\rm A} \right)}^{-0.7} & t\leqslant {10}^{7}\,\mathrm{yr}\\ {\tau }_{2}{\left(\lambda /5500\mathring{\rm A} \right)}^{-0.7} & t\gt {10}^{7}\,\mathrm{yr}\end{array}\right.,\end{eqnarray} \tag{ 3 }$$

where τ₁ and τ₂ are the attenuation around a young stellar population and in the ISM, respectively. See Charlot & Fall (2000), Figure 1, for a visual representation. In FSPS these two parameters are controlled via the dust1 and dust2 variables. For this analysis, dust1 was set to two times dust2, and dust2 was allowed to vary to match the observations. Conroy et al. (2009) claims good values of dust1 and dust2 are 1.0 and 0.3, respectively, while Charlot & Fall (2000) prefers values of dust1 and dust2 of 1.0 and 0.5, respectively. The allowed range for dust2 in this analysis is explained in Section 4.2 and is consistent with these recommendations.

A few of the host SEDs show an unusual r-band feature (e.g., SN4019). Adding nebular emission (setting add_neb_emission = True and cloudy_dust = True) adjusts the r-band magnitude for a young stellar population and allows the model to match this observed feature. This characteristic is shown to be achievable in the self-consistency validation test number 3, as explained in Section 5.1.

Finally, FSPS outputs the luminosity of 1 M_⊙, so an extra constant, δ, is used to scale the output SED of FSPS to match the observed SEDs.

FSPS has many inputs for describing the star formation history of a galaxy. The sfh parameter allows the user to select the functional form of star formation history. G11 used the simple τ-model: the star formation rate is proportional to e^−t/τ, with t being the time since the start of star formation and τ is a free parameter. G11 fit both τ and the length of star formation history. This is the simplest model, which is important when fitting a small number of free parameters. However, such a simple prescription makes it difficult to create a model with both old stars and recent star formation.

Simha et al. (2014) investigated the ability of several star formation history models to match simulated galaxies. This research looked at the simple τ-model, a linear-exponent model,⁶ and a four-parameter τ-model. A visual comparison is presented in Figures 3–5 of Simha et al. According to the calculations in Simha et al., the simple τ-model can overestimate the age by ∼1–2 Gyr, particularly for younger populations. Since we expect some supernovae to explode in young (≲2 Gyr) stellar populations, we decided to use a four-parameter τ-model.

The main feature of the four-parameter τ-model (4pτ-model) is that it separates the properties of early and late time star formation. This model can describe a wide range of star formation histories: an early burst, a history dominated by recent star formation, or both an early burst and recent star formation. The 4pτ-model is a piecewise combination of a linear-exponent star formation history, then a linearly rising or falling star formation. This model is used by FSPS when sfh = 5. Mathematically the 4pτ-model can be written as

$$\begin{eqnarray}{{\rm{\Psi }}}_{\star }(t)\propto \left\{\begin{array}{ll}(t-{t}_{0})\,{e}^{-(t-{t}_{0})/\tau } & {t}_{0}\leqslant t\leqslant {t}_{i}\\ {{\rm{\Psi }}}_{\star }({t}_{i})+{m}_{\mathrm{sf}}\,{\mathscr{R}}(t-{t}_{i}) & {t}_{{\rm{i}}}\lt t\leqslant {\mathscr{A}}(z)\\ 0 & \mathrm{else}\end{array}\right.,\end{eqnarray} \tag{ 4 }$$

where t₀ is when the star formation started, τ is the e-folding parameter, t_i is the star formation transition time, m_sf is the slope of the late time star formation, ${\mathscr{R}}$ is the ramp function, and ${\mathscr{A}}(z)$ is the redshift dependent time when the observed light was emitted. Note that the equation above allows negative values of Ψ, which is nonphysical. So we add an extra constraint that forces Ψ(t) to be 0 if calculated to be negative. The four variables in the equation (t₀, τ, t_i, m_sf) are the free parameters that give this model its name. These correspond to the FSPS parameters sf_start, tau, sf_trunc, and sf_slope, respectively. Also, this function takes t as the time from the start of the universe. A sample of various star formation histories calculated from the tau model can be seen in Figure 1.

Zoom In Zoom Out Reset image size

For any set of star formation parameters, we determine the average population age. The mass weighted average age is

$$\begin{eqnarray}&&\langle A{\rangle }_{\mathrm{mass}}={\mathscr{A}}(z)-{t}_{0}-\displaystyle \frac{{\int }_{{t}_{0}}^{{\mathscr{A}}(z)}(t-{t}_{0}){{\rm{\Psi }}}_{\star }(t){dt}}{{\int }_{{t}_{0}}^{{\mathscr{A}}(z)}{{\rm{\Psi }}}_{\star }(t){dt}},\end{eqnarray} \tag{ 5 }$$

where all of the variables are the same as described in Equation (4), so t − t₀ is simply the length of star formation. In the integral, a variable substitution of $t-{t}_{0}\to t$ transforms the time zero-point from the beginning of the universe to the start of star formation. We would also need to transform ${{\rm{\Psi }}}_{\star }(t)\to {{\rm{\Psi }}}_{\star }^{{\prime} }(t)$ so that ${{\rm{\Psi }}}_{\star }^{{\prime} }$ assumes t = 0 is the start of star formation. This makes the equation for the mass weighted average age to be

$$\begin{eqnarray*}&&\langle A{\rangle }_{\mathrm{mass}}={\mathscr{A}}(z)-{t}_{0}-\displaystyle \frac{{\int }_{0}^{{\mathscr{A}}(z)-{t}_{0}}t{{\rm{\Psi }}}_{\star }^{{\prime} }(t){dt}}{{\int }_{0}^{{\mathscr{A}}(z)-{t}_{0}}{{\rm{\Psi }}}_{\star }^{{\prime} }(t){dt}}\end{eqnarray*}$$

or the equation used in G11. In this paper, the age of a stellar population will refer to the mass weighted average described here.

The model using the star formation parameters t₀ = 8.0, τ = 0.1, t_i = 12, and m_sf = 20 produces a population with an average age of 437 Myr. This demonstrates that our SFH prescription can generate a dominant young population when SN Ia are expected to start exploding. A small τ is needed to keep the number of old stars from building up over the cosmic time and dominating over the recent linear star formation period. An example of this can be seen by the "old and young burst" star formation history in Figure 1.

This method uses a standard log-likelihood function for data with Gaussian uncertainties. Summing over each filter, we compare the observed apparent magnitude (m_i) with the resulting apparent magnitude from FSPS (m_FSPS,i), plus a scaling factor (δ) to account for FSPS's 1 M_⊙ output. Mathematically this is written as

$$\begin{eqnarray}&&\mathrm{ln}({\mathscr{L}})\propto \displaystyle \sum _{i}\displaystyle \frac{{\left({m}_{i}-\left({m}_{\mathrm{FSPS},i}+\delta \right)\right)}^{2}}{{\sigma }_{i}^{2}}+\mathrm{ln}\left(2\pi {\sigma }_{i}^{2}\right),\end{eqnarray} \tag{ 8 }$$

with σ_i as the uncertainty in each m_i measurement.

For five of the variables, we use bounded flat tops:

$$\begin{eqnarray}\begin{array}{rcl}2.5 & \lt & {t}_{i}\lt {\mathscr{A}}(z)\\ 0.1 & \lt & \tau \lt 10\\ -1.520838 & \lt & \phi \lt 1.520838\\ 0.5 & \lt & {t}_{0}\lt {t}_{i}-2.0\\ -45.0 & \lt & \delta \lt -5.0.\end{array}\end{eqnarray} \tag{ 9 }$$

Since a flat distribution in a slope parameter preferentially searches the high values space,⁷ the MCMC was performed over ϕ, the angle the ramp function makes with respect to the x-axis; therefore, ${m}_{\mathrm{sf}}=\arctan (\phi )$. The prior bound above, −1.520838 < ϕ < 1.520838, corresponds to −20 < m_sf < 20.

In addition, $\mathrm{log}(z/{z}_{\odot })$ uses a Gaussian distribution with μ = −0.5 dex and σ = 0.5 dex limited to the range of $-2.5\lt \mathrm{log}(z/{z}_{\odot })\lt 0.5$. This is a common assumption as seen in Belczynski et al. (2016). Our model reaches a lower metallicity than the grid search used in G11.

For the ISM dust parameter, τ₂, we assume a Gaussian prior on the top bounds. The Gaussian distribution is defined by μ = 0.3 and σ = 0.2, but only values between 0 and 0.9 are accepted. This allows for some variability but keeps the values near the 0.3 and 0.5 as recommended by Conroy et al. (2009) and Charlot & Fall (2000), respectively.

The first validation of this newly developed age estimator was to verify that it was self-consistent (i.e., it could correctly estimate the star formation parameters from an SED generated by FSPS).

Eight different models were used in this test. The initial metallicity, dust, and star formation parameters can be seen in Table 5. These values were put into FSPS, at a redshift of z = 0.05, to generate observed SEDs. The resulting SEDs, Table 6, were then analyzed with our Bayesian estimator. This set of models produced old populations (∼10.5 Gyr) and very young populations (∼0.5 Gyr). They explored the effect of metallicity (Model 5) and dust (Model 6). Model 2 also looked at a "mixed" population with an old burst of star formation and a strong increase of star formation to the present epoch, a stellar population that is not possible with a simpler star formation history.

Table 5.  SFH Parameters Used in the Self-consistency Test

ID	$\mathrm{log}(z/{z}_{\odot })$	τ2	τ	t0	ti	ϕ	Age
			(Gyr−1)	(Gyr)	(Gyr)	(rad)	(Gyr)
1	−0.5	0.1	0.5	1.5	9.0	−0.785	10.68
2	−0.5	0.1	0.5	1.5	9.0	1.504	1.41
3	−0.5	0.1	7.0	3.0	10.0	1.504	1.75
4	−0.5	0.1	7.0	3.0	13.0	0.0	4.28
5	−1.5	0.1	0.5	1.5	9.0	−0.785	10.68
6	−0.5	0.8	7.0	3.0	10.0	1.504	1.75
7	−0.5	0.1	0.5	1.5	6.0	1.504	2.40
8	−0.5	0.1	0.1	8.0	12.0	1.52	0.437

Note. All models are at a z = 0.05.

Download table as:  ASCIITypeset image

Table 6.  SEDs for the Self-consistency Test

ID	u	g	r	i	z	$\langle A\rangle$
	(mag)	(mag)	(mag)	(mag)	(mag)	(Gyr)
1	20.36	18.76	17.99	17.67	17.39	8.5 ± 1.5
2	20.31	18.74	17.98	17.66	17.39	7.7 ± 1.5
3	16.15	15.43	15.40	15.19	15.21	1.4 ± 0.5
4	17.65	16.74	16.49	16.26	16.16	4.2 ± 1.0
5	19.69	18.29	17.70	17.45	17.29	7.2 ± 1.6
6	17.66	16.58	16.25	16.01	15.86	2.7 ± 0.8
7	17.62	16.80	16.57	16.34	16.26	4.5 ± 1.4
8	19.72	18.37	17.88	17.68	17.56	4.4 ± 1.1

Note. SEDs were scaled with a δ = −25 mag.

Download table as:  ASCIITypeset image

This method can model a young stellar population, (like Model 3), but is unable to recover a starburst or mixed populations (Models 8 and 2, respectively) based on SED fitting. The "old and young burst" of Model 2 is not particularly blue, and as such, we identify an age of 4.4 Gyr. We do better with Model 3, where we estimate the correct age of 1.4 Gyr with a small uncertainty.

Figure 2 shows the FSPS output SEDs of the initial input star formation parameters and the best-fit parameters from our analysis. The SEDs are fit very well across all the models used in this self-consistency test. In addition to fitting the SEDs, our analysis was able to approximate the underlying parameters. An example corner plot of the posterior probabilities is displayed in Figure 3, and the the full set of figures is available online.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

View figure set (8 panels)

Tarball of figure set images

The final validation was to recalculate the global host ages originally presented in G11. A direct comparison between ages from G11 and the results from our analysis can be seen in Figure 4. Most points fit along the one-to-one line with a scatter of around 2 Gyr, implying that our analysis is consistent with the previous work. However, there are six hosts that our method estimates to be ≲2 Gyr, whereas none of the ages in G11 were that young. At ≲4 Gyr, our method systematically estimates a younger age. This is because our star formation history model allows for more recent star formation than the simple τ models permits. In their discussion on this topic, Simha et al. (2014) claimed that the star formation model used in G11 can overestimate young populations by ∼2 Gyr. This overestimation can be seen in Figure 4.

Zoom In Zoom Out Reset image size

First, we compare the Hubble residuals with ages derived from the global photometry of the hosts using C13 (Table 3) sample. The comparison is presented in Figure 5. This data set has a low 2.1σ correlation, as defined by the Spearman's value of −0.23. In addition, the distribution in Hubble residual-age space appears to show a distinctive "step" between 7 and 8 Gyr.

Zoom In Zoom Out Reset image size

Finally, we compare the Hubble residuals versus average local environment age for the data set derived from C13 (Table 1). These results are presented in Figure 6. Using an age derived from the local environment only slightly changes the Spearman's correlation between these two parameters, But this correlation, −0.21, only has a 1.8σ significance. The overall age distribution and apparent "step" at ∼8 Gyr are not significantly changed by switching to a local environment analysis. At first glance, the local age does not appear to contain any additional information not already present in the global age. Jones et al. (2018) found a marginally significant difference between their global and local analyses.

Zoom In Zoom Out Reset image size

Because a fraction of our sample is photometrically classified, there may be some CC contamination that would be found preferentially in the upper left of Figure 5 and 6. This contamination might contribute to the observed correlations, but at a level that is small compared with the 2.1σ and 1.8σ trends.

The global age estimated for spiral galaxies is an average of several stellar populations. The prompt component of SN Ia are expected to be found in the youngest regions of a galaxy. So we may see that the population age at the supernova location would tend to be younger than the global age of the galaxy. In Figure 7 we show the difference between the estimated global and local ages. For massive galaxies with local ages less than 4 Gyr, there is a tendency for the supernova to be located in a younger than average spot in the host. The effect is less apparent for low-mass galaxies, but that is likely because the size of the region measured for the local environment is a large fraction of the size of a small galaxy.

Zoom In Zoom Out Reset image size

Between local ages of 4 and 8 Gyr the difference between global and local age estimates show a large scatter with no apparent trend. Beyond 8 Gyr the scatter between the local and global ages is reduced, probably because the stellar population in ellipticals is fairly uniform. In these older hosts the global age estimates tend to be 1–2 Gyr younger than the local ages. The reason for this difference is not clear, but it may be due to activity at the center of ellipticals contaminating the stellar colors.

Although we see no statistical difference between the use of local and global ages when comparing with Hubble residuals, our results suggest that using local photometry to characterize the supernova environment has some benefit over the global average. For example, the SED measured local to the supernova can indicate a population that is a factor of two younger than the global average in large star-forming galaxies. When feasible, the measurement of the local environment, particularly of younger populations, provides a more accurate representation of the progenitor age than simple averaging the light from the host.

The Hubble residuals in the C13 sample have a monotonically decreasing trend (at 2.1σ) that appears to be a break or step at an age of ∼8 Gyr. Figure 8 plots the same data as Figure 6, but this time splits the data into two age bins: ≤8 Gyr and >8 Gyr. Both age ranges show very small Spearman's correlations (−0.07 and −0.09, respectively) that are consistent with a flat distribution. The correlation decreasing in significance when the data set is split in two implies that the monotonic function seen in Figure 6 is really a step-like function with a transition at ∼8 Gyr. This function may have a transition width making it more like a continuous sigmoid function with a transition faster than our age resolution. In both the local and global age measurements, we find a significant age step in Hubble residuals.

Zoom In Zoom Out Reset image size

The amplitude of the jump between these two populations split at ∼8 Gyr appears to be more than 0.1 mag. The younger population has a mean Hubble residual just over zero (0.005 ± 0.018 mag), whereas the older population has a mean Hubble residual of −0.109 ± 0.035 mag. This makes the resulting difference in the means equal to 0.114 ± 0.039 mag (2.9σ). This is almost two times larger than the commonly used mass step of 0.06 mag. A step of this size may affect precision cosmological measurements since the fraction of each subpopulation is expected to change with redshift. For example, one would not expect to find stellar populations as old a 8 Gyr at redshifts greater than ∼1, given the current standard cosmology. Such an age step may also impact local measurements of the Hubble constant, as the peak luminosity of SN Ia has so far been calibrated with Cepheid variables found exclusively in star-forming galaxies. With the age step being detected at a 3σ significance, Figure 8 suggests a need for an additional SN Ia luminosity correction based on their host's local or global stellar age.

Childress et al. (2014) argues for a link between host stellar mass and the delay time between stellar formation and SN Ia explosion, with results summarized in their Figure 4 where the SN Ia progenitor age distribution is divided into host-galaxy stellar mass bins. They find, with reasonable assumptions of star formation histories and SN Ia delay times, that there is a natural division in host mass and age between prompt and "tardy" SN Ia. Prompt SN Ia occur in lower mass galaxies with ages ≲6 Gyr, while tardy events continue in high mass galaxies with ages ≳6 Gyr. Projecting their model onto the host mass axis results in an overlapping distribution of prompt/tardy explosions with a transition near 10^10.5 M_⊙. This transition in progenitor age may correspond to the mass step observed at 10¹⁰ M_⊙. Projection of the Childress et al. (2014) model onto the age axis provides a clean separation between the prompt and tardy SN Ia with a dearth of events between 4 and 8 Gyr. Our age estimates are not consistent with a deficit of supernovae exploding in that range, but we do see a shift in the SN Ia light curve properties around a stellar age of about 8 Gyr.

We have estimated the host-galaxy masses in our sample to test if our measurements agree with the prediction of Figure 4 in Childress et al. (2014). We use the kcorrect code (v4_3)⁸ described in Blanton & Roweis (2007). This code utilizes spectral fitting templates based on the stellar population synthesis models of Bruzual & Charlot (2003), which are calculated using the Chabrier (2003) initial mass function. We input the SDSS model magnitudes and the pipeline redshifts for each galaxy (presented in Table 3) to calculate the k-corrections and stellar mass-to-light ratios. The stellar masses are output in units of M_⊙h⁻¹, and we convert them to units of M_⊙ using the C13 cosmology described in Section 2.4. The uncertainties on the stellar masses are approximately ±0.3 dex. The resulting estimated stellar masses are reported in Table 7.

Table 7.  Results of Median Age, Stellar Mass, and PC₁

SN	Local Age	σa	Global Age	σa	$\mathrm{log}(M/{M}_{\odot })$	PC1
762	5.1	2.7	8.0	3.9	11.1	0.07
1032	5.8	2.6	8.6	4.3	10.5	−1.78
1371	8.9	4.3	8.1	1.3	10.7	−0.82
1794	4.0	1.4	2.0	0.9	8.8	2.20
2372	5.9	1.5	4.8	2.9	10.2	−0.07

Note. Stellar masses have a 0.3 dex uncertainty, and PC₁ was calculated with local environment ages.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

Figure 9 shows the distribution of the local SN Ia age versus stellar mass for our sample. For our sample, the distribution of hosts in age-stellar mass space shows similarities and differences with that predicted in Childress et al. (2014). The observations do show young hosts extending to low stellar masses reproducing the backward "L" seen in the Childress et al. (2014) simulation. This indication of "cosmic downsizing" is not as pronounced in our data, as we are probably still overestimating the ages of the extremely young populations. Most notably missing is that the predicted bimodal age feature is not present. An island of young galaxies is expected between 0.5–1.0 Gyr and that is not seen in Figure 9. Our age distribution is relatively flat from 2 to 5 Gyr and then declines out to 11 Gyr. This does not match the predicted distribution in Figure 4 of Childress et al. (2014), and there is no clear peak of old hosts. The fact that our sample of supernovae extends to a redshift of 0.2 may contribute to the lack of a clear peak in old age. Childress et al. does predict that the old-age peak decreases in height and age, as the data are taken at higher redshifts simply due to the finite age of the universe.

Zoom In Zoom Out Reset image size

Light curve shape and color of SN Ia have been shown to be strongly correlated with peak luminosity. Other modest trends with population age, host mass, and gas metallicity have been reported, and here we have identified a rather significant jump in Hubble residuals with local and global stellar age. Because these environmental observables are highly correlated with each other, it is interesting to ask if there is a linear combination of observables that have a clear and significant correlation with SN Ia Hubble residual. As an initial test, we performed a principal component analysis (PCA) on our data set.

PCA is a linear algebra tool that transforms the basis of a matrix of data to orthogonal axes where the new axes are maximally aligned and sorted with the intrinsic scatter of the data.⁹ Therefore the first principal component contains the most amount of information and the last principal component contains the least. Thus, it is possible to reduce the dimensions of a problem by retaining only the principal components that comprise most of the variance (or relative information) in the original data set. In addition, the first principal component will identify a linear combination of the original variables that accounts for most of the variance in the data. Interpretation of PCA is difficult because the results are sensitive to noise, specific parameterizations, normalization, and other implementation details. However, PCA can still be useful to understand how multiple parameters work together within a data set. In the search for a trend in Hubble residuals versus SN Ia parameters, we applied PCA on the parameters of SALT2 stretch (x₁) and color (c), as well as host mass and age. We calculated PCA coefficients first using local ages and then with global ages. Hubble residuals were not included as a PCA parameter. Only after the PCA did we search for correlations between the resulting principal components and the Hubble residuals.

Before running PCA, we normalized all parameters by removing the mean and scaling to unit variance resulting in normalized input parameters (${x}_{1}^{{\prime} }$, c', m', and a'). An example, for x₁, is defined as

$$\begin{eqnarray}&&{x}_{1}^{{\prime} }=({x}_{1}-{\mu }_{{x}_{1}})/{\sigma }_{{x}_{1}}.\end{eqnarray} \tag{ 10 }$$

The means and standard deviations used in the normalization process can be seen in Table 8.

Table 8.  Normalization Parameters Applied before PCA

	μ	σ
x1	−0.177	1.015
c	0.0100	0.0829
$\mathrm{log}(M/{M}_{\odot })$	10.15 dex	0.69 dex
age	5.22 Gyr	2.11 Gyr

Download table as:  ASCIITypeset image

From these four observables, PCA yields four principal components (PC_i). Table 9 shows the linear combination of the observable variables that make up each principal component, as well as the explained variance. Substituting global age for local age resulted in some minor differences between the PCA coefficients and variance, but the overall conclusions are very similar between the two analyses. PC₂ accounts for a quarter of the variance and is dominated by the SALT2 color parameter. The other two PCA components contribute only a small portion of the variance.

Table 9.  PCA Coefficients Using Local and Global Ages

	x1	c	$\mathrm{log}(M/{M}_{\odot })$	Agea	% Variance	x1	c	$\mathrm{log}(M/{M}_{\odot })$	ageb	% Variance
PC1	0.557	−0.103	−0.535	−0.627	44	0.465	−0.134	−0.596	−0.641	47
PC2	−0.159	0.956	−0.213	−0.118	25	−0.206	0.937	−0.263	−0.101	25
PC3	−0.651	−0.258	−0.710	0.071	18	−0.846	−0.315	−0.385	−0.190	19
PC4	0.490	0.086	−0.404	0.767	11	−0.158	0.070	0.654	−0.737	8

Note. All observables are normalized via Equation (10).

^aThe median age of the local environment posterior. ^bThe median global age posterior.

Download table as:  ASCIITypeset image

As a first analysis, it is important to see if there are any correlations between Hubble residual and these SN Ia-host-galaxy principal components. Looking at the PCA done with local age, PC₃ and PC₄ versus Hubble residual have extremely low Spearman's coefficients of 0.21 and −0.15, respectively. PC₂, dominated by the SALT2 color term, when compared with Hubble residual is also a scatter plot with a low Spearman's correlation of −0.06. Surprisingly there is an increase in the scatter at the high PC₂ (high color) domain. Using the PCA done with the global age, PC₂, PC₃, and PC₄ have similar Spearman's correlation coefficients. Figures showing the relationship between the Hubble residual and each principal component (four per analysis method, eight total figures) are available as a figure set in the online journal.

Figure 10 shows the relationship between Hubble residuals and PC₁ using the estimated local age. The correlation when PCA is applied using the global host age can be seen in Figure 11(a). PC₁ with either local or global ages shows a very strong correlations with Hubble residuals. The first principle component, using local age estimates, is defined as

$$\begin{eqnarray}&&{\mathrm{PC}}_{1}=0.56{x}_{1}^{{\prime} }-0.10c^{\prime} -0.54m^{\prime} -0.63a^{\prime} ,\end{eqnarray} \tag{ 11 }$$

or it can be approximated by ignoring the color term because it barely contributes to PC₁.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

View figure set (8 panels)

Tarball of figure set images

Zoom In Zoom Out Reset image size

Uncertainties for these coefficients can be obtained via bootstrap resampling (Wall & Jenkins 2012, Section 6.6). Here we run a PCA on 100,000 data sets that were created randomly (with replacement) from our original data set using local ages. Since principal components are invariant to being multiplied by −1, a constraint was made that the bootstrap eigenvector needed to have a positive dot product with the original eigenvector. If this was not true, the bootstrap eigenvector's direction was reversed. The resulting distribution of coefficients for PC₁ can be seen in Figure 12 and a statistical summary can be seen in Table 10. As expected, the color coefficient is consistent with zero. Interestingly, the coefficient for local age is more constrained than for stellar mass.

Zoom In Zoom Out Reset image size

In PC₁, the SALT2 color coefficient is very small, implying supernova color is not a strong contributor to this Hubble residual correlation. The mass, age, and stretch parameters in our PCA analysis are similar in amplitude and likely are of similar importance in any further standardization of SN Ia distances. The SALT2 stretch parameter, x₁, has a surprising large contribution to PC₁, given that the SN Ia have already been corrected for light-curve shape. Figure 10 suggests a correlation between stretch and Hubble residual, but this correlation was removed at the start of the analysis. Instead, PC₁ shows that the value of x₁ is related to mass and age. Attempting to correct for the stretch-luminosity relationship requires including host properties, as all three have an influence on Hubble residuals. Our results suggest that the α parameter derived from the SALT2 fit is not ideal, because the effects of host mass and population age are not distributed uniformly with stretch.

The Hubble residual-PC₁ correlation in Figure 10 has a Spearman's correlation of 0.44. This trend is highly significant, with a p-value of 3.86 × 10⁻⁶ corresponding to a 4.7σ significance. The trend in Hubble residual versus age has a significance of only 2.1σ (Figure 5), while including mass and stretch greatly increases the significance. This strong correlation suggests that host properties influence the luminosities of SN Ia beyond the currently applied corrections.

This data set did not have any significant correlations between Hubble residual and x₁ or c. But when x₁ is combined with stellar mass and age, there is a significant correlation with Hubble residual.

When using global ages, the underling trend has a Spearman's correlation of 0.40 or a 4.0σ significance. The Spearman's correlation with PC₁ from local ages is slightly more significant than for global ages, suggesting that measurements near the event provide some improvement in correlating supernova with environment.

The best-fit linear regression of this trend is

$$\begin{eqnarray}&&\mathrm{HR}=0.051\,\mathrm{mag}\times {\mathrm{PC}}_{1}-0.012\,\mathrm{mag},\end{eqnarray} \tag{ 12 }$$

with uncertainties in the slope and intercept as ± 0.011 mag and ± 0.016 mag, respectively. This trend reduces the 1σ scatter in Hubble residual from 0.17 to 0.15 mag. The best-fit linear regression using global ages values is HR = 0.045 × PC₁ −0.012, with uncertainties in the slope and intercept as ±0.011 mag and ±0.016 mag, respectively. This reduces the 1σ scatter in the Hubble residual to 0.16 mag. As a note, the PCA does not maximize this correlation or minimize the χ² parameter in a linear regression. For example, reversing the sign of the x₁ coefficient nearly doubles the slope of the correlation with Hubble residual as seen in Figure 11(b) and significantly reduces the χ² parameter of a linear fit. Understanding the best use of these coefficients will be part of future work already in preparation.

Since the age and mass coefficients make up nearly two-thirds of the correlation amplitude, we see that a change of ±2σ in the normalized mass and age values results in a 0.24 mag shift in SN Ia brightness. Cosmic downsizing suggests these two parameters are correlated since young galaxies tend to be small while typical old hosts are massive in the current epoch. Low-mass, young galaxies tend to be metal deficient, while old, massive galaxies can be metal-rich. Thus, this combination of age and mass may indicate progenitor metallicity influences the peak luminosities of SN Ia.

The correlation between Hubble residual and PC₁ is very significant, and it is unlikely that PCA could generate such a strong correlation from a random distribution. To test the probability that this correlation is caused by chance, we applied a bootstrap style method for hypothesis testing. We took our PCA input matrix (x₁, c, host mass, local environment age) and shuffled the order along each parameter, creating 103 "new" SN Ia. There was no cross shuffling, so a stretch always stays a stretch; it just corresponds to a different SN Ia. We applied PCA on each shuffled sample and tested for any correlations with Hubble residual. After 32,000 runs, a 4σ test, the maximum Spearman's correlation coefficient between Hubble residual and any principal component was 0.43, while the measured Spearman's test of the non-shuffled values was 0.44. The distribution of the bootstrapped Spearman's correlation coefficients is shown in Figure 13. This bootstrap style analysis shows that the false-positive Spearman's correlation follows a Gaussian distribution with a σ ≈ 0.1. This is the expected distribution for Spearman's correlation coefficients for a data set of N ≈ 100. It is exceedingly unlikely that the correlation between Hubble residual and PC₁ found in Figure 10 would appear at random. As a result, this luminosity-stretch-mass-age relationship is the most significant systematic seen between calibrated SN Ia and a host-galaxy environment.

Zoom In Zoom Out Reset image size

Since our PC₁ strongly correlates with Hubble residuals, it is reasonable to consider a modification to the Tripp formula that is used to correct SN Ia peak luminosities. The new correction coefficients would be the multiplication of the PCA coefficients and the slope of the Hubble residual-PC₁ trend. Except for the SALT color, the individual components making up PC₁ have significant weights and thus are included in this modified equation. The stretch parameter in PC₁ can be grouped with the original Tripp coefficient, leaving a new term with just host properties. Since PC₁ has a positive correlation with Hubble residual, the coefficient should come in with a negative sign. From these results, and following the example of others (e.g., Moreno-Raya et al. 2016b), we propose a change in the distance modulus corrections performed by SALT2 by modifying the equation to include PC₁. This new equation would be

$$\begin{eqnarray}&&\mu ={m}_{B}-{M}_{B}+(\alpha -\alpha ^{\prime} ){x}_{1}-\beta c+\gamma \,m^{\prime} +\gamma \,a^{\prime} ,\end{eqnarray} \tag{ 13 }$$

where α' ≈ γ ≈ 0.03 mag. This uses both the approximate form of PC₁ and encompasses the very slight renormalization of x₁ into ${x}_{1}^{{\prime} }$. More research is needed to accurately determine α' and γ.

This modified distance correction formula now includes a host-galaxy stellar mass correction and a host age term. Childress et al. (2014) attempted to explain the mass step as an age dependency, and we see that both appear to have an impact on SN Ia luminosities. Our result is bolstered by Jones et al. (2018), who showed that stellar population color effects are still present after stellar mass corrections. Stellar mass plus population age is similar to the Mannucci relationship (Mannucci et al. 2010) that connects galaxy mass and star formation rate with gas phase oxygen abundance. Our combined age plus mass parameters work in the same sense as the Mannucci relation and may be a proxy for metallicity. Both Hayden et al. (2013) and Moreno-Raya et al. (2016b) have indicated that metallicity is the underlying galactic variable influencing Hubble residuals, and our results support this suggestion.

Kasen et al. (2009) shows that varying the metallicity of SN Ia progenitors will shift the slope of the Phillips relationship. To investigate this, we assume our SN Ia have a ${M}_{\mathrm{peak}}=-19$. We then add in the Phillips relationship (with α = 0.16) and our Hubble residual-PC₁ relationship. The result is seen in Figure 14. As fast decliners are preferentially in larger and older hosts, we see the host-galaxy's effect on the Phillips relationship is stretch dependent. Put another way, the slope of the Phillips relationship, α, is dependent on the host-galaxy properties of the sample. These results show that these parameters are interrelated and that there is likely a multi-dimensional relationship between peak absolute magnitude, decline rate, and progenitor metallicity. Research into this multi-dimensional relationship will be a part of a future work already in preparation.

Zoom In Zoom Out Reset image size