THE CARNEGIE SUPERNOVA PROJECT: INTRINSIC COLORS OF TYPE Ia SUPERNOVAE

The left-hand graph of Figure 5 shows the time of first maximum of the NIR light-curves as a function of the decline-rate Δm₁₅(B). In general, slowly declining SNe Ia peak in the NIR 2–6 days earlier than B band, which presents an observational challenge to obtaining NIR coverage of the peak for objects whose observations are typically triggered based on optical band discoveries. However, it has been known for some time that the faster evolving events tend to peak later (after B-band maximum). Examining the right-hand graph in Figure 5, it would appear that we are seeing two groups of objects: a set with early NIR peak and a set with late NIR peak, which has been seen before (Krisciunas et al. 2009; Phillips 2012). The early NIR peak objects span the "normal" range 0.7 < Δm₁₅(B) < 1.7, while the late NIR peak objects are exclusively in the fast-declining region Δm₁₅(B) > 1.7. Furthermore, in the fast-declining region, there is a large spread of peak times with no correlation with Δm₁₅(B). On the other hand, when examining the right-hand graph of Figure 5, where we have used the s_BV parameter, the time of NIR peak is clearly correlated with s_BV. As a result, it can be argued that the bifurcation of the fast-declining sample into early- and late-risers is an artifact of the behavior of the light-curve parameter Δm₁₅(B).

Zoom In Zoom Out Reset image size

A second conspicuous feature of the NIR light-curves for normal SNe Ia is the presence of a secondary maximum. This is due to recombination of iron group elements in the supernova ejecta (Kasen 2006) and begins to be seen in r band as an inflection point, and develops into a secondary maximum in i band. In the z band and Y band, the secondary maximum can in fact be stronger than the primary (Stritzinger et al. 2002; Burns et al. 2010). It becomes less prominent in the J band, and then regains strength in the H band. For the fast-declining SNe Ia, however, both optical and NIR light-curves show either a very weak secondary or none at all. Following Krisciunas et al. (2001), we consider the strength of the i band secondary maximum as a function of decline rate. They define 〈f_λ(i)〉_{20 − 40} as the average flux (normalized to maximum) in the i-band from day 20 to day 40 after B maximum in the rest-frame of the SN Ia. This gives a measure of the prominence of the secondary peak. In the left panel of Figure 6 we plot this as a function of Δm₁₅(B). We see a clear correlation, albeit with significant scatter (see also Figure 8 of Folatelli et al. 2010, hereafter F10), and the fast decliners seem to separate into two groups. However, if we instead plot against s_BV (see right panel of Figure 6), there is again a more continuous transition between the slow and fast decliners.

Zoom In Zoom Out Reset image size

The transition (if one exists) from the slow to fast declining SNe Ia has been a major stumbling block for light-curve fitters. Clearly, a simple stretch in time cannot account for the NIR light-curves losing the secondary peak as we go from slow decliners to fast decliners. To our knowledge, SNooPy is the only light-curve fitter that currently attempts to fit NIR light-curves for both slow and fast decliners. Up until recently, SNooPy's NIR templates for the fast decliners were unreliable due to two problems: (1) small numbers of fast-declining SNe Ia with well-observed light-curves, and (2) the apparent failure of the Δm₁₅(B) parameter to distinguish between the fastest decliners.

Figure 7 shows the sorting of the 11 fastest decliners in the SNooPy training sample into bins of Δm₁₅(B). One can clearly see the problem: while the B-band light-curves seem to be sorted out correctly (with the possible exception of SN 2006mr), the Y-band light-curves clearly do not follow any kind of continuous progression from Δm₁₅(B) ≃ 1.5 to Δm₁₅(B) ≃ 2.0. The same holds true for the other NIR bands. The situation, however, is completely different when one uses the s_BV parameter instead. In Figure 8, we can see that SN 2006mr, which is clearly the fastest declining of the SNe Ia, is now the lowest color stretch. SN 2007ba and SN 2007on are now categorized as slower events. This has the very noticeable effect that the Y-band light-curves now show a more continuous progression from s_BV ≃ 0.2 to s_BV ≃ 0.6. Another notable feature is that SN 2005ke now looks very much like a transition object, having no secondary peak, but clearly an inflection point.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To further investigate this, we examined the spectra of four representative objects in this range of s_BV: SN 2006ax (s_BV = 0.985), SN 2007on (s_BV = 0.574), SN 2005ke (s_BV = 0.419), and SN 2006mr (s_BV = 0.260). Figure 9 shows the spectra of these four objects near maximum light. SN 2005ke shows several spectroscopic features common to the "Cool" SNe Ia: prominent Ti ii and Si ii λ5972. On the other hand, it shows an intermediate strength in the S ii "W" feature. The presence of Ti ii is highly sensitive to effective temperature and will therefore turn on rather abruptly, making any transition with s_BV equally abrupt. Another feature known to be correlated with light-curve width is Si ii λ5972 (Folatelli et al. 2013). In Figure 10 we plot the pseudo-equivalent width of Si ii λ5972 from Folatelli et al. (2013) as a function of s_BV. Interestingly, there is a more complicated relation than the linear trend seen with Δm₁₅. In particular, there seems to be very little trend for s_BV < 0.5.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

While it remains unclear whether s_BV is in fact a "better" parameter than Δm₁₅ in predicting the physics of the SNe Ia, the fact remains that it does produce a much more continuous family of light-curve templates for the faster-declining objects. For this reason, the SN Ia template-fitting package SNooPy will incorporate templates parameterized by s_BV and we will adopt it as our shape parameter for the remainder of the paper. It is also unclear whether s_BV will present a superior alternative when analyzing high-redshift SNe Ia for the purposes of cosmology, as the fast-declining events are less numerous and significantly fainter than normal SNe Ia.

The first thing we determined was whether or not the intrinsic color-stretch relations, $P_{j}^{N}(s_{{\rm BV}})$, require a linear or quadratic relation. To accomplish this, we introduce a selection function S(p₁, p₂) into our likelihood which selects a linear color relation, (N = 1) with probability p₁ and the quadratic color relation (N = 2), with probability p₂. In this way, the Markov Chain Monte Carlo (MCMC) chains will spend some time trying the linear model and some time trying the quadratic. We can compute the posterior probability P(N) by counting nodes in the chain. From that we can compute the odds ratio:

$$\begin{equation} O=\frac{P(N=2)}{P(N=1)}\frac{p_{1}}{p_{2}}. \end{equation} \tag{ 8 }$$

The odds ratio gives a measure of the likelihood of one model versus the other, giving us the ability to select one if the data so warrant. Note that the probabilities P in Equation (8) are computed by marginalizing over all parameters and so there is a built-in "Occam factor" that penalizes the N = 2 model as it has one more free parameter (Gregory 2011). Several pilot runs of our MCMC gave odds factors between 16 and 88 in favor of the quadratic model. The lower values of the odds ratio occur when we do not use the LRS prior and restrict the data to s_BV > 0.5. This can easily be seen in Figure 13 where fits to the data without the "low reddening" prior tend to have less curvature than when the LRS are constrained to have zero reddening. In either case, an odds ratio of greater than 10 is considered reasonably strong evidence (Jeffreys 1961) to favor the quadratic model over the linear one and we continue throughout with quadratic solutions for the SN Ia intrinsic colors.

When solving for the intrinsic colors, the factor that has the largest effect is whether or not we use the LRS prior. In most cases, $P_{j}^{N}$ ends up having a higher curvature when we use the LRS. This is most notable in the optical minus NIR colors V − Y, V − J, and V − H, where there are a significant number of objects at high color stretch that are bluer than the LRS and therefore tend to straighten out the intrinsic color curves. This is easily understood as the LRS prior, being associated with specific objects, implicitly defines a color-stretch dependence on the intrinsic colors, whereas the Cauchy and exponential priors apply to the sample as a whole and are "color-stretch blind." At the high- and low-color-stretch ends of the distribution, there are relatively few objects and so the solution for the Cauchy and exponential prior has more freedom at those ends.

The resulting E(B − V)_host color excesses also differ the most when comparing the LRS prior to the analytic priors. The left panel of Figure 15 shows the difference between the color excesses from the exponential and Cauchy priors. There is a systematic offset of about 0.02 mag in the sense that the exponential prior produces higher extinctions and there is a systematic difference with s_BV, though both effects are less than the 1-σ errors in the extinctions themselves. The right panel of Figure 15 shows the difference between the color excesses from the Cauchy and LRS priors. One can see there is a much larger offset in the sense that the LRS prior produces lower color excesses and there is a more pronounced systematic with respect to s_BV. This is also easily understood as Figure 13 shows that the intrinsic colors from the LRS are significantly redder than those from the Cauchy or exponential priors. In both cases, the systematic trend in the differences with s_BV is due to low- and high-s_BV ends having fewer objects so that the choice of prior has more of an effect.

Zoom In Zoom Out Reset image size

Along with varying the priors, we also varied the sample used to fit the intrinsic colors. In particular, we were interested in whether the fast-declining objects could be fit together with the normal Ia's. Figure 13 shows the best-fit intrinsic colors obtained when excluding (solid line) and including (dashed line) the fast-declining SNe Ia. The fits are remarkably close. Indeed, when using the color-stretch parameter s_BV instead of Δm₁₅, the fast-declining objects seem to lie on a smooth extension to the color loci of the normal objects.

For reasons stated above, we prefer to be agnostic about the reddening and go forward with the results obtained without using the LRS. We also adopt the Cauchy prior for the reddening as it allows a higher probability for larger reddenings. Table 2 lists the fits to the intrinsic color coefficients when using different subsamples and assuming different priors. It also tabulates the intrinsic scatter in the colors. These vary between 0.05 and 0.15 mag. The optical colors B − V and g − r have the lowest intrinsic widths, followed by the optical-NIR colors. The u − B color stands out as having a particularly large scatter. We consider the results for the reddening law in the next sections.

If we use a simple uniform prior for R_V that is universal for the entire sample of SNe Ia, we find that the recovered value depends on the subsample of objects we fit and the form of the reddening law. Table 3 lists the resulting values of R_V for these different cases. The sensitivity to subsample was seen in our earlier work (F10) and continues with our larger sample of objects. Generally speaking, when we include those SNe Ia with pseudo-colors redder than B_max − V_max ≃ 0.5, the resulting value for R_V tends to be lower.

In all cases, the F99 reddening law results in higher R_V than CCM+O. This is simply due to the different behaviors of either reddening law as we extrapolate to R_V < 2.6. All we can say is that when the reddest SNe Ia are removed from the analysis, the resulting value of R_V = 2.15 ± 0.16 for the F99 reddening law is well below the typical value for the MW. For the remainder, we will primarily use F99 to describe ISM dust, however we also include results using CCM+O in our tables for comparison and in order to be consistent with previous publications.

The dependence of R_V on subsample was seen in F10 and it was argued that a possible explanation was that the reddest SNe Ia could have increased extinction due to the proximity of dust local to the progenitor system (G08). Indeed, a detailed analysis in F10 of the extinction as a function of wavelength seemed to favor the G08 power-law over the standard CCM+O law due to its ability to better fit the u-band extinction. We therefore also attempt to fit the G08 reddening law:

$$\begin{equation} \frac{A_{j}}{A_{V}}=1-a+a\left(\frac{\lambda _{j}}{\lambda _{V}}\right)^{-p}, \end{equation} \tag{ 9 }$$

where λ_j is the effective wavelength of filter j. The differential extinction is therefore

$$\begin{equation} \Delta A_{j}=A_{V}a\left(\left(\frac{\lambda _{B}}{\lambda _{V}}\right)^{-p}-\left(\frac{\lambda _{j}}{\lambda _{V}}\right)^{-p}\right), \end{equation} \tag{ 10 }$$

where we see that A_V and a are completely degenerate. We can therefore only make inferences on the product of these two parameters. Column 4 of Table 3 lists the best-fit values for the power index p.

As with fitting for R_V, the resulting values of p depend quite a bit on the data. We fit with the same four scenarios as before. As with CCM+O and F99, removing the reddest objects results in a less steep reddening law, and there is a slight dependence on inclusion of the u band.

While it remains true that the u-band color excess for the two reddest objects (SN 2005A and SN 2006X) are better fit by the G08 power-law, these seem to be the only two and it may simply be a statistical fluke. Figure 16 shows the residuals in the color model for u − B. Clearly, SN 2005A and SN 2006X show a better fit, but very few others show a similar improvement. SN 2006X is known to have CSM (Wang et al. 2008), which was the motivation for the G08 power-law, so it is interesting that it fits better. However, SN 2007le is also known to have CSM (Simon et al. 2009), and yet it is fit equally well by CCM+O and F99. It is also conspicuous that the largest outliers are all spectroscopically classified as broad line (BL) and high-velocity (HV). This could indicate that it is spectral features that are responsible for the non-standard u − B colors. To investigate this, we plot the u − B residuals as a function of the pseudo-equivalent width of the Ca ii H & K lines, pW1 from Folatelli et al. (2013). This is shown in Figure 17. Clearly, the objects with high pW1 have the highest residuals. We therefore conclude that the anomalous u − B colors of these BL objects are due to their prominent Ca ii lines rather than an anomalous reddening law due to CSM.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

As a check, we can easily estimate the change in u − B color due to a change in pW1. We start by measuring synthetic u and B flux from the Hsiao et al. (2007) SED. We then artificially remove the Ca ii H & K lines by interpolating the continuum. The flux is then re-measured and a flux decrement in each filter is computed. This flux decrement is scaled from the Hsiao et al. (2007) pW1 = 108 Å to pW1 = 200 Å and removed from the original u and B fluxes. The resulting change in u is 0.33 mag fainter, while the change in B is 0.03 mag fainter, so that the change in u − B is 0.3 mag redder, consistent with the effect seen.

In a recent paper, Mandel et al. (2014) presented evidence that BL–HV events show intrinsic color differences in B − V and B − R of ∼0.05–0.10 mag with respect to normal SNe Ia, and that these differences originate in the B band. The intrinsic colors for the CSP sample plotted in Figure 13 show similar deviations in B − V for the same BL–HV events in our sample that display anomalous u–B colors. Foley & Kasen (2011) have argued that this is a natural expectation for HV objects since for saturated lines, as ejecta velocity increases, the width of the line will also increase. This, in turn, leads to larger equivalent widths and generally higher opacity in wavelength regions where there are many strong lines. This effect is particularly strong in u due to the strength of the Ca ii H & K absorption.

This trend with Ca ii can also be seen in the work by Chotard et al. (2011), who modeled a reddening law based on spectroscopic observations of SNe Ia. After correcting for decline rate, the reddening law had sharp features due to variations in the strength of Ca ii H & K, Ca ii NIR triplet, and Si ii λ6355. They apply a color correction based on the Ca ii H & K strength, which is precisely what we see in Figure 17.

Given that the derived value of R_V depends on the sample used, it is natural to ask whether certain objects are driving the solution, particularly those with larger color excesses which have a greater "pull." In Figure 18, we show a montage of the extinctions in different NIR filters relative to E(B − V)_host. In these extinction-extinction plots, a fixed reddening law, R_V, should produce a correlation with fixed slope equal to (R_X − R_Y) where X and Y correspond to the particular filter combination and the R_λ values are dependent on the value of R_V through the reddening law. Clearly, there is a large spread in the correlation with a tendency to small values of R_V. Because of this, we allow the reddening coefficient R_V to vary for each SN Ia in the sample.

Zoom In Zoom Out Reset image size

As a first step, we again use uniform priors on the individual R_V. For objects with very little extinction or no NIR photometry, R_V will be very poorly constrained. Nevertheless, the inferred E(B − V)_host and intrinsic color loci will be well-defined. We are primarily interested in those objects for which R_V can be constrained and require only uniform priors.

The results of fitting E(B − V)_host and R_V are shown in Figure 19 and tabulated in Columns 6, 7, and 10 of Table 1. As expected, the value of R_V becomes unconstrained as the amount of reddening tends to zero. In general, it would seem that the objects with moderate color excesses (0.2 < E(B − V)_host < 0.5) have a range of R_V from the rather low values that have been endemic of SN Ia studies to values more in line with the canonical MW value of R_V = 3.1, though these seem to be in the minority. All the more reddened objects tend to the low value of R_V ∼ 1.7. Figure 20 shows a plot similar to Figure 14 in F10 where we plot the color excess in different filters as a function of inverse wavelength. We have normalized the color excesses by E(B − V)_host in order to focus on the shape of the extinction curve. CCM+O curves with a range of values of R_V are shown for reference. Clearly, the colors of some objects like SN 2006gj are consistent with values closer to the MW value of R_V = 3.1 whereas colors of objects like SN 2007le would imply a lower value of R_V ∼ 1.7. Others have attempted to model a relationship between R_V and the observed colors (Mandel et al. 2011), however in the absence of any theoretical model to motivate such a relationship, and significantly fewer objects in our sample, we will not do so at this time. It is fair to say, however, that a simple one-parameter description of the extinction suffered by SNe Ia is insufficient to explain the diversity of inferred values of R_V.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

There are four objects which prefer very low values of R_V, such that they give us upper limits. These objects, SN 2004eo, SN 2005ke, SN 2007ba, and SN 2008gp, have moderately red optical colors, but very blue optical-NIR colors. These cannot be explained with any reasonable reddening law and must be due to intrinsically peculiar colors or systematic errors in the NIR photometry. We have already argued that SN 2005ke is a transitional object between fast and normal SNe Ia and SN 2007ba has a similar s_BV. They, and SN 2004eo, are classified as spectroscopic-type CL, while SN 2008gp is unclassified. However, there are other CL-types that have ranges of R_V that are perfectly normal. Unfortunately, we do not possess NIR spectra of these objects and cannot determine whether there are spectroscopic peculiarities in the NIR.

Throughout this paper, we have assumed foreground MW color excesses from Schlafly & Finkbeiner (2011). While we include the error for these measurements in our analysis, one might be concerned that there may be a systematic difference between the true color excess and the Schlafly & Finkbeiner (2011) values. In order to test this, we can look at objects whose host-galaxy extinctions are expected to be low (due to host-type, projected distance, and lack of Na i absorption), yet have significant MW extinction according to Schlafly & Finkbeiner (2011).

There are three such objects in our sample: SN 2006kf (E(B − V)_MW = 0.210), SN 2008bc (E(B − V)_MW = 0.225), and SN 2008ia (E(B − V)_MW = 0.195). In a trial run, these objects' foreground extinctions were artificially set to zero, allowing us to determine the foreground extinction based on the observed colors of the SNe rather than extinction maps. The results are given in Table 4. The values of E(B − V)_host are quite close and show no sign of significant systematic difference, though it is difficult to be certain with only three objects. The derived values of R_V are consistent with the spread in values associated with MW sight-lines (Cardelli et al. 1989; Fitzpatrick 1999).

Table 4. Inferred Milky Way Color Excesses

SN	E(B − V)S11	E(B − V)	RV
2006kf	0.210(0.002)	0.19(0.02)	$4.1^{+1.1}_{-0.5}$
2008bc	0.225(0.004)	0.25(0.03)	$2.1^{+0.5}_{-0.3}$
2008ia	0.195(0.005)	0.24(0.03)	$3.9^{+0.6}_{-0.4}$

Notes. Column 1: IAU name; Column 2: Milky Way E(B − V) color excesses from Schlafly & Finkbeiner (2011); Column 3: Milky Way E(B − V) color excesses inferred from SN colors; Column 4: Milky Way R_V inferred from SN colors.

Download table as:  ASCIITypeset image

We finish with an attempt to construct a prior for R_V that can be used when either (1) the extinction is too low or (2) there is insufficient filter coverage. We do this in two ways: by constructing a Gaussian mixture model as done in Kelly (2007) for all objects and splitting the objects into bins of E(B − V)_host and constructing an independent Gaussian distribution for each bin.

The Gaussian mixture model is simply a sum of N Gaussians, forming a composite prior that can have larger wings or multiple peaks. Details of the construction of this prior can be found in the Appendix and Kelly (2007). We run our MCMC analysis using the same sub-samples and filter sets as in Section 4.4.3. We find that two Gaussian components are sufficient, though the second is barely significant. Table 5 summarizes the values of the mixture model's hyper-parameters. It is evident that the parameters are quite insensitive to which subsample and filter set we use. We can therefore use any to represent a prior on R_V for future studies.

Our second approach is to bin the data by E(B − V)_host and solve for an independent single Gaussian prior for each bin. The hope is that the larger numbers of objects at low color excesses will balance their relatively weak pull on the data and that the reddest objects will not overly influence their prior. We split the data into five bins. The results are given in Table 6 and plotted in Figure 19 as red circles denoting the mean and error bars denoting the standard deviation of each prior. Clearly, the average R_V diminishes with E(B − V)_host, starting with typical R_V ∼ 3–4 found in the MW for the least reddened objects. It is hard to imagine a scenario in which such a correlation should exist and could simply be due to low numbers highly reddened objects and biases in follow-up selection (e.g., highly reddened, fainter objects close to their host galaxy are selected against for spectroscopic follow-up and typing). One should also note that all the BL and HV objects have low R_V, while the other objects show more of a spread. Deciding whether this is due to differences in ISM environment or spectroscopic diversity (SN physics) will require more objects with moderate to high color excesses in a variety of locations within their hosts.

The values of R_V for each object that result from the mixture model prior and binned prior are tabulated in Columns 8–12 of Table 1.

With the probability of the data defined, all that remains is to specify the priors. Wherever possible, we employ uniform priors, relying on the data to constrain the parameters.¹¹ This is done for the N + 1 coefficients of the polynomial $P_{j}^{N}$ for each intrinsic color, the true values of the color-stretch s_{BV, i}, and the MW color excess E(B − V)_{MW, i}. For the host-galaxy color excess, we use one of three possible priors. The first is a conditional prior such that objects in the LRS are assigned zero color excess and the rest are drawn from a uniform prior:

$$\begin{eqnarray} p\left(E(B-V)_{{\rm host},i}\right) & = & \delta (0),\quad i\in \mathrm{LRS} \nonumber \\ & = & U\left(0,\infty \right)\quad i\notin \mathrm{LRS}, \end{eqnarray} \tag{ A6 }$$

where δ() is the Dirac delta function. The second is the same prior as used by Jha et al. (2007), namely an exponential distribution with scale τ:

$$\begin{equation} p\left(E(B-V)_{{\rm host},i}\right)\propto \exp \left(-\frac{E(B-V)_{{\rm host},i}}{\tau }\right). \end{equation} \tag{ A7 }$$

Lastly, we investigate using a prior that is similar to the exponential prior, but with longer tail, the truncated Cauchy distribution:

$$\begin{eqnarray} p(E(B&-&V)_{{\rm host},i}|\tau,E_{{\rm max}}) = \left[\tau \arctan \left(\frac{E_{{\rm max}}}{\tau }\right)\right. \nonumber\\ &&\left. \!\!\! \times \left(1+\left(\frac{E(B-V)_{{\rm host},i}}{\tau }\right)^{2}\right)\right]^{-1}, \nonumber\\ && \! 0<E(B-V)_{{\rm host},i} <E(B-V)_{{\rm max}} \, {=}\, 0,\ \mathrm{otherwise}, \nonumber \\ \end{eqnarray} \tag{ A8 }$$

where again τ represents a scale length and we impose a maximum E(B − V)_max to ensure finite probability. In both the exponential and Cauchy priors, an inverse-gamma distribution is used for the hyper-parameter τ:

$$\begin{equation} p\left(\tau \right)=\Gamma ^{-1}\left(\alpha =1,\beta =0.2\right), \end{equation} \tag{ A9 }$$

where α and β are chosen to give a relatively broad prior peaked near the typical color excesses found in SNe Ia.

In this paper, we investigate several priors for R_V that are motivated by previous work done in the field. The first and simplest is the assumption that there is a single universal R_V for all SNe Ia. We only insist that it be strictly positive:

$$\begin{equation} p\left(R_{V}\right)=U\left(0,\infty \right). \end{equation} \tag{ A10 }$$

We therefore solve for the most likely value for R_V and its error. This is the prior typically used in most light-curve fitters, SNooPy included. A more realistic assumption would be that there is some intrinsic distribution of R_V and we instead try to infer the properties of this distribution. Inspired by the work of Kelly (2007), we use "Gaussian mixture" prior, which is a sum of Gaussians, each with its own mean, standard deviation, and normalization:

$$\begin{equation} p\left(R_{V,i}\left|\pi _{k},\mu _{k},\sigma _{k}^{2}\right.\right)\propto \sum _{k}\pi _{k}\exp \left(-\frac{1}{2\sigma _{k}^{2}}(R_{V,i}-\mu _{k})^{2}\right). \end{equation} \tag{ A11 }$$

The hyper-parameters μ_k, $\sigma _{k}^{2}$, and π_k are themselves given the following hierarchical priors:

$$\begin{eqnarray} p\left(\pi _{k}\right) & = & \mathrm{Dirichlet\left(1,\ldots,1\right),} \nonumber \\ p\left(\mu _{k}\left|\mu _{0},\sigma _{0}^{2}\right.\right) & = & N\left(\mu _{0},\sigma _{o}^{2}\right), \nonumber \\ p\left(\sigma _{k}^{2}\left|w^{2}\right.\right)=p\left(\sigma _{0}^{2}\left|w^{2}\right.\right) & = & \mathrm{Inv}\chi ^{2}(1,w^{2}), \nonumber \\ p\left(\mu _{0}\right) & = & U\left(-\infty,\infty \right), \nonumber \\ p(w^{2}) & = & U\left(0,\infty \right), \end{eqnarray} \tag{ A12 }$$

where the Dirichlet prior simply ensures ∑_kπ_k = 1 and N is the normal distribution. Though complicated, this set of priors simply ensures that the difference between the means of the Gaussians μ_k are on the order of their widths $\sigma _{k}^{2}$. The parameter μ₀ controls the overall mean of the composite distribution while w₂ controls its overall width and separation of local maxima.

The final prior we use reflects the findings of Mandel et al. (2011), namely that R_V tends to be smaller on average for larger E(B − V)_host. We therefore construct a prior in which we bin the data based on the value of E(B − V)_host and assign an independent prior for each bin. The prior for each bin is taken to be a single Gaussian with independent mean and standard deviation:

$$\begin{eqnarray} p\left(R_{V,i}\left|E(B-V)_{{\rm host},i},\mu _{k},\sigma _{k}^{2}\right.\right) &=& N\left(\mu _{k},\sigma _{k}^{2}\right),\; \nonumber\\ x_{k-1} \;&{<} & E(B-V)_{{\rm host},i}<x_{k}, \nonumber\\ \end{eqnarray} \tag{ A13 }$$

where the x_k denote the bin boundaries.