Galactic Extinction: How Many Novae Does It Hide and How Does It Affect the Galactic Nova Rate?

A. Kawash; L. Chomiuk; J. A. Rodriguez; J. Strader; K. V. Sokolovsky; E. Aydi; C. S. Kochanek; K. Z. Stanek; K. Mukai; K. De; B. Shappee; T. W.-S. Holoien; J. L. Prieto; T. A. Thompson

doi:10.3847/1538-4357/ac1f1a

1. Introduction

A classical nova occurs in an interacting binary system with a white dwarf primary, referred to as cataclysmic variable (CV; see Warner 1995). The white dwarf accretes material from the secondary, usually a low-mass main-sequence star, until a critical pressure and temperature are reached, leading to a thermonuclear runaway at the bottom of the hydrogen-rich shell accreted by the white dwarf (see Bode & Evans 2008 for a review). Expulsion of the accreted envelope occurs, causing the system to brighten significantly, by 5 to >19 mag (Kawash et al. 2021). Studies of classical novae in M31 have constrained the peak luminosities to a range from M_V ≈ −4 to −10 mag (Shafter 2017).

Historically, amateur astronomers have been the driving force in finding classical novae, with discoveries dating back thousands of years (Patterson et al. 2013; Shara et al. 2017). Novae began to be more systematically discovered in the mid-20th century, when, on average, ∼three per year were visually discovered. When film photography became commonly used in the 1980s and 1990s, there were ∼four discovered novae per year and then ∼eight per year in the 2000s and 2010s when digital cameras became widely available. Then, in 2017, the All-Sky Automated Survey for Supernovae (ASAS-SN) became the first survey to systematically observe the entire night sky with nearly daily cadence down to g ≈ 18.5 mag (Shappee et al. 2014), significantly deeper than most amateur observations. Since 2017, there have been roughly 10 classical novae discovered per year, on average.

The first prediction for the total frequency of Galactic nova eruptions was by Lundmark (1935), who estimated there should be about 50 novae yr^–1 (see Della Valle & Izzo 2020 for a review of Galactic nova rate predictions). Estimates from the early 1990s predicted much lower rates ranging from 11 to 20 yr^–1, derived from observations of other galaxies (Ciardullo et al. 1990; van den Bergh 1991; della Valle & Livio 1994). More recent surveys of M31 have increased these predictions for the Milky Way rate to ${34}_{-12}^{+15}$ yr⁻¹ (Darnley et al. 2006) and as high as ∼50 to ∼70 yr⁻¹ when accounting for incompleteness for faint and fast novae (Shafter 2017). Recent work modeling novae in our Galaxy predicted rates of ${50}_{-23}^{+31}$ yr⁻¹ from a sample of bright novae (Shafter 2017) and ${43.7}_{-8.7}^{+19.5}$ yr⁻¹ from a sample of novae detected in IR observations (De et al. 2021). If these recent, higher estimates are correct, novae could be key contributors to various isotopes present in the Galaxy (José et al. 2006), like ²⁶Al (José & Hernanz 1998; Bennett et al. 2013) and ⁷Li (Starrfield et al. 1978; Hernanz et al. 1996), but there must be a reason the majority of novae go undetected.

So, what is the cause of the discrepancy between the predicted and observed rates? One idea put forward is that the majority of observable classical nova events go undetected due to insufficient sky coverage of observations. However, the emergence of large sky surveys, including ASAS-SN, should solve this issue. The most common Galactic transient ASAS-SN discovers is a dwarf nova outburst, and Kawash et al. (2021) explored the possibility that some classical novae were being mistaken for dwarf novae. Though it is possible that a small number of novae can be misclassified, there are too few to significantly alter the discovery rate of classical novae.

Another possibility is that interstellar dust obscures a majority of classical novae that erupt in the Galaxy. This prospect is supported by the recently discovered sample of highly reddened and optically missed novae discovered by the Palomar Gattini IR (PGIR) survey (De et al. 2020a, 2021). These results suggest that dust could cause a substantial fraction of Galactic nova events to go undetected by optical observations, but how many remains an open question.

Shafter (2017) used an exponential disk to model extinction in the Galaxy, with 1 mag of extinction per kiloparsec in the midplane in the V band. This predicts that over 90% of all Galactic novae should get brighter than V = 18 mag, inconsistent with the results of De et al. (2021). The primary goal of this work is to investigate how utilizing a more robust Galactic extinction model changes these conclusions, and the availability of three-dimensional Galactic extinction maps now make it possible to explore this exact question.

Recent large photometric surveys combined with models of stellar spectral energy distributions have resulted in several three-dimensional extinction models of the Galaxy. Green et al. (2019) used a combination of Gaia, Pan-STARRS, and Two Micron All Sky Survey (2MASS) observations to model the extinction north of δ = −30°, and Marshall et al. (2006) used 2MASS data to model extinction in the inner Galactic plane. Neither of these maps cover the entire sky, but Bovy et al. (2016) combined these maps along with analytic models from Drimmel et al. (2003) to build an all-sky three-dimensional model of extinction in the Galaxy.

The Green et al. (2019) three-dimensional dust map was used in De et al. (2021) to build a model of Galactic novae in the PGIR field of view (δ > −28 fdg 9). This model was then compared to the PGIR sample of novae discovered at the J band and optically discovered novae. De et al. (2021) found that the amount of extinction predicted by the model was much more consistent with the PGIR sample of novae, implying that IR observations are sensitive to a significant fraction of novae that optical observations are not, but this analysis did not cover a large portion of the bulge and inner disk. Here we extrapolate to the entire sky using the Bovy et al. (2016) all-sky dust map to build upon these recent findings.

The goal of this work is to model the distribution of novae and dust within the Galaxy and explore how that combination affects how bright novae are when observed in the optical. In Section 2, we discuss the implementation and assumptions of the model, including the resulting apparent magnitude distribution of novae, how many optically missed novae we should expect, and where to find them with IR or optical observations in redder bands. Then, in Section 3, we use our model along with observational constraints from various surveys to estimate the global frequency of nova eruptions in the Galaxy. Lastly, in Section 4, we explore how sensitive our results are to our model assumptions.

2. Nova Model/Results

Here we discuss the components and assumptions that go into our Galactic nova model, the resulting Galactic apparent magnitude distribution of novae, and the spatial distribution of optically observable versus unobservable novae.

2.1. Stellar Density Profile

To analyze the effects of extinction on novae, we first must assume some distribution of novae within the Galaxy. For our primary model, we simply assume that the distribution of novae follows the distribution of stellar mass, but in Section 4.3, we also consider additional models with a higher rate of nova production per unit mass in the bulge (as compared with the disk). The stellar mass distribution is inferred by implementing a version of the Besançon Galactic mass model with a two-component bulge from Simion et al. (2017) and a thin disk, thick disk, and halo component from Robin et al. (2003). An explanation of the stellar density model, including the model parameters assumed, can be found in Appendix A.1. The distribution of mass is calculated on a three-dimensional Cartesian grid with size (x, y, z) = (30, 30, 30) kpc and a resolution of 0.1 kpc. The total mass of each component is calculated and shown in Table A1. In Section 4.2, we also distribute novae based on the contracted halo version of the mass profile presented in Cautun et al. (2020) to investigate how our results depend on the assumed distribution of novae, and this does slightly affect the resulting apparent magnitude distribution.

2.2. Extinction Model

To date, there is no single three-dimensional dust map that models Galactic extinction across the entire sky. However, the maps of Drimmel et al. (2003), Marshall et al. (2006), and Green et al. (2015) were combined and made publicly available at http://github.com/jobovy/mwdust to provide a stitched-together map over the entire sky (Bovy et al. 2016). This model, hereafter referred to as mwdust, can use several different map combinations, and here we use the combined19 version. This uses the updated map of the sky north of decl. δ = −30° from Green et al. (2019), the Marshall et al. (2006) maps covering the sky around the Galactic center −100° ≤ l ≤ 100° and −10° ≤ b ≤ 10°, and the Drimmel et al. (2003) map for the rest of the sky not covered by the first two. The Marshall et al. (2006) map takes precedence over the Green et al. (2019) map where they overlap; because the former relies more heavily on near-infrared data, it is able to model extinction out to larger distances at low latitude. A detailed explanation of the default model can be found in the Appendix of Bovy et al. (2016).

The mwdust model can estimate the extinction in a wide array of observing bands by assuming an extinction law with ${A}_{V}/{A}_{{K}_{s}}=8.65$ . However, Nataf et al. (2016) argued that toward the inner Galaxy, the extinction law varies, and they found a median value of ${A}_{V}/{A}_{{K}_{s}}=13.44$ , 55% larger than the default mwdust value. This roughly corresponds to a 25% increase in the extinction in the ASAS-SN g-band filter, A_g, so we also explore how our results change from this steeper extinction curve. This change introduces additional systematic uncertainties in the model, because extinction in the various maps that make up mwdust are measured in different filters and with different techniques. In addition, not all novae in our model are in the inner Galaxy. Nonetheless, this alternate model roughly captures the effects of a different Galactic extinction curve on our ability to find novae.

2.3. Positions and Distances

We ran Monte Carlo simulations of 1000 Galactic novae by probabilistically distributing them following the stellar mass of the Besançon Galactic model from Robin et al. (2003). The resulting distribution of novae in Galactic cylindrical coordinates for the primary model and a secondary stellar density model from Cautun et al. (2020) are shown in Figure 1. In Figure 2, these positions are transformed to sky coordinates in the reference frame of the Sun at R_⊙ = 8.122 kpc (Gravity Collaboration et al. 2018; Cautun et al. 2020) in the top panel, and a face-on view of the Galaxy is shown in the bottom panel. As expected, novae hug the disk plane, and there is a large increase in density toward the Galactic center.

**Figure 1.** Distributions of Galactic positions in cylindrical radii (top) and height above the disk (bottom) for 1000 randomly sampled novae distributed by randomly sampling from the Robin et al. (2003) model, shown as the blue dashed line, and the Cautun et al. (2020) model, shown as the red solid line.
Download figure:
Standard image High-resolution image

**Figure 2.** Top: positions of N = 1000 simulated novae in Galactic coordinates, distributed by randomly sampling the Robin et al. (2003) stellar density model. Bottom: same as top panel but in an external face-on view of our Galaxy. Novae are plotted in blue, and the position of the Sun is shown by an orange plus sign.
Download figure:
Standard image High-resolution image

The distances to Galactic novae are often hard to constrain even with Gaia parallaxes (Schaefer 2018). Figure 3 shows the expected distribution of distances based on our primary model. We expect the median distance to a nova to be 8.7 kpc, 68% of novae to have distances between 6 and 13 kpc, and 95% of novae to be within 16 kpc of the Sun. Also shown in Figure 3 is the distribution of distances from a magnitude-limited sample of novae that we expect to more closely resemble the sample of optically discovered novae. The median distance of this distribution is 7.9 kpc, very similar to the global population, and we expect 95% of this magnitude-limited sample to be within 13 kpc. Overall, the distribution of distances for a magnitude-limited sample is not significantly different than the global population, suggesting that distance alone is not the determining factor for missed novae.

2.4. Brightness Distribution

The peak absolute magnitude of each nova is randomly sampled from a normal distribution with a mean and standard deviation of μ = −7.2 and σ = 0.8 mag, respectively (Shafter 2017). This distribution was derived for M31 novae in the V band, and we assume there is no magnitude difference between this and Milky Way novae at the g-band peak (van den Bergh & Younger 1987; Miroshnichenko 1988; Hachisu & Kato 2014). In deriving Galactic nova rates, both Shafter (2017) and De et al. (2021) explored altering the luminosity function for the bulge and disk novae together and separately. We only assume the above luminosity function for all of the novae in our model, as De et al. (2021) found that it was not a significant factor for their results; however, we plan to investigate the effects of the luminosity function in A. Kawash et al. (2022, in preparation). After a Galactic position was randomly assigned to the nova, the accompanying extinction for that line of sight and distance was estimated using the mwdust package. These values were combined with the distance and randomly assigned absolute magnitude to estimate the peak apparent magnitude for each nova.

The cumulative distribution of the peak apparent magnitudes of Galactic novae is shown in blue and red for an assumed extinction curve of ${A}_{V}/{A}_{{K}_{s}}=8.65$ and 13.44, respectively, in Figure 4 and compared to a distribution excluding dust and a distribution implementing the disk extinction model of Shafter (2017). It is clear that the disk model of extinction utilized in Shafter (2017) vastly underestimates the effects of dust relative to the estimates from three-dimensional dust maps. Specifically, the extinction from the three-dimensional dust maps predicts that only 53%, to as low as 46% for the steeper extinction curve, of novae in the Galaxy will get brighter than g = 15 mag, while Shafter (2017) predicted that 82% of novae will be brighter than g = 15 mag. This could explain why ASAS-SN observations have not resulted in a significant increase in the nova discovery rate, and it is consistent with the scenario that a large fraction of nova eruptions are too faint to be detected in blue optical bands but are detectable in the IR, like the recently discovered PGIR sample (De et al. 2021). Also, this likely means that the deeper observations of the Vera C. Rubin Observatory Legacy Survey of Space and Time (LSST; Tyson 2002) will discover many more Galactic novae than previously thought if the bulge and plane are observed at a moderate cadence. Also shown in Figure 4 is the distribution of peak apparent magnitudes for novae discovered since 2000. As expected, observations are heavily biased toward detection of the brightest novae, so most faint and highly extinguished novae are undetectable by optical observations.

**Figure 4.** Cumulative distributions of peak apparent magnitudes of N = 1000 simulated novae compared to the distribution of observed novae discovered since 2000 (orange dashed–dotted line). Note that this observed distribution of peak brightness is measured in various observing bands, and no color correction has been made. Four models are shown: one that includes no dust (pink dotted line), one with an exponential disk as implemented in Shafter (2017; brown dashed line), and the `mwdust` model using the default reddening law (blue solid line) and a higher value of ${A}_{V}/{A}_{{K}_{s}}=13.44$ (red dotted line). The different models yield significantly different distributions, highlighting the importance of accurately modeling dust to estimate the nova rate.
Download figure:
Standard image High-resolution image

**Figure 4.** Cumulative distributions of peak apparent magnitudes of N = 1000 simulated novae compared to the distribution of observed novae discovered since 2000 (orange dashed–dotted line). Note that this observed distribution of peak brightness is measured in various observing bands, and no color correction has been made. Four models are shown: one that includes no dust (pink dotted line), one with an exponential disk as implemented in Shafter (2017; brown dashed line), and the `mwdust` model using the default reddening law (blue solid line) and a higher value of ${A}_{V}/{A}_{{K}_{s}}=13.44$ (red dotted line). The different models yield significantly different distributions, highlighting the importance of accurately modeling dust to estimate the nova rate.
Download figure:
Standard image High-resolution image

The accuracy of this model distribution relies heavily on the ability of the three-dimensional dust maps to estimate high extinctions at low latitudes out to large distances. A majority of the novae (94%) are in the area of the sky that the Marshall et al. (2006) map covers, and for regions at high column densities, this map only has information out to ∼7 kpc. Only 5% of novae fall within the Green et al. (2019) region, and this model only extends to a few kiloparsecs. The remaining few novae lie in the Drimmel et al. (2003) region, where the analytic model extends out to a galactocentric radius of R = 15 kpc, or the entire size of the grid. So, a large fraction of the novae in our model could have underestimated extinction values, but we suspect that this is only for the severely extinguished novae, and thus they are already unobservable in the optical but not necessarily in the IR. For example, in the inner Galactic plane (∣l∣ < 30° and ∣b∣ < 2°) out to 7 kpc, the mwdust map predicts, on average, A_g ∼ 10 mag of extinction in the g band but only A_J ∼ 2 mag of extinction in the J band. So, a nova with a luminosity of M = −7.2 mag would peak at m_g = 17 mag in the ASAS-SN g-band filter but as bright as m_J = 9 mag in the PGIR J-band filter. Therefore, we do not expect that our prediction that only 46%–53% of novae in the Galaxy get brighter than g = 15 mag would change if the extinction model was complete for the entire Galaxy, although it would change predictions for IR surveys and observations carried out in redder bands.

2.5. Reddened Novae

Figure 5 shows the positions of our model novae in Galactic coordinates around the Galactic center, distinguished by whether the peak apparent magnitude reached g = 15 mag. As expected, almost all of the heavily obscured, and therefore faint, novae lie within a couple degrees of the Galactic plane. This implies that optical observations will struggle to discover novae in regions within a couple degrees of the plane and especially toward the Galactic center.

fdg — **Figure 5.** Galactic coordinate positions of N = 1000 simulated novae around the Galactic center. The blue points indicate novae that reach a peak brightness of g = 15 mag or brighter and would likely be discovered by optical observations. The red points indicate novae that never reach g = 15 mag and have a much lower chance to be discovered by optical observations. The amount of extinction integrated out to 15 kpc from the `mwdust` model is shown as a gray-scale color map. The resolution of the map is 025 with a maximum extinction value of 28 mag for visualization purposes.
Download figure:
Standard image High-resolution image

To explore how our model predictions compare to the known sample of optically discovered novae, we have compiled a list of known novae by combining the sources from the CBAT list of novae in the Milky Way¹² and Koji's List of Recent Galactic Novae.¹³ The CBAT list consists of objects from 1612 to 2010, and we only include objects with eruptions after 1900. A literature search was then performed on the entire list to investigate if any contaminating sources were present and how many objects have been spectroscopically confirmed as novae. We found that 10 objects from the CBAT list are not classical novae and removed them from our list. Of the objects, 351 have spectroscopic or photometric observations suggesting that they are indeed classical novae, but we found no information about the remaining 47 objects. We assume that these objects are classical novae, but the possibility of contamination by other sources still remains.

Figure 6 shows the positions of optically discovered novae from our list in Galactic coordinates toward the Galactic center. Consistent with predictions from our model, there appears to be a significant lack of novae near the Galactic plane where most of the obscuring dust resides (∣b∣ ≲ 2°). However, there also appears to be a bias against discovery of novae at lower decl. This is likely due to a historic lack of observations in the Southern Hemisphere, although this issue has been recently addressed by ASAS-SN's Southern Hemisphere facilities and increasing numbers of amateur observers in Australia and Brazil (e.g., the Brazilian Transient Search).

Also shown in Figure 6 are nova candidates discovered by the Vista Variables in the Via Lactea (VVV; Minniti et al. 2010) and the Optical Gravitational Lensing Experiment (OGLE; Udalski et al. 2015). The VVV deep near-IR observations and the OGLE I-band observations of the Galactic bulge and nearby disk are better suited than most optical observations to discover novae in dustier fields at low Galactic latitudes, although the cadence of the OGLE observations is much lower in high-extinction regions at low latitude (Udalski et al. 2015). The 20 candidates from VVV (Saito et al. 2012, 2013; Beamin et al. 2013; Saito et al. 2014; Montenegro et al. 2015; Saito et al. 2015; Contreras Pena et al. 2016; Gutierrez et al. 2016; Saito et al. 2016, 2017) and 19 from OGLE (Kozlowski et al. 2012; Wyrzykowski et al. 2014; Mroz & Udalski 2014, 2016) were either discovered in the data after eruption or not followed up spectroscopically. Many are likely classical nova eruptions, though the sample could be contaminated by a few dwarf novae and young stellar objects. As seen in Figure 6, it does appear to be the case that there are more VVV and OGLE nova candidates closer to the plane, but there are still regions where few to no novae or nova candidates have been discovered.

For example, there has never been a nova or nova candidate discovered in the 20 deg² patch of sky with a Galactic longitude and latitude of −10° < l < 0° and −1° < b < 1°, respectively. Our model predicts that ∼10% of Galactic novae should be in this region, but the average extinction for this region is A_g = 23, according to the Marshall et al. (2006) dust map. This is a region that OGLE observed less frequently than other bulge fields (Udalski et al. 2015), but ASAS-SN has observed this region at a high cadence.

The absence of novae at low Galactic latitude is perhaps shown more clearly in Figure 7, which compares the Galactic latitude distribution of simulated novae with known optically discovered novae. Our model predicts that ∼60% of novae erupt within $\left|b\right|\lt 2^\circ ;$ however, only ∼20% of the optically discovered sample resides within this region. Also shown in Figure 7 is the distribution of a model magnitude-limited (g < 15 mag) sample of novae. This distribution peaks at a Galactic latitude of 2°, similar to the observed distribution, consistent with a historic magnitude-limited sample with dust as the determining factor.

3. Global Nova Rate Estimates

Next, we discuss how our nova model, described in Section 2, can be used to explore what fraction of simulated novae would be observable for various surveys and what that implies for the global Galactic nova rate.

The all-sky and nearly 1 day cadence observations of ASAS-SN provide an unprecedented opportunity to better constrain the Galactic nova rate. Even though the limiting magnitude of ASAS-SN is as deep as g ∼ 18 mag, nova searches have some unique challenges, and we believe the current infrastructure of the transient candidate pipeline is best suited to discover novae brighter than g ∼ 15 mag for various reasons. First, ASAS-SN usually flags portions of the sky known to contain variable stars and avoids searching for transients in these regions, but classical novae often have variable hosts. To address this issue, a special "ASAS-SN Nova Alert" email is generated for essentially all transients g ≲ 15 mag, regardless of previous variability, and immediately sent to alert several of the coauthors upon detection. Second, confusion from neighboring sources—many of which are variable—is a much larger issue in the Galactic plane and essentially translates to a shallower detection limit. Additionally, the candidate pipeline dedicated to discovering novae has a cutoff at this threshold. And finally, the number of contaminants greatly increases closer to the detection limit, so candidates fainter than g = 15 mag are not always checked. The ASAS-SN can, and does, find Galactic transients fainter than this threshold, but the detection efficiency likely falls off quickly at g > 15 mag. This detection threshold is bright, especially when compared to the number of fainter CV candidates and extragalactic supernova candidates discovered by the survey,¹⁴ and the detection efficiency of crowded fields at deeper thresholds will be explored in future work.

Currently, the detection efficiency of Galactic transients brighter than g ≈ 15 mag is not known for ASAS-SN. For extragalactic supernovae, Holoien et al. (2019) found that ASAS-SN is essentially complete down to m = 16.2 mag, but we do not expect the completeness to be as deep in the plane where almost all novae reside. To estimate this value, a fake transient recovery analysis performed on ASAS-SN data is needed, but it is beyond the scope of this work. Even with perfect recovery in observable fields, an optical transient survey is limited to detection rates ≲80% due to solar conjunction (Mróz et al. 2015), and preliminary estimates of detecting fainter but longer-lived Type Ia supernovae in ASAS-SN suggest detection capabilities between 70% and 80% (D. Desai et al. 2021, in preparation). A fake transient recovery analysis was performed on PGIR data in De et al. (2021), and they estimated that 36% of all Galactic novae that reach J = 14 mag in their field of view could be detected (δ > −28 fdg 9 at a cadence of ≈2 nights). This lower detection fraction is largely due to crowding/blending from an 8'' pixel scale, contamination from nearby bright stars, and the Galactic center being unobservable for a large fraction of the year from PGIR's Mt. Palomar location. We expect the first two issues to be present in ASAS-SN data, since the two surveys have the same pixel scale and a majority of novae should be found within a couple degrees of the crowded plane, but the last issue does not affect ASAS-SN because it has facilities in both the Northern and Southern Hemisphere. Taking all of the above information, we estimate that ASAS-SN can detect 60% of Galactic novae that reach g = 15 mag. To be clear, this estimated detection efficiency is only applicable to Galactic novae, and the detection efficiency is known to be higher for extragalactic transients in less crowded fields off of the plane.

Between 2018 and 2020, there were 31 known Galactic novae that peaked brighter than g = 15 mag. A majority of these were clearly detected and flagged as transients in ASAS-SN data, but there are at least a few examples of novae that were not detected or flagged as nova candidates. Object V3731 Oph (De et al. 2020b) was detected as a transient candidate in the ASAS-SN pipeline but was confused with a coincident variable within a pixel of the nova. Object V6567 Sgr (De et al. 2020c) was detected on the rise in ASAS-SN data and initially reported as a CV candidate. And though V1709 Sco (Kawash et al. 2020) was detected as bright as V = 12.7 mag by other observers, it was never detected brighter than g = 15 mag in ASAS-SN data, likely due to facilities being shut down for a large portion of 2020 due to the pandemic. It is likely that more observable nova events were missed by all transient surveys and observers and even more due to solar conjunction, so this is consistent with our 60% detection efficiency estimate for this time period. This estimation is very crude, and it will be one of the major goals of A. Kawash et al. (2022, in preparation) to better understand and constrain it.

We use our model to distribute N = 1000 novae in a mock galaxy, estimate what fraction would be detectable by ASAS-SN by assuming detections of 60% ± 6% that reach g < 15 mag, and extrapolate to a global rate from an annual discovery rate of R = 10.3 ± 1.9. This analysis is carried out for 1000 iterations, each time randomly distributing novae according to the stellar density model and sampling a normal distribution with a mean and standard deviation equal to the estimated value and uncertainty, respectively, for the discovery rate, detection efficiency, and luminosity function of novae in order to evaluate the most likely Galactic nova rate based on ASAS-SN observations.

Because this is the first nova rate estimate from ASAS-SN, and because the detection efficiency is not well constrained, we compare our results to those derived from other transient surveys. First, it was estimated that between 2010 and 2013, OGLE-IV observations discovered up to 80% of novae brighter than I = 17 mag in the most frequently visited fields in their field of view (−10° < l < 10° and −7° < b < 5°, with a cadence varying from 20 minutes to a few days; Mróz et al. 2015). The OGLE discovery rate was R = 4.8 ± 1.1 yr⁻¹ over this time period. We carry out the same analysis as we did for ASAS-SN to derive a Galactic nova rate from OGLE-IV observations using a detection efficiency as a function of sky position estimated from Figure 9 of Mróz et al. (2015), but we ignore novae directly in the plane ( $\left|b\right|\lt 1^\circ$ ) and fields with a detection efficiency less than 25% to account for the different extinction model used in our analysis. We believe this is a safe assumption, as only one nova from Table 1 of Mróz et al. (2015) lies within this ignored region. The OGLE-IV observations have a much lower cadence in the highest-extinction regions in the plane (Udalski et al. 2015), so a large fraction of novae are likely undetected despite the better pixel scale and redder filter compared to ASAS-SN. However, roughly 40%–50% of these novae are too highly extinguished to be detectable even with improved monitoring of the field.

Table 1. Galactic Nova Rates for Various Parameters

		Parameters			Implied Galactic Rate (yr⁻¹)
θ	Dust Model	Mass Model	N_d/N_b	ASAS-SN	OGLE	PGIR	Average ( $\bar{\chi }$ )
1.0	`mwdust`: ${A}_{V}/{A}_{{K}_{s}}=8.65$	Besançon	1.7	33 ± 7	36 ± 9	44 ± 15	35
1.0	`mwdust`: ${A}_{V}/{A}_{{K}_{s}}=8.65$	Cautun	5.0	38 ± 8	52 ± 13	45 ± 15	43
1.0	Exponential disk	Besançon	1.7	24 ± 5	36 ± 9	43 ± 14	28
0.4	`mwdust`: ${A}_{V}/{A}_{{K}_{s}}=8.65$	Besançon	0.7	31 ± 6	27 ± 6	43 ± 14	31
0.4	`mwdust`: ${A}_{V}/{A}_{{K}_{s}}=8.65$	Cautun	1.7	38 ± 8	41 ± 10	44 ± 15	40
1.0	`mwdust`: ${A}_{V}/{A}_{{K}_{s}}=13.44$	Besançon	1.7	38±9	36 ± 9	55 ± 19	40
1.0	`mwdust`: ${A}_{V}/{A}_{{K}_{s}}=13.44$	Cautun	5.0	44±10	54 ± 13	44 ± 15	47

Note. Model parameters and derived Galactic nova rates based on ASAS-SN, OGLE, and PGIR detections of novae. The default extinction curve for the mwdust model is ${A}_{V}/{A}_{{K}_{s}}=8.65$ , but the last two rows show results for a steeper value. Here θ is the ratio of disk to bulge novae per unit mass, N_d/N_b is the resulting disk-to-bulge ratio of novae for a given mass model, and $\bar{\chi }$ is the weighted average of the Galactic nova rate from the three surveys.

Download table as: ASCII Typeset image

Another survey with a published rate of nova discovery is the PGIR survey (De et al. 2021). Over the first 17 months of observations, they discovered 7.8 novae yr⁻¹ and estimated that they could detect 36% of all novae brighter than J ∼ 14 mag in their field of view (δ > −28 fdg 9 at a cadence of ≈2 nights). We run our analysis on a PGIR detection rate of r = 7.8 ± 2.3 yr⁻¹ and detection efficiency of 36% ± 3.6% in the field of view.

Lastly, Shafter (2017) derived a Galactic nova rate with a bright nearby sample of novae. Between 1900 and 2020, there were only seven novae that reached an apparent magnitude of m = 2. We carry out our analysis assuming that 90% of all novae that reached this brightness were discovered.

The distributions of nova rates derived from these various observational constraints are shown in Figure 8. All of these distributions except for the m < 2 mag constraint are well fit by a normal distribution, and for those, we derive the mean and standard deviation (listed in Table 1). The ASAS-SN (32.5 ± 7.3 yr⁻¹), OGLE (35.6 ± 8.8 yr⁻¹), and PGIR (44.2 ± 14.7 yr⁻¹) derived rates are all consistent at the 1σ level. The PGIR rate has a higher variance because of the limited field of view and because they are unable to observe the Galactic center for a large portion of the year. The bright nova constraint (m < 2 mag) results in a distribution with high variance due to a low observed rate (seven novae over 120 yr). The ASAS-SN derived distribution results in the distribution with the lowest variance, but this distribution is derived with only a rough estimate of the detection efficiency of novae. Once this value is better constrained and more observations are accumulated, ASAS-SN could provide the best constraint on the Galactic nova rate to date.

The observational constraints from these various surveys are almost completely independent. The OGLE rate is derived from OGLE-IV observations occurring between 2010 and 2013. The ASAS-SN discoveries occurred between 2018 and 2020, overlapping with the PGIR discoveries from 2019 July to 2020 November. There have been no m = 2 mag novae that have erupted since any of these surveys started observing. A simple weighted average of ASAS-SN, OGLE, and PGIR derived rates results in a Galactic nova rate of R ≈ 35 yr⁻¹.

4. How Sensitive Are the Results to Our Assumptions?

Here we explore how our results change as we vary certain assumptions in our model. In Section 4.1, we change the extinction model from the three-dimensional dust maps to the simple exponential used in Shafter (2017) and also explore implementing a different reddening law. Then, in Section 4.2, we see how our results depend on the mass model of the Galaxy. Finally, we briefly explore how assuming different populations of bulge and disk novae affects our results in Section 4.3. The derived rates for various sets of parameters are shown in Table 1.

4.1. Extinction Models

The dust maps of Green et al. (2019), Marshall et al. (2006), and Drimmel et al. (2003)—stitched together by Bovy et al. (2016)—form the best all-sky three-dimensional dust map to date. It is almost certainly superior to simply assuming a disk of dust (i.e., as in Shafter 2017), but it is only able to model extinction out to a few to ∼10 kpc, depending on the direction. For this reason, we investigate how the rates change if we implement the exponential disk model of extinction used in Shafter (2017). The results are shown in the top right panel of Figure 9.

**Figure 9.** Same as Figure 8 but for various combinations of stellar density models, extinction models, and ratio of disk to bulge novae per unit mass (θ). Top left: Galactic nova rate distributions resulting from using the Cautun et al. (2020) mass model. Top right: Galactic nova rate distributions resulting from changing the `mwdust` model to an exponential disk; this results in rate estimates inconsistent at the 1σ level. Bottom left: Galactic nova rate distributions from an elevated rate of nova production in the bulge using the Cautun et al. (2020) mass model. Bottom right: same as bottom left but for the Robin et al. (2003) model. Overall, the ASAS-SN and OGLE derived rates are sensitive to the model assumptions, but the PGIR rate is not.
Download figure:
Standard image High-resolution image

As expected, the rate estimated from the g-band observations of ASAS-SN is extremely sensitive to the dust model used. The exponential model of extinction underestimates the amount of dust in the plane relative to the mwdust model, therefore yielding a higher detection fraction and, ultimately, a much lower rate. The OGLE derived rate is not sensitive to the dust model, since we assume that they do not detect novae in the fields with the highest extinction, and the IR observations of PGIR are also not sensitive to the dust model. The rates derived from the three surveys are no longer consistent at the 1σ level when using this dust model, and it is clear that an under- or overestimation of Galactic extinction will cause a significant error in the derived nova rate from ASAS-SN observations. This conclusion is consistent with predictions for observing the next Galactic supernova. Adams et al. (2013) found that different dust models yield different likelihoods of observing a Galactic supernova in the optical but were less important for near-IR observations.

The mwdust model assumes an extinction ratio of ${A}_{V}/{A}_{{K}_{s}}=8.65$ across the entire sky, but this is not a safe assumption along all lines of sight, especially toward the inner Galaxy (Nataf et al. 2016). Therefore, we also investigate how our results change assuming an extinction ratio of ${A}_{V}/{A}_{{K}_{s}}=13.44$ , causing extinction in the g band to be 25% higher than the default mwdust value. As expected, the g-band observations of ASAS-SN are sensitive to variations in the extinction curve, with steeper reddening laws leading to higher inferred rates of nova production (see last two rows of Table 1). However, not all novae in the model are in the inner Galaxy; therefore, the extinction is likely overestimated for a fraction of the lines of sight. Hence, this alternate model should be considered closer to a limit on the effects of a different extinction law on the Galactic nova rate, rather than a best estimate of the effect.

4.2. Mass Model

The Robin et al. (2003) stellar density model is just one of many widely used Galactic mass models, so we explored whether our results change if we implement another model. A stellar density model inferred from Gaia DR2 data was recently presented in Cautun et al. (2020), and we implement the bulge, thin disk, and thick disk components from the contracted halo version of this model. The grid was calculated using the same resolution and boundaries as our primary model, and the mass of each component was found to be consistent with the derived mass in Table 2 of Cautun et al. (2020).

The Besançon model of the Galaxy has a more massive, bar-like bulge component, while Cautun et al. (2020) assumed an axisymmetric bulge for simplicity and stated that they are unable to constrain the bulge mass or its radial profile. However, by implementing the Cautun et al. (2020) stellar density model, we are able to see how the derived nova rate changes with a different assumed distribution of novae in the model galaxy. As seen in Figure 1, the Cautun et al. (2020) model places fewer novae at smaller Galactic radii but more novae at a lower Galactic height from the plane. So, the Cautun et al. (2020) model results in fewer bulge novae but more in the highest-extinction regions in the plane. As seen in Figure 9, this predicts a slightly higher rate from ASAS-SN observations and a significantly higher rate from OGLE observations. The PGIR rate does not appear to be sensitive to the stellar distribution model. Overall, using the Cautun et al. (2020) model results in a higher prediction of the Galactic nova rate, but it is still consistent at the 1σ level with using the Robin et al. (2003) model.

4.3. Differing Bulge and Disk Populations

It has been posited that novae that erupt in the bulge have different properties than those that erupt in the disk because of different progenitor populations hailing from differing star formation histories (Della Valle & Izzo 2020). Darnley et al. (2006) found that the favored model of novae in M31 supported separate disk and bulge populations that erupted at different rates per unit r-band flux. They found that, per unit r-band flux, the ratio of disk novae to bulge novae was 0.18. Shafter & Irby (2001) also studied the spatial distribution of novae in M31 and estimated this ratio to be 0.4. In a similar fashion, we define θ as the ratio of disk novae to bulge novae per unit mass in our model. So far, we have assumed one population of novae that traced the overall stellar mass of the Galaxy (θ = 1), resulting in a ratio of disk-to-bulge novae of N_disk/N_bulge ≈ 5 for the Cautun et al. (2020) mass model and N_disk/N_bulge ≈ 1.7 for the Robin et al. (2003) mass model.

Does the the Milky Way produce more novae per unit mass in the bulge than in the disk? Because distances are often hard to constrain, this is not an easy question to answer, but to first order, the higher the nova rate in the bulge, the more novae we should expect to find near l = 0°. From our model, where novae simply trace the stellar mass of the Galaxy (θ = 1), we expect 40% of bright (g < 15 mag) novae to be located within $\left|l\right|\leqslant 10^\circ$ . Of all the known Galactic novae from our list, 45% have erupted within $\left|l\right|\leqslant 10^\circ$ . This could support the idea of bulge enhancement, or the bulge producing more novae per unit mass relative to the disk, especially since the observed sample of novae is likely biased toward nearby disk novae.

For this reason, we explore how our results changed when using an elevated bulge rate of θ = 0.4. This model predicts that 49% of bright (g < 15 mag) novae are within $\left|l\right|\leqslant 10^\circ$ based on the Cautun et al. (2020) mass model and 67% for the Robin et al. (2003) model. The value of θ is difficult to constrain due to the unknown number of foreground disk novae, but the number of known novae around the Galactic center suggests it is larger than θ = 0.4. The nova rate results for this elevated bulge distribution of novae are shown in the bottom panels of Figure 9 for both the Cautun et al. (2020) and Robin et al. (2003) mass models.

If the production of novae in the bulge is elevated relative to the disk per unit mass, the global rate based on OGLE-IV observations decreases. This is expected, as OGLE observations are heavily biased toward finding bulge novae. The ASAS-SN and PGIR observations are less sensitive to the ratio of disk to bulge novae, as the fields of view of these surveys are less biased toward the bulge or disk.

5. Conclusions

We have used an all-sky three-dimensional dust map to explore the effects of Galactic extinction on the discovered nova rate. This model predicts that roughly half of nova eruptions will be too faint for current ASAS-SN discovery abilities, much higher than previous estimates that used a simpler dust model, and likely explaining much of the discrepancy between observed and predicted rates. Many of the highly extinguished, reddened novae lie within 2° of the plane. Our model predicts that ∼60% of all novae should erupt within 2° of the plane, compared to only ∼20% of the discovered sample being found in this region. This further highlights the necessity of optical surveys to observe in redder bands and IR transient surveys like PGIR to detect these highly reddened novae in the plane, although these fields have lower recovery rates for surveys with large pixel scales.

For the first time, we have estimated a Galactic nova rate based on an all-sky survey with nightly cadence. The ASAS-SN observations between 2018 and 2020 suggest a Galactic nova rate of 33 ± 7 yr⁻¹. This derived rate relies heavily on the detection efficiency (assumed to be 60% in this work for novae brighter than g < 15 mag), a value that will need to be better understood in the future. However, the derived rate from ASAS-SN is consistent with rates derived from OGLE (R = 36 ± 9 yr⁻¹) and PGIR (R = 44 ± 15 yr⁻¹) observations. Our results are consistent with the recent higher rates derived in De et al. (2021) and Shafter (2017), though we cannot rule out lower rates at this time due to the large uncertainties in the model parameters.

The derived rate from ASAS-SN's blue g-band filter is sensitive to the extinction model implemented, so a precise nova rate estimation will depend on how accurately dust is modeled close to the plane. Similarly, the OGLE rate is sensitive to the level of bulge enhancement and, along with the ASAS-SN rate, the model used to place novae within the mock Galaxy. The PGIR rate does not appear to be sensitive to altering any of the assumptions of our model. For any combination of stellar density model, extinction model, and level of bulge enhancement, the observations of ASAS-SN, OGLE, and PGIR suggest a Galactic nova rate of ∼30–40 yr^–1.

Overall, this work makes significant progress in constraining the Galactic nova rate, but it can still be greatly improved. In A. Kawash et al. (2022, in preparation), we plan to estimate the detection efficiency of ASAS-SN through fake transient recovery and incorporating various decline rates into our simulated sample of novae. Knowing this, along with continued observations from ASAS-SN and PGIR, will allow us to further constrain the rate of novae in the Galaxy. We can further quantify the effects of solar constraint on nova discovery rates and make predictions for next-generation transient facilities like LSST.

We thank P. Mroz for providing data contributing to this work. We also thank Jo Bovy, D. J. Marshall, A. C. Robin, and I. T. Simion for helpful discussions and the referee for suggestions that substantially improved the paper. A.K., L.C., E.A., and K.V.S. acknowledge the financial support of NSF award AST-1751874 and a Cottrell fellowship of the Research Corporation. J.S. acknowledges support from the Packard Foundation. B.J.S., C.S.K., and K.Z.S. are supported by NSF grant AST-1907570. C.S.K. and K.Z.S. are supported by NSF grant AST-181440.

We thank the Las Cumbres Observatory and its staff for their continuing support of the ASAS-SN project. The ASAS-SN is supported by the Gordon and Betty Moore Foundation through grant GBMF5490 to the Ohio State University and NSF grants AST-1515927 and AST-1908570. Development of ASAS-SN has been supported by NSF grant AST-0908816, the Mt. Cuba Astronomical Foundation, the Center for Cosmology and AstroParticle Physics at the Ohio State University, the Chinese Academy of Sciences South America Center for Astronomy (CAS-SACA), and the Villum Foundation.

The analysis for this work was performed primarily in ipython (Perez & Granger 2007) using numpy (Oliphant 2006; Van Der Walt et al. 2011), Astropy (Price-Whelan et al. 2018), Matplotlib (Hunter 2007), and scipy (Virtanen et al. 2020).

Appendix

A.1. Besançon Mass Model

Here we discuss the form of each component of the Besançon mass model used to distribute novae for our primary model. The total mass and normalization of each component are shown in Table A1.

Table A1. Mass and Normalization Values for the Various Components of the Galactic Model

Component	Normalization	Total Mass
	M_⊙ pc⁻³	10⁹ M_⊙
Thin disk	1.45	35.0
Thick disk	0.002^a	4.67
Bulge (model S)	2.37	22.1
Bulge (model E)	1.17	1.20
Halo	0.00005	0.55

Notes. Normalization values and the total mass of the various components of the Galactic model utilized in this work. We set the normalization values to achieve a total mass consistent with that derived in Robin et al. (2003), Simion et al. (2017), or Bland-Hawthorn & Gerhard (2016).

^aDensity at the solar position, where the other normalization values refer to the density at the Galactic center.

Download table as: ASCII Typeset image

We use a Cartesian grid with resolution 0.1 kpc, R = $\sqrt{{x}^{2}+{y}^{2}}$ is the radial distance from the Galactic center, and z is the distance perpendicular to the Galactic plane. The assumed solar radius from Robin et al. (2003) is R_⊙ = 8.5 kpc.

A.1.1. Thin Disk

The form of the thin disk density is from Robin et al. (2003),

$\begin{eqnarray}\begin{array}{rclrcl}\rho & = & {\rho }_{0} & & \times & \left\{\exp \left[-\left({0.5}^{2}+\displaystyle \frac{{a}^{2}}{{h}_{{R}_{+}}^{2}}\right)\right]-\exp \left[-\left({0.5}^{2}+\displaystyle \frac{{a}^{2}}{{h}_{{R}_{-}}^{2}}\right)\right]\right\},\end{array}\end{eqnarray} \tag{ A1 }$

where a² = R² + (z/)², ${h}_{{R}_{+}}$ = 2.5 kpc is the scale length of the disk, ${h}_{{R}_{+}}$ = 0.9 kpc is the scale length of the hole, and = 0.0791.

A.1.2. Thick Disk

A piecewise thick disk density distribution is utilized from Robin et al. (2003),

$\begin{eqnarray}\rho =\left\{\begin{array}{ll}{\rho }_{0}\exp \left(-\displaystyle \frac{R-{R}_{\odot }}{{h}_{R}}\right)\times \left(1-\displaystyle \frac{{z}^{2}/{h}_{z}}{\xi \times \left(2+\xi /{h}_{z}\right)}\right) & \mathrm{if}\,z\leqslant \xi \\ {\rho }_{0}\exp \left(-\displaystyle \frac{R-{R}_{\odot }}{{h}_{R}}\right)\times \exp \left(-\displaystyle \frac{\left|z-{z}_{\odot }\right|}{{h}_{z}}\right)\times \displaystyle \frac{2\exp \left(\xi /{h}_{z}\right)}{2+\xi /{h}_{z}} & \mathrm{if}\,z\gt \xi \end{array}\right.,\end{eqnarray} \tag{ A2 }$

where h_R = 2.5 kpc is the radial scale length, h_z = 0.8 kpc is the vertical scale length, and ξ = 0.4 kpc. The local density ρ₀ = 0.002 M_⊙ pc⁻³ is set to be 4% of the local thin disk density (Bland-Hawthorn & Gerhard 2016).

A.1.3. Bulge/Bar

For the bulge, we use an updated fit to VVV data from Simion et al. (2017). The best-fit model combines a hyperbolic secant density distribution

$\begin{eqnarray}&&\rho ={\rho }_{0}{{\rm{sech}} }^{2}\left({r}_{s}\right)\,(\mathrm{model}\,{\rm{S}})\end{eqnarray} \tag{ A3 }$

and an exponential distribution

$\begin{eqnarray}&&\rho ={\rho }_{0}\exp \left(-0.5{r}_{s}^{n}\right)(\mathrm{model}\,{\rm{E}}),\end{eqnarray} \tag{ A4 }$

where

$\begin{eqnarray}&&{r}_{s}={\left\{{\left[{\left(\displaystyle \frac{x}{{x}_{0}}\right)}^{{c}_{\perp }}+{\left(\displaystyle \frac{y}{{y}_{0}}\right)}^{{c}_{\perp }}\right]}^{\displaystyle \frac{{c}_{\parallel }}{{c}_{\perp }}}+{\left(\displaystyle \frac{z}{{z}_{0}}\right)}^{{c}_{\parallel }}\right\}}^{1/{c}_{\parallel }}.\end{eqnarray} \tag{ A5 }$

Here c_∥ and c_⊥ control the face-on and edge-on shape of the bulge, respectively, and x₀, y₀, and z₀ are the scale lengths in each respective direction. We use the best-fit parameters using the Besançon disks presented in Simion et al. (2017). For the ${\rm{sech}}$ component (model S), the best-fit parameters are c_∥ = 2.89, c_⊥ = 1.49, and (x₀, y₀, z₀) = (1.65, 0.71, 0.50) kpc, and for the exponential component, the best-fit parameters are c_∥ = 3.64, c_⊥ = 3.54, and (x₀, y₀, z₀) = (1.52, 0.24, 0.27) kpc, and n = 2.87. The elongated direction of the bar is offset from the Sun–Galactic center line by 21 fdg 1 and 2 fdg 1 for the hyperbolic secant and exponential component, respectively. The bulge density has a cutoff radius R_c = 6.96 kpc implemented by multiplying the bulge density ρ by the function

$\begin{eqnarray}\begin{array}{ll}f(R)=1 & R\lt {R}_{c}\\ f(R)=\exp \left[-2{\left(R-{R}_{c}\right)}^{2}\right] & R\gt {R}_{c}.\end{array}\end{eqnarray} \tag{ A6 }$

A.1.4. Stellar Halo

We use a power-law form of the halo similar to the one presented in Robin et al. (2003),

$\begin{eqnarray}\rho =\left\{\begin{array}{ll}{\rho }_{0}{\left({a}_{c}/{R}_{\odot }\right)}^{n} & a\lt {a}_{c}\\ {\rho }_{0}{\left(a/{R}_{\odot }\right)}^{n} & a\gt {a}_{c}\end{array}\right.,\end{eqnarray} \tag{ A7 }$

where $a=\sqrt{{x}^{2}+{y}^{2}+{\left(z/\epsilon \right)}^{2}}$ , a_c = 0.5 kpc is the cutoff radius, and n = −2.44.

Galactic Extinction: How Many Novae Does It Hide and How Does It Affect the Galactic Nova Rate?

Article metrics

Permissions

Author e-mails

Author affiliations

ORCID iDs

Dates

Abstract

1. Introduction