## On standardizing kilonovae and their use as standard candles to measure the Hubble constant

Michael W. Coughlin, 1, 2 Tim Dietrich, 3 Jack Heinzel, 4, 5 Nandita Khetan, 6 Sarah Antier, 7 Mattia Bulla, 8 Nelson Christensen, 4, 5 David A. Coulter, 9 and Ryan J. Foley 9

1

School of Physics and Astronomy, University of Minnesota, Minneapolis, Minnesota 55455, USA 2 California Institute of Technology, 1200 East California Blvd, MC 249-17, Pasadena, CA 91125, USA 3 Nikhef, Science Park, 1098 XG Amsterdam, The Netherlands 4 Artemis, Universit´ e Cˆ ote d'Azur, Observatoire Cˆ ote d'Azur, CNRS, CS 34229, F-06304 Nice Cedex 4, France 5 Carleton College, Northfield, MN 55057, USA 6 Gran Sasso Science Institute (GSSI), I-67100 L'Aquila, Italy 7 APC, UMR 7164, 10 rue Alice Domon et Lonie Duquet, 75205 Paris, France 8 Nordita, KTH Royal Institute of Technology and Stockholm University, Roslagstullsbacken 23, SE-106 91 Stockholm, Sweden 9 Department of Astronomy and Astrophysics, University of California, Santa Cruz, CA 95064, USA

The detection of GW170817 is revolutionizing many areas of astrophysics with the joint observation of gravitational waves and electromagnetic emissions. These multi-messenger events provide a new approach to determine the Hubble constant, thus, they are a promising candidate for mitigating the tension between measurements of Type Ia supernovae via the local distance ladder and the Cosmic Microwave Background. In addition to the 'standard siren' provided by the gravitational-wave measurement, the kilonova itself has characteristics that allow to improve existing measurements or to perform yet another, independent measurement of the Hubble constant without gravitationalwave information. Here, we employ standardization techniques borrowed from the type-Ia community and apply them to kilonovae, not using any information from the gravitational-wave signal. We use two versions of this technique, one derived from direct observables measured from the lightcurve, and the other based on inferred ejecta parameters, e.g., mass, velocity, and composition, for two different models. These lead to Hubble Constant measurements of H 0 = 109 +49 -35 km s -1 Mpc -1 for the measured analysis, and H 0 = 85 +22 -17 km s -1 Mpc -1 and H 0 = 79 +23 -15 km s -1 Mpc -1 for the inferred analyses. This measurement has error bars within ∼ 2 to the gravitational-wave measurements ( H 0 = 74 +16 -8 km s -1 Mpc -1 ), showing its promise as an independent constraint on H 0 .

PACS numbers: 95.75.-z,04.30.-w

A precise knowledge of the Hubble constant ( H 0 ), to determine the expansion rate of the Universe, is one of the most important measurements driving the study of cosmology [1, 2]. It has been known for a long time that the combined detection of gravitational waves (GWs) and their potential electromagnetic counterparts are useful for measuring the expansion rate of the universe [3]. These measurements are interesting since the GW standard siren measurements of H 0 do not rely on a cosmic distance ladder and do not assume any cosmological model as a prior.

This measurement has been made possible by the detection of GW170817 [4] and AT2017gfo, a 'kilonova,' which is thermal emission produced by the radioactive decay of neutron-rich matter synthesized from the ejecta of the compact binary coalescence at optical, nearinfrared, and ultraviolet wavelengths [5-21]. The analysis of GW170817 and the redshift of its host galaxy led to a measurement of H 0 = 74 +16 -8 km/s/Mpc (median and symmetric 68% credible interval), where degeneracy in the GW signal between the source distance and the weakly constrained angle of inclination between the total angular momentum of the binary and the line of sight dominated the H 0 measurement uncertainty [22].

It has been estimated that ∼ 50-100 GW events with identified optical counterparts would be required to have a H 0 precision measurement of ∼ 2% [23]. Of course, these searches are a significant observational challenge because one has to cover, over a short time interval, a large localization region, typically larger than the ∼ 20 square degrees for GW170817 [24, 25]. The resulting H 0 measurements can be improved with, for example, high angular resolution imaging of the radio counterpart. Hotokezaka et al. [26] applied this technique for GW170817 and obtained H 0 = 68 9 . +4 7 . -4 6 . km/s/Mpc.

In this letter, we employ techniques borrowed from the type-Ia supernova community to measure distance moduli based on kilonova lightcurves. When constraining distances with type-Ia supernovae, it is required to anchor them to some other 'primary' distance indicators in order to calibrate their absolute magnitude. In this work, we are using modelled light curves and here, the absolute magnitude is given by these models, and so we do not require to know their distances apriori by some other technique. The method relies on differences in the modeled lightcurves due to ejecta parameters such as the ejecta mass ( M ej ), ejecta velocity ( v ej ), and lanthanide fraction ( X lan ). Recently, it was also shown that

FIG. 1: On the left is the color-magnitude diagram for the models in [27]. On the right is ∆ m 7 (measured in K-band over 7 days) - magnitude diagram for the same. Lanthanide fraction values X lan = [10 -9 , 10 -5 , 10 -4 , 10 -3 , 10 -2 , 10 -1 ] are shown by the color bar, while for v ej = [0.03, 0.05, 0.1, 0.2, 0.3], the markers are circle, triangle, upside-down triangle, triangle pointing to the right, and triangle pointing to the left.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00889v4-figure-4a3e1a78-99e5-4718-9459-e727ce1fcd86.jpg)

a similar approach can be used for the bolometric luminosity to standardize kilonovae [28]. Qualitatively, the imprint of the ejecta properties can be modeled by the semi-analytic methods of Arnett [29], and more broadly by these models' success at predicting the lightcurves for GW170817 (e.g. Ref. [19]). Arnett-based models for kilonovae require a power term, which has generally been taken to be P ∝ t β as appropriate for radioactivity models [30], in addition to having other parameters such as the ejecta mass, the energy (or equivalently the velocity), and the opacity. Most important to the overall luminosity is the heating rate per mass of the ejecta, determined by the product of intrinsic decay power and thermalization efficiency. Estimates of the thermalization efficiency based on simulations exist [31], although they still are the largest systematic error budget, as the mass scales roughly inversely with the powering level. From Ref. [29], the diffusion timescale is τ ∝ ( κM v ) 1 / 2 , and similarly, the density is ρ ∝ M V ∝ M ( vτ ) 3 ; in this way, all of these quantities affect the observables. For the simplest model where all of the energy is injected at t = 0, corresponding to the time of peak luminosity, then the luminosity as a function of time L t ( ) ∝ L e 0 -t/τ , which implies that log( L t /L ( ) 0 ) ∝ -t/τ . This argues that the change in magnitude will be proportional to τ .

Based on this, we explore color-magnitude diagrams for kilonovae, with the idea that measuring time constants and colors may be useful for determining the underlying luminosities. In the left-hand panel of Fig. 1, we show these quantities for all of the spherical models made available in Ref. [27] plotting i -band minus K -band [52]. A few trends stand out: As expected, the simulations with lower X lan have lower absolute K -band magnitudes than those with higher X lan , with 1-2 mag differences seen depending on the lanthanide fraction. The much larger effect is on the color. From the lowest to highest X lan , the i -band minus K -band color can vary by up to 56 mag. In addition, the velocities predominantly change the color, but at a much lower level, changing the color glyph[lessorsimilar] 1 mag. The trend with M ej is a clear increase in peak magnitude, which is true of all lanthanide fraction and ejecta velocity pairs. Moreover, the overall K-band peak luminosity increases as X lan increases (or as one looks to the right in the grid). There is a similar but smaller trend with velocity.

In the right hand panel, we plot the change in luminosity, ∆ m 7 , between peak and 7 days later in K -band. As X lan decreases, the effects on ∆ m 7 increase, with ∆ m 7 glyph[greaterorsimilar] 1 mag at low X lan and ∆ m 7 glyph[lessorsimilar] 1 mag at higher X lan . We plot the peak K-band magnitudes vs. ejecta mass for the available lanthanide fractions and ejecta velocities of the employed simulation set in the Supplementary material (this is essentially the same plot as Fig. 1, but with the points separated out by ejecta velocity and lanthanide fraction).

The clear linear structure in Fig. 1 motivates the potential for their standardization, similar to SNe Ia lightcurves (see e.g. Ref. [32]). In the case of SN Ia lightcurves, they typically reach peak 17 days after explosion and then decay on a timescale of a few months. This motivates the use of the peak brightness, the time of the peak, and the 'width' of the lightcurve as characteristic variables that can be compared, as realized early on [1, 2]. Similar to the decline-rate parameter used in the SN 1a community (typically ∆ m 15 ), we will define a K-band decay parameter over 7 days as discussed above. This has the benefit of being measured from the observed light curve, with a downside of being tied to a particular filter and photometric system. It also requires the peak in

FIG. 2: The top panel shows the 'Fundamental plane' plot for the fit of Eqs. 1 and 3 and to the color-magnitude diagram shown in Figure 1 for the models in [27]. The 'measured' values are derived from direct observables measured from the lightcurve, and the 'inferred' values are from modeldependent ejecta parameters (mass, velocity, and lanthanide fraction). The bottom panel shows the difference between the computed values from the model and the fits for the simulations analyzed here.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00889v4-figure-8426dda0-eb48-48bb-8b69-c7780b03c506.jpg)

this passband to be well-measured, which is perhaps more straight-forward in the near-infrared where the peaks occur after a few days and therefore may be identified more easily. One downside might be that this band is less likely to be imaged in typical follow-up observations (i.e., before a kilonova has been confidently identified) because of the lack of infrared imagers on typical telescopes. Therefore, while our analysis makes one choice, there could be others more suitable for the particular observational situation.

For the moment, we will ignore the so-called K corrections that arise from the fact that the observed spectral energy distribution is redshifted by a factor (1 + z) and are effectively observing with filters that have been blueshifted in the rest frame of the kilonova. Given the local sensitivity volume for these kilonovae, this is reasonable, although their inclusion could be straight-forward by adopting a spectral energy distribution. A similar concern is that we implicitly ignore dust reddening in the analysis. In general, the changes in color due to the evolution of the kilonova is much faster than the evolution of the dust reddening timescales. This color term can be disentangled by reasonably high-cadence, multi-color imaging to measure the dust component that follows a typical reddening law and is constant in time.

We adopt two models where we do regressions in three dimensions in order to fit a distance for the kilonovae. The first is based purely on observable quantities, while the second will be on quantities inferred from the modelbased lightcurve fitting. We note that these points in parameter space are given equal weights in the model, and therefore implicitly assume that all portions of the parameter space are equally likely. To do this, we use a Gaussian Process Regression (GPR) based interpolation [33] to create a fit to the peak K-band magnitude for arbitrary parameters. This takes the form:

$$M _ { K = K _ { \max } } = & f ( \log _ { 1 0 } \Delta m _ { K = 7 }, \log _ { 1 0 } \Delta m _ { i = 7 }, \quad ( 1 ) \\ & m _ { i = K _ { \max } } - m _ { K = K _ { \max } } )$$

We employ a GPR based interpolation instead of a linear fit due to the significant covariance between parameters. As discussed above, the quantities we fit are ∆ m 7 in Kband, ∆ m 7 in i-band, and the i minus K -band color at peak in K-band. In order to compare with the apparent magnitudes, and therefore compute a distance modulus

$$M _ { \text{abs} } = m _ { \text{app} } - \mu$$

where m app is the apparent magnitude and µ = 5 log 10 ( D 10pc ) is the desired distance modulus (here, D is the distance). In principle, after the ∆ m 7 and colordependent corrections have been applied, there will remain an intrinsic dispersion in the lightcurves, perhaps arising from location in the galaxy or perhaps dependent on the galaxy type.

We compare this 'measured' fit to a fit based on 'inferred' quantities of ejecta mass, ejecta velocity, and lanthanide fraction:

$$M _ { K = K _ { \max } } = f ( \log _ { 1 0 } M _ { e j }, v _ { e j }, \log _ { 1 0 } X _ { \text{lan} } ) \quad ( 3 )$$

We note that this is not specific to these observables; for example, Ref. [34][53] use ejecta mass M ej , the temperature at 1 day after the merger T 0 [54], the halfopening angle of the lanthanide-rich component Φ (with Φ = 0 and Φ = 90 ◦ corresponding to one-component lanthanide-free and lanthanide-rich models, respectively) and the observer viewing angle θ obs (with cos θ obs = 0 and cos θ obs = 1 corresponding to a system viewed edgeon and face-on, respectively). We will compare these two models in the following. For consistency, we also employ a GPR based interpolation here, although we have found that a linear fit based on these variables gives comparable results. When performing the analysis, we include an overall error during the fit that represents scatter from intrinsic variability in the kilonovae models. These errors are ∼ 0.7 mag for the 'measured' case and ∼ 0.4 mag for the inferred case, which are derived from the median error reported by the GPR across the parameter space. Fig. 2 shows the performance of the fits of Eqs. (1) and (3) compared to the models in [27]. In general, there is some systematic scatter in the 'measured' case, indicating that the ∆ m 7 in K-band, ∆ m 7 in i-band, and the i minus K -band color at peak in K-band are not sufficient to completely standardize across the parameter space. In all likelihood, the proper choice of variables for this purpose will become more apparent as detections are made. The fits broken up by lanthanide fraction and velocity for the inferred and measured cases are shown in the Supplementary material.

FIG. 3: Posterior distributions for distance to GW170817, where the results of the GW-only analyses (high spin, dashed and low spin, solid) [35], the Kasen et al. [27] kilonova-only analyses (measured, dashed and inferred, solid), and the Bulla [34] kilonova analysis based on the fit of Eq. 3 are shown.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00889v4-figure-ba92ea61-3cc8-4435-b3af-4a74e273d619.jpg)

We use the fit of Eqs. (1) and (3) and apply them to the posteriors on ejecta mass, velocity, and lanthanide fraction derived previously [36-38]. In these fits, we assume a 1.0 mag systematic uncertainty on the model; in this sense, we do not require that the models are 'perfect' but instead encode some of the systematic uncertainties associated with them. We sample from the posteriors, in addition to the distributions for the measured parameters, to derive a constraint of D = 31 +17 -11 Mpc for the Kasen et al. [27] measured analysis. We find D = 37 +8 -7 Mpc and D = 40 +9 -8 Mpc for the Kasen et al. [27] and Bulla [34] inferred analyses respectively. This is consistent (given the relatively broad error bars) with other measurements of the host galaxy for GW170817, e.g. [39-41]. We show this constraint in Fig. 3 along with the high and low spin posteriors presented in [35].

Following the analysis of [22], we compute the corresponding values of the local Hubble constant for the kilonova analyses. For the measured analysis, we show the Kasen et al. [27] kilonova-only Hubble constant measurement of H 0 = 109 +49 -35 km s -1 Mpc -1 , median and symmetric 68% credible interval, in Fig. 4. For the inferred analyses, we also show the Kasen et al. [27] and Bulla [34] analyses, giving H 0 = 85 +22 -17 km s -1 Mpc -1 and H 0 = 79 +23 -15 km s -1 Mpc -1 . Models from Kasen et al. [27] and Bulla et al. [34] differ in several aspects, with the latter assuming parameterized opacities and temperatures, adopting different density profiles and ejecta geometries and taking into account the interplay between the two ejecta components. While the use of two different kilonovae models is not a truly independent analysis, the fact that they are different but the corresponding distributions consistent with one another gives confidence that the analysis is not particularly model dependent.

We also compute the corresponding values of the local Hubble constant for an analysis which combines the kilonova and GW derived distances (using the high-spin

FIG. 4: Posterior distributions for H 0 for GW170817, where the results of the GW-only analysis, the Kasen et al. [27] and Bulla [34] models kilonova-only analyses, and the Kasen et al. [27] combined GW-EM analysis are shown. The 1- and 2σ regions determined by the 'superluminal' motion measurement from the radio counterpart (blue) [26], Planck CMB (TT,TE,EE+lowP+lensing) (green) [42] and SHoES CepheidSN distance ladder surveys (orange) [43] are also depicted as vertical bands.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00889v4-figure-88707df2-8cfa-4455-9123-f2d103857210.jpg)

posteriors, as done in [22]). Because the GW and kilonova data are independent, the posterior probabilities for the distances can simply be multiplied. For the Kasen et al. [27] inferred analysis, we find a combined measurement of H 0 = 78 +14 -9 km s -1 Mpc -1 , while for the measured analysis, the results are the same as for the GW analysis. Altogether, this proves that the kilonova measurement is competitive with the GW measurements to obtain an independent constraint on H 0 .

In this article, we have demonstrated how to use parameters derived from radiative transfer simulations to give distance measurements using kilonovae, see also [28] for an alternative approach. We have adopted the peak luminosity in K-band, the decay in K-band over 7 days, ∆ m 7 , and the i minus K -band colors. We have shown that these distance estimates are consistent with other measurements of GW170817's host galaxy directly and provide competitive measurements of H 0 to GW distance measurements alone.

These techniques will play an important role in an era of both large-scale optical surveys, and radio-to-X-ray follow-up efforts, in identifying electromagnetic counterparts of compact binary coalescences. The Large Synoptic Survey Telescope (LSST) [44] and the Wide-Field InfraRed Survey Telescope (WFIRST) [45] in particular would be able to identify candidates with optical/NIR emission similar to GW170817 beyond 300 Mpc (i.e., the limit of current Advanced LIGO GW surveys for compact binary coalescences). Techniques such as those described here will play a key role in making these detections useful for studies of the Hubble constant and, thus, Dark Energy.

To explore this further, we examine the potential constraints over a sample of similar events. First of all, we assume that the 26% measurement here is representative of the typical width of the H 0 measurement from an individual event; this might be optimistic as AT2017gfo has a very well-sampled lightcurve. We also assume that the constraints converge as 1 √ N , where N is the number of kilonovae. Ref. [23] note that for the GW counterpart measurement, the typical tails in the distributions 'average out' to make a smaller effective width than is usual for any given single event, so this may be a conservative assumption. Under these assumptions, the number of events required to make a 6% and 2% measurement is ∼ 20 and ∼ 170 respectively. For the sake of comparison, Ref. [23] found that 50-100 binary neutron stars with counterparts are required to make a 2% measurement. They also found that for binary neutron stars without counterparts, where a galaxy-catalog based statistical method is required, ∼ 50 binary neutron stars are required to make a 6% measurement. In this sense, it is less powerful than binary neutron stars with counterparts but more so than binary neutron stars without counterparts. It is also likely that in practice, the error bars on a kilonova only analysis will become smaller with time as the models improve. Further detections of kilonovae will also make it possible to perform fits for Eq. (1) directly using data rather than using models, which may also yield improved constraints.

There are other techniques that could be employed that may yield further constraints. Ref. [28] uses a combination of semi-analytic kilonovae models with ejecta mass and velocity fits to numerical relativity simulations to show correlations between the maximum bolometric luminosity and decline from peak luminosity for a simulated population of binary neutron star mergers, showing that it may be possible to incorporate either models or future observations into the standardization of bolometric lightcurves. In addition, Ref. [46] laid the groundwork for techniques for performing joint standard candle-standard siren measurements that bring together the analyses of GWs and their electromagnetic counterparts, demonstrating in particular their importance in case where the selection effects from kilonova and GW observations can be significant.

In addition, these techniques can be used to augment GW detector calibration, which are currently limited to a few percent in amplitude and a few degrees in phase in the latest observing runs and are reaching their fundamental limits in performance [47, 48]. The match between the waveform predicted by general relativity and the GW strain signal alone can calibrate a single detector's relative amplitude and phase responses to a GW as a function of frequency [49]. This GW strain signal unaccompanied by other messenger signals can also calibrate relative responses between two GW detectors. However, when augmented by independent measurements of the event's distance and inclination angle, the detectors' absolute amplitude and phase responses to GWs can be calculated. Therefore, our techniques could be a critical method in the effort to calibrate responses of GW detectors using astrophysical signals. Although there will need to be many more GW detections with electromagnetic counterparts to improve the precision of the astrophysical calibration measurements over the current existing in situ measurements, single events can also be used to improve and verify in situ measurements. Ref. [49] showed that GW data alone can constrain the relative amplitude calibration uncertainty to less than 10% with 10-20 events (less than 1% with 1000-2000 events); when combining conservative constraints with broad distance and inclination constraints from electromagnetic data, it takes 400-600 events. Given the similar distance constraints assumed to those derived here, this is the right ballpark for this level of calibration uncertainty.

Further work is needed to understand how our restriction to a spherical geometry with a single-component would change for multiple components including possibly reprocessing between ejecta components, e.g. [50]. Discovering more kilonovae will improve the understanding of the effect of ejecta geometries, supporting the use of multiple components and informing how they should be included. In addition, the standardization assumes that the radiative-transfer models are consistent to produce proper absolute magnitudes and colors (at least within the assumed error bars), which motivates continued work to improve the accuracy of the models and the grids that they are simulated on. The potential for more sources, at the price of higher systematics, with this method leads to a trade-off between the use of kilonovae or purely GW measurements; GW measurements will have relatively lower levels of systematics at the price of fewer objects to use. This method can be used in particular for any kilonovae detection with measurements of the host galaxy redshift, by way of fits of the ejecta parameters (e.g. Ref. [51]). It may also be possible to constrain the inclination using associated GRB detections or upperlimits on potential sub-threshold gamma-ray transients. As both of these analyses will be directly testable by future kilonovae observations, the utility of analyses of this type will be dependent on the coming comparisons.

Acknowledgements. The authors would like to thank Daniel Kasen for making his models publicly available. MC is supported by the David and Ellen Lee Postdoctoral Fellowship at the California Institute of Technology. TDacknowledges support by the European Union's Horizon 2020 research and innovation program under grant agreement No 749145, BNSmergers. NC and JH acknowledge support from the National Science Foundation with grant number PHY-1806990. SA is supported by the CNES Postdoctoral Fellowship at Laboratoire Astroparticle et Cosmologie. The UCSC team is supported in part by NASA grant NNG17PX03C, NSF grants AST1518052 and AST-1911206, the Gordon &amp; Betty Moore Foundation, the Heising-Simons Foundation, and by a fellowship from the David and Lucile Packard Foundation to R.J.F. DC acknowledges support from the Na-

- [1] Adam G. Riess et al. Observational evidence from supernovae for an accelerating universe and a cosmological constant. Astron. J. , 116:1009-1038, 1998.
- [2] S. Perlmutter, G. Aldering, G. Goldhaber, et al. Measurements of Ω and Λ from 42 High-Redshift Supernovae. Astrophys. J. , 517:565-586, June 1999.
- [3] B. F. Schutz. Determining the Hubble constant from gravitational wave observations. Nature , 323:310, September 1986.
- [4] Abbott et al. Gw170817: Observation of gravitational waves from a binary neutron star inspiral. Phys. Rev. Lett. , 119:161101, Oct 2017.
- [5] K. D. Alexander, E. Berger, W. Fong, et al. The Electromagnetic Counterpart of the Binary Neutron Star Merger LIGO/Virgo GW170817. VI. Radio Constraints on a Relativistic Jet and Predictions for Late-time Emission from the Kilonova Ejecta. ApJL , 848:L21, October 2017.
- [6] Chornock et al. The electromagnetic counterpart of the binary neutron star merger ligo/virgo gw170817. iv. detection of near-infrared signatures of r -process nucleosynthesis with gemini-south. The Astrophysical Journal Letters , 848(2):L19, 2017.
- [7] P. S. Cowperthwaite, E. Berger, V. A. Villar, et al. The Electromagnetic Counterpart of the Binary Neutron Star Merger LIGO/Virgo GW170817. II. UV, Optical, and Near-infrared Light Curves and Comparison to Kilonova Models. ApJL , 848:L17, October 2017.
- [8] M. R. Drout, A. L. Piro, B. J. Shappee, et al. Light curves of the neutron star merger GW170817/SSS17a: Implications for r-process nucleosynthesis. Science , 358:15701574, December 2017.
- [9] P. A. Evans, S. B. Cenko, J. A. Kennea, et al. Swift and NuSTAR observations of GW170817: Detection of a blue kilonova. Science , 358:1565-1570, December 2017.
- [10] D. Haggard, M. Nynka, J. J. Ruan, et al. A Deep Chandra X-Ray Study of Neutron Star Coalescence GW170817. ApJL , 848:L25, October 2017.
- [11] G. Hallinan, A. Corsi, K. P. Mooley, et al. A radio counterpart to a neutron star merger. Science , 358:1579-1583, December 2017.
- [12] M. M. Kasliwal, E. Nakar, L. P. Singer, et al. Illuminating gravitational waves: A concordant picture of photons from a neutron star merger. Science , 358(6370):15591565, 2017.
- [13] Kilpatrick et al. Electromagnetic evidence that SSS17a is the result of a binary neutron star merger. Science , 358(6370):1583-1587, 2017.
- [14] R. Margutti, E. Berger, W. Fong, et al. The Electromagnetic Counterpart of the Binary Neutron Star Merger LIGO/Virgo GW170817. V. Rising X-Ray Emission from an Off-axis Jet. ApJL , 848:L20, October 2017.
- [15] C. McCully, D. Hiramatsu, D. A. Howell, et al. The Rapid Reddening and Featureless Optical Spectra of the Optical Counterpart of GW170817, AT 2017gfo, during the First Four Days. ApJL , 848:L32, October 2017.
- [16] Nicholl et al. The electromagnetic counterpart of the binary neutron star merger LIGO/Virgo GW170817. iii. optical and UV spectra of a blue kilonova from fast polar
- ejecta. The Astrophysical Journal Letters , 848(2):L18, 2017.
- [17] B. J. Shappee, J. D. Simon, M. R. Drout, et al. Early spectra of the gravitational wave source GW170817: Evolution of a neutron star merger. Science , 358:1574-1578, December 2017.
- [18] E. Pian, P. D'Avanzo, S. Benetti, et al. Spectroscopic identification of r-process nucleosynthesis in a double neutron-star merger. Nature (London) , 551:67-70, November 2017.
- [19] Smartt et al. A kilonova as the electromagnetic counterpart to a gravitational-wave source. Nature , 551:75 EP -, 10 2017.
- [20] E. Troja, L. Piro, H. van Eerten, et al. The X-ray counterpart to the gravitational-wave event GW170817. Nature (London) , 551:71-74, November 2017.
- [21] Y. Utsumi, M. Tanaka, N. Tominaga, et al. J-GEM observations of an electromagnetic counterpart to the neutron star merger GW170817. Proceedings of the Astronomical Society of Japan , 69:101, December 2017.
- [22] B. P. Abbott, R. Abbott, T. D. Abbott, et al. A gravitational-wave standard siren measurement of the Hubble constant. Nature (London) , 551:85-88, November 2017.
- [23] Hsin-Yu Chen, Maya Fishbach, and Daniel E. Holz. A two per cent Hubble constant measurement from standard sirens within five years. Nature , 562(7728):545-547, 2018.
- [24] Michael W. Coughlin, Tom´ as Ahumada, S. Bradley Cenko, et al. 2900 square degree search for the optical counterpart of short gamma-ray burst GRB 180523b with the zwicky transient facility. Publications of the Astronomical Society of the Pacific , 131(998):048001, feb 2019.
- [25] Coughlin et al. GROWTH on s190425z: Searching thousands of square degrees to identify an optical or infrared counterpart to a binary neutron star merger with the zwicky transient facility and palomar gattini-IR. The Astrophysical Journal , 885(1):L19, oct 2019.
- [26] K. Hotokezaka, E. Nakar, O. Gottlieb, et al. A hubble constant measurement from superluminal motion of the jet in GW170817. Nature Astronomy , 2019.
- [27] Daniel Kasen, Brian Metzger, Jennifer Barnes, Eliot Quataert, and Enrico Ramirez-Ruiz. Origin of the heavy elements in binary neutron-star mergers from a gravitational-wave event. Nature , 551:80 EP -, 10 2017.
- [28] Rahul Kashyap, Gayathri Raman, and Parameswaran Ajith. Can kilonova light curves be standardized? The Astrophysical Journal , 886(1):L19, nov 2019.
- [29] W. D. Arnett. Type I supernovae. I -Analytic solutions for the early part of the light curve. Astrophys. J. , 253:785-797, February 1982.
- [30] Brian D. Metzger. Kilonovae. Living Reviews in Relativity , 20(1):3, May 2017.
- [31] J. Barnes and D. Kasen. Effect of a High Opacity on the Light Curves of Radioactively Powered Transients from Compact Object Mergers. Astrophys. J. , 775:18, 2013.
- [32] Christopher R. Burns et al. The Carnegie Supernova

Project: Absolute Calibration and the Hubble Constant. Astrophys. J. , 869(1):56, 2018.

- [33] Zoheyr Doctor, Ben Farr, Daniel E. Holz, and Michael Prrer. Statistical Gravitational Waveform Models: What to Simulate Next? Phys. Rev. , D96(12):123011, 2017.
- [34] M. Bulla. POSSIS: predicting spectra, light curves, and polarization for multidimensional models of supernovae and kilonovae. Monthly Notices of the Royal Astronomical Society , 489(4):5037-5045, Nov 2019.
- [35] B. P. Abbott et al. Properties of the binary neutron star merger GW170817. Phys. Rev. , X9(1):011001, 2019.
- [36] M. Coughlin, T. Dietrich, K. Kawaguchi, et al. Toward Rapid Transient Identification and Characterization of Kilonovae. Astrophys. J. , 849:12, November 2017.
- [37] Coughlin et al. Constraints on the neutron star equation of state from at2017gfo using radiative transfer simulations. Monthly Notices of the Royal Astronomical Society , 480(3):3871-3878, 2018.
- [38] Michael W. Coughlin, Tim Dietrich, Ben Margalit, and Brian D. Metzger. Multimessenger Bayesian parameter inference of a binary neutron star merger. Mon. Not. Roy. Astron. Soc. , 489(1):L91-L96, 2019.
- [39] Jens Hjorth, Andrew J. Levan, Nial R. Tanvir, et al. The distance to NGC 4993: The host galaxy of the gravitational-wave event GW170817. The Astrophysical Journal , 848(2):L31, oct 2017.
- [40] A. J. Levan et al. The environment of the binary neutron star merger GW170817. Astrophys. J. , 848(2):L28, 2017.
- [41] Michele Cantiello et al. A Precise Distance to the Host Galaxy of the Binary Neutron Star Merger GW170817 Using Surface Brightness Fluctuations. Astrophys. J. , 854(2):L31, 2018.
- [42] P. A. R. Ade et al. Planck 2015 results. XIII. Cosmological parameters. Astron. Astrophys. , 594:A13, 2016.
- [43] A. G. Riess, L. M. Macri, S. L. Hoffmann, et al. A 2.4% Determination of the Local Value of the Hubble Constant. Astrophys. J. , 826:56, July 2016.
- [44] Z. Ivezic, J. A. Tyson, R. Allsman, J. Andrew, and R. Angel. LSST: from Science Drivers to Reference Design and Anticipated Data Products. 2008.
- [45] D. Spergel, N. Gehrels, C. Baltay, et al. Wide-Field InfrarRed Survey Telescope-Astrophysics Focused Telescope Assets WFIRST-AFTA 2015 Report. arXiv eprints , page arXiv:1503.03757, Mar 2015.
- [46] Zoheyr Doctor. Thunder and Lightning: Using NeutronStar Mergers as Simultaneous Standard Candles and Sirens to Measure Cosmological Parameters. arXiv eprints , page arXiv:1912.12218, Dec 2019.
- [47] Craig Cahillane, Joe Betzwieser, Duncan A. Brown, et al. Calibration uncertainty for advanced ligo's first and second observing runs. Phys. Rev. D , 96:102001, Nov 2017.
- [48] F Acernese, T Adams, K Agatsuma, et al. Calibration of advanced virgo and reconstruction of the gravitational wave signal h(t) during the observing run o2. Classical and Quantum Gravity , 35(20):205004, sep 2018.
- [49] Reed Essick and Daniel E Holz. Calibrating gravitational-wave detectors with GW170817. Classical and Quantum Gravity , 36(12):125002, may 2019.
- [50] Kyohei Kawaguchi, Masaru Shibata, and Masaomi Tanaka. Radiative transfer simulation for the optical and near-infrared electromagnetic counterparts to GW170817. The Astrophysical Journal , 865(2):L21, sep 2018.
- [51] Stefano Ascenzi et al. A luminosity distribution for kilo-
- novae based on short gamma-ray burst afterglows. Mon. Not. Roy. Astron. Soc. , 486(1):672-690, 2019.
- [52] https://github.com/dnkasen/Kasen Kilonova Models 2017
- [53] https://github.com/mbulla/kilonova models
- [54] In this analysis, we use simulations from [34] where the temperature is taken as a parameter rather than calculated from the adopted densities and nuclear heating rates. We will address this point and use self-consistent models in future analyses.

Ejecta Mass [solar mass]

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00889v4-figure-ddbf1979-132d-4496-9957-41c6a1b83b7d.jpg)

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00889v4-figure-65e70677-9a48-4bd4-9e4c-f06701390158.jpg)

FIG. 5: The fit of Eqs. (1) and (3) and to the color-magnitude diagram shown in Figure 1 for the models in Ref. [27], varying the available lanthanide fractions and ejecta velocities of the employed simulation set. The black points are the fits (with the measured error bar from the chi-squared) while the blue points are computed directly from the models in Ref. [27].
Michael W. Coughlin Tim Dietrich Jack Heinzel Nandita Khetan Sarah Antier Mattia Bulla Nelson Christensen David A. Coulter Ryan J. Foley
On standardizing kilonovae and their use as standard candles to measure the Hubble constant
The detection of GW170817 is revolutionizing many areas of astrophysics with
the joint observation of gravitational waves and electromagnetic emissions.
These multi-messenger events provide a new approach to determine the Hubble
constant, thus, they are a promising candidate for mitigating the tension
between measurements of Type Ia supernovae via the local distance ladder and
the Cosmic Microwave Background. In addition to the "standard siren" provided
by the gravitational-wave measurement, the kilonova itself has characteristics
that allow to improve existing measurements or to perform yet another,
independent measurement of the Hubble constant without gravitational-wave
information. Here, we employ standardization techniques borrowed from the
type-Ia community and apply them to kilonovae, not using any information from
the gravitational-wave signal. We use two versions of this technique, one
derived from direct observables measured from the lightcurve, and the other
based on inferred ejecta parameters, e.g., mass, velocity, and composition, for
two different models. These lead to Hubble Constant measurements of $H_0 =
109^{+49}_{-35}$\,km $\mathrm{s}^{-1}$ $\mathrm{Mpc}^{-1}$ for the measured
analysis, and $H_0 = 85^{+22}_{-17}$\,km $\mathrm{s}^{-1}$ $\mathrm{Mpc}^{-1}$
and $H_0 = 79^{+23}_{-15}$\,km $\mathrm{s}^{-1}$ $\mathrm{Mpc}^{-1}$ for the
inferred analyses. This measurement has error bars within ~$\sim$\,2 to the
gravitational-wave measurements ($H_0=74^{+16}_{-8}$\,km $\mathrm{s}^{-1}$
$\mathrm{Mpc}^{-1}$), showing its promise as an independent constraint on
$H_0$.