CLEAR: Paschen-$\beta$ Star Formation Rates and Dust Attenuation of Low Redshift Galaxies

We use \Pab\ (1282~nm) observations from the Hubble Space Telescope ($\HST$) G141 grism to study the star-formation and dust attenuation properties of a sample of 29 low-redshift ($z<0.287$) galaxies in the CANDELS Ly$\alpha$ Emission at Reionization (CLEAR) survey. We first compare the nebular attenuation from $\Pab/\Ha$ with the stellar attenuation inferred from the spectral energy distribution, finding that the galaxies in our sample are consistent with an average ratio of the continuum attenuation to the nebular gas of 0.44, but with a large amount of excess scatter beyond the observational uncertainties. Much of this scatter is linked to a large variation between the nebular dust attenuation as measured by (space-based) $\Pab$ to (ground-based) $\Ha$ to that from (ground-based) $\Ha/\Hb$. This implies there are important differences between attenuation measured from grism-based / wide-aperture $\Pab$ fluxes and the ground-based / slit-measured Balmer decrement. We next compare star-formation rates (SFRs) from $\Pab$ to those from dust-corrected UV. We perform a survival analysis to infer a census of \Pab\ emission implied by both detections and non-detections. We find evidence that galaxies with lower stellar mass have more scatter in their ratio of \Pab\ to attenuation-corrected UV SFRs. When considering our \Pab\ detection limits, this observation supports the idea that lower mass galaxies experience"burstier"star-formation histories. Together, these results show that \Pab\ is a valuable tracer of a galaxy's SFR, probing different timescales of star-formation and potentially revealing star-formation that is otherwise missed by UV and optical tracers.


INTRODUCTION
Star-formation rates (SFRs) are a critical quantity in the understanding of galaxy evolution. There exist several different methods of estimating SFR in a given galaxy, the most direct of which is counting identifiable stars of a specific age (Kennicutt & Evans 2012). With current instrumentation, this method of star counting is limited to the most immediate of Milky Way satellites. In more distant galaxies, the primary methods of measuring SFRs are using continuum and emission-line tracers (e.g., Kennicutt & Evans 2012).
Near-ultraviolet (UV) continuum observations of a galaxy measure the photospheric emission of massive young stars, and so the UV continuum acts as a direct tracer of recent starformation, timescales of order hundreds of Myr (Kennicutt & Evans 2012;Reddy et al. 2012). However, UV continuum observations are highly sensitive to attenuation by interstellar dust. In principle, this attenuation can be corrected using the UV slope β, but in practice the unknown intrinsic UV slope and differences in the UV shape of different attenuation laws complicate this approach (e.g., Salim & Narayanan 2020). Another approach is to add the reprocessed IR emission to the observed UV continuum for a "ladder" SFR (Wuyts et al. 2011;Whitaker et al. 2014), but this is similarly complicated by UV optical depth effects and/or potential anisotropy of the IR emission (Kennicutt & Evans 2012;Barro et al. 2019).
Optical and near-infrared (IR) emission lines from ionized gas around massive stars are also widely used as SFR indicators. These emission lines receive peak contribution from stars of mass 30-40 M , and as such are tracers of stars with lifetimes of 3-10 Myr. Recombination lines of hydrogen are especially useful to trace star-formation since they are insensitive to metallicity and relatively insensitive to gas temperature and density (Osterbrock 1989). Since the continuum and emission-line tracers correspond to star-formation on differ-2 ent timescales, their ratio can be used to measure the burstiness of the star-formation (e.g., Guo et al. 2016;Weisz et al. 2012). Still, optical emission lines are susceptible to dust attenuation. For example, Hα flux is reduced by a factor of ∼2 at a modest attenuation of A V = 1, and reduced by a factor of ∼10 in a dusty galaxy with A V = 3. The Balmer decrement (Hα/Hβ) can be used to correct for the attenuation, but this correction is inaccurate in regions of high optical depth to Balmer emission and can be inaccurate if the emission and/or attenuation scales are smaller than the spatial resolution (Kennicutt & Evans 2012). Uncertainties in correcting for dust attenuation fundamentally limit measurements of star-formation burstiness from the UV continuum and optical emission-line SFR tracers (Broussard et al. 2019).
Near-IR recombination lines of hydrogen offer a solution to the problem of dust attenuation in measuring SFR. Just like the more commonly used Balmer series, the Paschen lines of hydrogen are highly sensitive to the ionizing (E > 13.6 eV) radiation of OB stars formed within the last 10 Myr, while remaining relatively insensitive to nuisance parameters like the temperature and density of the star-forming gas (Osterbrock 1989). Unlike the optical Balmer lines, the near-IR Paschen lines are far less affected by interstellar dust extinction, so the Paschen lines can reveal otherwise hidden star-forming regions that are shrouded in gas and dust that is optically thick to Balmer emission.
In previous work, Paβ and Hα have been studied in a 2 galaxy sample (Kessler et al. 2020). Previous studies have also used the Paα (18750Å) emission line to calibrate mid-IR SFR indicators in nearby starburst and luminous IR galaxies (Alonso-Herrero et al. 2006;Calzetti et al. 2007) and in rare lensed galaxies at higher redshift (Papovich et al. 2009;Finkelstein et al. 2011;Shipley et al. 2016).
In this work we study Paβ (1282 nm), the n = 5 → 3 hydrogen recombination line, as an SFR indicator. We use Paβ fluxes measured from near-IR spectroscopy from the HST/WFC3 grisms taken as part of the 3D-HST (Momcheva et al. 2016) and CLEAR (Simons et al. 2020) surveys, as described in Section 2. Section 3 discusses the viability of Paβ/Hα as an attenuation indicator compared to the Balmer decrement and V-band continuum attenuation. In Section 4, we compare Paβ to other SFR indicators, demonstrating that Paβ includes star-formation missed by UV and optical tracers and showing evidence for burstier star formation at low mass. We summarize our results and discuss future applications with the James Webb Space Telescope (JWST) in Section 5.

CLEAR Parent Sample G102 and G141 Spectroscopy, Redshifts and Line Fluxes
The grizli (grism redshift and line analysis) pipeline (Brammer 2019) serves as the primary method of data reduction for the CLEAR dataset. In contrast to traditional methods of extracting one-dimensional (1D) spectra from slit observations, grizli directly fits the two-dimensional (2D) spectra with model spectra convolved to the galaxy image and for multiple position angles of grism observations. This process yields complete and uniform characterization of the suite of spectral line features of all objects observed in each of the G102 and G141 grisms. The most relevant of these spectral properties for our analysis are redshifts, line fluxes, and emission-line maps. The Paβ line is not included in the grizli fits by default, but was included in the CLEAR reductions for this work.
Our parent sample represents all CLEAR galaxies within the redshift range for detectable Paβ in the G141 spectrum (z < 0.287). The CLEAR spectral extractions are limited to galaxies with m F105W < 25.

Sample Selection
We select a sample of Paβ-emitting galaxies from the CLEAR parent catalog using the following steps: • Require z < 0.287, such that Paβ is within the observed-frame spectral range and blueward of the G141 sensitivity decline at 1.65 µm.
-Primary sample: A Paβ signal-to-noise ratio (SNR) of SNR > 3. (20 objects) -Secondary sample: A marginal Paβ SNR > 1 and a reliable spectroscopic redshift from either ground-based optical spectroscopy or from Hα emission in the G102 spectrum. (9 objects) The primary sample ensures reliable Paβ-detections, and the secondary sample ensures reliable redshifts from other (brighter) emission lines, while remaining as inclusive as possible to non-spurious Paβ emission. The combined sample is constructed to include all non-spurious Paβ-detections, even when the Paβ SNR is only marginal. For clarity in all of the following figures, objects in our primary sample will be plotted with larger symbol sizes than those in our secondary sample.
Our total (primary+secondary) sample includes 29 Paβemitting galaxies: 20 from the primary sample and 9 from the secondary sample. This total sample comprises approximately 19% of all CLEAR galaxies in this redshift range. The median Paβ SNR of the total sample is 3.9 with the median SNR of the >3σ primary sample detections of 5.1. The grism redshifts and Paβ line fluxes for the objects in our sample are taken from the CLEAR release v2.1.0 (Simons et al. 2020). We note that 17 of the objects in our sample have matching spectroscopic redshifts from ground-based programs (as compiled in the 3D-HST catalog), and another 6 have a redshift from well-detected Hα in the G102 grism. Figure 1 shows the color-mass relation for the galaxies in our Paβ-detected sample (in red) with respect to all CLEAR galaxies in the redshift range z < 0.287 (in gray). The F435W and F775W magnitudes and stellar masses are taken from the CANDELS/SHARDS multiwavelength catalog (Barro et al. 2019). Our galaxy sample is broadly representative of the larger galaxy population in CLEAR, with some preference for blue galaxies with log(M * /M ) 8. Figure 2 shows the redshift distribution of the galaxies in our sample. The redshift range of objects in our sample is 0.11 < z < 0.28 with a median redshift z = 0.23. Figure 3 shows RGB images (i, Y , and H band) of all 29 Paβ-detected galaxies in our primary and secondary samples, binned by stellar mass. Our sample includes a broad range of galaxy morphologies, including high-mass extended disks, bright compact sources, and diffuse irregular objects. Figure 4 shows the stacked one dimensional spectrum for all objects in the sample with Hα SNR > 2. We also show a suite of several other emission features in the near-IR, including relatively strong emission from doubly ionized sulfur (S III). The stacked emission features are visibly consistent with the intrinsic Paβ/Hα ratio of 1/17.6. We note that the stacked emission-line fluxes and flux ratios do not affect the main conclusions of the paper (see Figure 15), which rely on the scatter of measurements between individual galaxies. Figure 5 shows one-dimensional (1D) and twodimensional (2D) spectra for a galaxy in our Paβ-detected sample. We show one (GS 27549) object with strong enough Paβ SNR to extract spatially-resolved attenuation maps. The points of the 1D spectra show only the median points of the each individual grism exposures for this object. We note that the intrinsic Paβ/Hα ratio is a relatively weak 1/17.6. The spatially-resolved line maps of Hα and Paβ show nebular Figure 1. The relation between F435W-F775W color and stellar mass for galaxies of redshift z < 0.287. Our sample is shown as red stars, with the rest of the CLEAR galaxies in this redshift range shown as gray points. Larger symbols represent objects in our primary sample with Paβ SNR>3, and smaller symbols represent objects in our secondary sample with reliable redshifts and marginal Paβ 1<SNR<3. The galaxies in our sample are broadly representative of the population of z < 0.287 star-forming galaxies with log(M * /M ) 8. The sample of Paβ-detected galaxies also includes a few red galaxies that are likely dust-obscured (see Section 3). emission visually consistent with the intrinsic Paβ/Hα ratio (Osterbrock 1989).
Our CLEAR Paβ-detected sample may include a bias to high Paβ fluxes that are not representative of the full galaxy population, especially considering the line flux limit of 1.5 × 10 −17 erg s −1 cm −2 for the G141 grism observations (Momcheva et al. 2016). We use the Paβ flux limits, given as the line flux uncertainty from grizli when it tries to fit a Paβ line, for the rest of the 152 CLEAR galaxies in the same z < 0.287 redshift range. We consider this sample of galaxies with Paβ flux limits in a survival analysis in Section 4.2 that considers the relationship between Paβ and UV SFRs for both Paβ detections and non-detections.

Photometry and Derived Quantities
We take stellar masses for objects in our sample from the 3D-HST catalog (Skelton et al. 2014), derived from the CANDELS photometry (Grogin et al. 2011;Koekemoer et al. 2011). The stellar masses are calculated with FAST (Kriek et al. 2009), using a Bruzual & Charlot (2003) stellar population synthesis model library, a Chabrier (2003) IMF, solar metallicity, and assuming exponentially declining starformation histories. The stellar masses of our z < 0.287 galaxies are generally robust to these assumptions because Figure 2. Histogram of our sample of Paβ-selected objects binned by grism redshift. The G141 grism wavelength range limits the detection of Paβ to z < 0.287. Our sample has a redshift range of 0.11 < z < 0.28, with a median grism redshift of 0.23. the peak of the stellar emission is well-constrained by the high-quality CANDELS near-IR imaging. We additionally use the V -band attenuation (A V ) measured from the same FAST fit to the spectral energy distribution.
We use UV continuum SFRs from the catalog of Barro et al. (2019), which supplements the CANDELS multiwavelength data with SHARDS photometry (Pérez-González et al. 2013) in GOODS-N. Attenuation-corrected UV SFRs are calculated using the Kennicutt (1998)  (1) L 280 and A 280 are the UV luminosity and dust attenuation at rest-frame λ = 280 nm, respectively. The UV luminosity L 280 ≡ νL ν (280 nm) is calculated from EAZY with a best-fit spectral energy distribution (Brammer et al. 2008;Wuyts et al. 2011). The UV attenuation is inferred iteratively, measured from the best-fit SED while ensuring consistency with the IR (non)detection and the star-formation mass sequence (see Appendix D of Barro et al. 2019). The shortest wavelength filter in the Barro et al. (2019) catalog is the U band covering observed-frame λ ≈ 320 nm at its bluest end. This iterative approach is designed to produce attenuationcorrected UV SFRs that are robust to poorly measured UV photometry, which is especially useful given the limited restframe UV coverage of the low-redshift galaxies in our sample.
The conversion factor for the UV luminosity is derived in Bell et al. (2005). The attenuation-corrected UV SFR cali-brations are metallicity-dependent, with a systematic uncertainty of 0.05 dex from Solar to 20% Solar based on the Bruzual & Charlot (2003) models.
The peak timescale probed by the near-UV-derived SFRs at 2800Å is of order hundreds of Myr. (Reddy et al. 2012). SFRs derived from bluer luminosities probe slightly shorter timescales, of order 10-100 Myr for SFRs from L 1500 (Reddy et al. 2012).
UV + IR "ladder" SFRs are calculated for objects with mid/far-IR detections following Wuyts et al. (2011): The relative scale of the UV and IR contribution is based on local universe calibrations (Kennicutt & Evans 2012), and the overall scale assumes a Chabrier (2003) initial mass function. The UV + IR ladder SFR ultimately measures the UV continuum emission that is not attenuated by dust plus the reprocessed IR continuum emission from the UV which is absorbed by the dust. For objects which do not have IR detections, the ladder SFRs are defined to be equal to the attenuation-corrected UV SFRs of Eq 1 ( Barro et al. 2019).
In Figure 6, we compare the attenuation-corrected UV SFRs and the UV + IR ladder SFRs for CLEAR galaxies with mid-IR detections. For the CLEAR sample, the attenuationcorrected UV and UV+IR SFR indicators agree with each other with a median absolute deviation of ∼0.09, similar to the scatter reported in Barro et al. (2019). This indicates that the attenuation-corrected UV SFRs are likely to be reliable for our sample of galaxies.
Our galaxies have morphology measurements from van der Wel et al. (2012). We use effective (50% light) radii and Sérsic (1968) indices for galaxies with "good" GALFIT (Peng et al. 2010) fits with flag = 0 (see van der Wel et al. 2012 for details). From these effective radii and Sérsic indices, we calculate the central density within 1 kpc, Σ 1kpc .
To correct for dust attenuation of Paβ, Hα, and Hβ, we assume a (Calzetti et al. 2000) attenuation model. We use Calzetti et al. (2000) over other attenuation models of the Milky Way or the Small Magellanic Cloud (Fitzpatrick 1999;Gordon et al. 2003) to maintain consistency with the attenuation-corrected UV SFRs from Barro et al. (2019). The choice of attenuation model has little impact for this work since the attenuation models are very similar at optical and near-IR wavelengths (i.e., for the Balmer and Paschen lines). We use the Calzetti et al. (2000) attenuation model in accordance with previous studies showing a nebular-to-stellar attenuation ratio of ≈ 2 (Salim & Narayanan 2020;Calzetti et al. 2000).

Optical Spectra
A subsample of 11 galaxies in GOODS-N match to publicly available optical spectra from the Team Keck Treasury Redshift Survey (TKRS, Wirth et al. 2004), from which . Stacked spectra of objects in the primary and secondary sample with Hα signal to noise ratio greater than 2. We also show several other emission-line features visible in the G102 (blue) and G141 (red) spectral coverage. We observe stacked fluxes corresponding to the relatively faint intrinsic ratio of Hα/Paβ = 17.6 assuming Case B recombination as described in Section 1. We show an example spectrum for a single object in Figure 5.
we use Hα and Hβ fluxes. The TKRS spectroscopic observations of GOODS-N were taken using DEIMOS on the Keck II telescope, with the spectra extracted using the DEEP2 Redshift Survey Team pipeline (Newman et al. 2013). Disk-integrated Balmer line fluxes are estimated from TKRS spectra in a way that accounts for slit losses under the assumption that emission line equivalent widths are invariant across the stellar disk; see Section 2 of Weiner et al. (2007) for details. We note that this assumption of invariant equivalent widths makes the TKRS measurements potentially susceptible to issues when measuring Balmer-line fluxes for large objects with significant color gradients. We discuss this further in Section 3.1.

Linear Regression Methods and Significance of Fits
Throughout this analysis, we use the linmix package in python (Kelly 2007) to calculate our linear regression fits. For the remainder of this work, we consider a correlation between two quantities to be significant if the linmix mean best fit slope is 3σ different from zero. The standard deviation σ is the standard deviation of the linmix best-fit slopes.

Paβ AND DUST ATTENUATION
Because the ratios of the fluxes of recombination lines of hydrogen are relatively insensitive to metallicity, temperature, and density, their ratios can be used to estimate dust attenuation. The most commonly used emission-line indicator of dust attenuation is the Balmer decrement, Hα/Hβ, which has an intrinsic ratio of 2.86 assuming Case B recombination with T = 10 4 K and n e = 10 4 cm −3 (Osterbrock 1989). This rest-frame optical ratio only works for modestly attenuated Figure 5. Top left: Observed-frame one-dimensional (1D) spectrum for a galaxy in our sample. The G102 (blue) and G141 (red) spectra show the median points from all exposures for this object. The inset shows the region around the Paβ line. The gray shaded regions show the Paβ and Hα lines. Hα is available in the G102 in a small redshift window where Paβ is simultaneously available in the G141 (0.22 < z < 0.287). Bottom left: Observed-frame two-dimensional (2D) spectra for the same galaxy in our sample. We indicate both Paβ and Hα (where available) by annotated regions outlined in black lines. Right: Spatially-Resolved emission-line maps of Hα (top) and Paβ (bottom) for the same object from the grizli extractions. Paβ is a relatively weak line, with the intrinsic ratio of Hα/Paβ ≈ 17.6 (Osterbrock 1989). galaxies, since at A V ∼ 1 − 2 Balmer decrement attenuation measurements will entirely miss regions of the ISM that are optically thick to Hβ emission.

Paβ and Nebular Attenuation Indicators
We compare the Paβ and Hα fluxes and SFRs to investigate if the near-IR Paβ emission line reveals star formation that is otherwise hidden in optical emission. Figure 7 shows the Paβ and Hα fluxes for galaxies in our sample that have optical spectroscopy from TKRS in GOODS-N (Wirth et al. 2004). Open symbols show the observed fluxes, while filled symbols show the attenuation-corrected fluxes using the Balmer decrement and a Calzetti et al. 2000 attenuation model (assuming an intrinsic Case B Balmer decrement of Hα/Hβ = 2.86). The Paβ dust corrections are generally a factor of two or less, while the Hα dust corrections are often a factor of several or more. Even after correcting for dust attenuation according to the Balmer decrement, many of the galaxies have dust-corrected Paβ fluxes that are significantly greater than the expected ratio for Hα/Paβ (17.6/1, Osterbrock 1989). This suggests that the Balmer decrement is likely to underestimate the dust attenuation affecting Hα in many of our galaxies. In regions of high optical depth to Hβ in particular, we have highly uncertain attenuation correc-tions, which leads to highly uncertain attenuation-corrected Hα fluxes. Figure 8 shows the log ratio of the Paβ and Hα SFRs. As in Figure 7, open symbols show the uncorrected SFRs and filled symbols are dust-corrected using the Balmer decrement and a Calzetti et al. (2000) attenuation model. Linear regression suggests there are no significant correlations between the Paβ/Hα ratio and the stellar mass or observed Balmer decrement, with regression slopes consistent with zero with a constant offset of ∼0.3 dex. High ratios of the attenuationcorrected Paβ/Hα are likely to occur if the attenuation is underestimated, i.e. if the Balmer decrement does not measure all of the attenuation. Figure 8 suggests that, at least within our small sample, Balmer decrements underestimate the attenuation in galaxies spanning a broad range of stellar mass and Balmer decrement.
There is a potential selection effect that could bias our sample toward dustier galaxies because we have selected galaxies based on their Paβ flux. One potential remedy to this is to study the SFRs of the Hα non-detections (similar to the analysis in Figure 15), which would select less dusty galaxies than our Paβ-selected sample where the Balmer decrement is more reliably measured. We propose to augment this work with future studies of these Paβ/Hα ratios with much larger Figure 6. Left: The relation between attenuation-corrected UV continuum SFRs and UV+IR "ladder" SFRs for IR-detected galaxies in our CLEAR sample (large hexagons) and the Barro et al. (2019) CANDELS galaxies (small circles) for the z < 0.287 regime. Right: The relation between IR SFRs and UV+IR "ladder" SFRs for the same objects. We exclude the IR non-detections because these objects are defined to have equivalent attenuation-corrected UV and "ladder" SFRs. The attenuation-corrected UV SFRs more accurately models the true SFR than the IR SFR, which underestimates the total SFR by missing star-formation only visible in the UV.
datasets, which will be much better suited to perform the statistical analyses necessary to resolve this degeneracy.
We also draw attention the nuance of the different arguments presented in Figure 9 (where we compare the total SFRs derived from different means) and Figure 6 (where we compare the Pa-beta/H-alpha ratios against galaxy dust attenuation estimates). The continuum SFR indicators (e.g., the rest-UV and far-IR) presented in Figure 6 represent a longer timescale (the UV traces the light from B stars with ages of order ∼100 Myr, while the far-IR responds to light from longer lived A and even F stars with ages of order ∼1 Gyr; see e.g., Salim & Narayanan 2020) than the timescale probed by Hydrogen recombination lines (which are primarily produced from O stars with ages of ∼10 Myr). We show in Figure 6 that the attenuation-corrected UV SFRs agree with the UV+IR ladder SFRs, yet in Figure 9 we argue that the differences in SFR derived from the near-IR Paβ line and the optical Hα arise from attenuation. We argue these observations are not discrepancies but rather indicate variability in the timescale of star-formation as traced by the different SFR indicators. This leads us to conclude that the variation between the SFRs derived from the Hydrogen recombination lines (e.g., from Paβ) and those from the UV trace variability in the star-formation histories. We return to this point in Section 4.2. Figure 9 directly compares the observed Paβ/Hα and Hα/Hβ ratios for the 11 galaxies in our CLEAR sample with optical spectra from TKRS (Wirth et al. 2004). The lines indicate the expected ratios for Calzetti et al. (2000), SMC (Gordon et al. 2003), and Milky Way (Fitzpatrick 1999) at-tenuation models, using intrinsic Case B ratios of Hα/Hβ = 2.86 and Paβ/Hα = 1/17.6. Most (8/11) galaxies have line ratios that are broadly (within <3σ) consistent with the expectation from Calzetti et al. (2000), but three have significantly larger Paβ/Hα ratios than expected from their Hα/Hβ ratios. As in Figures 7 and 8, we interpret these galaxies as having Hα/Hβ ratios that underestimate the true attenuation. In the galaxy with the highest Hα/Hβ (upper right of Figure  9), this is somewhat expected since the Balmer decrement is likely to be inaccurate due to gas which is optically thick to Hβ emission. This object has a significant dust lane, visible in the images in Figure 3. The other two galaxies might similarly have Balmer decrements that probe only optically thin gas, with an additional ISM component that is optically thick to Hβ (and possibly Hα) but optically thin to Paβ emission.

Nebular Paβ/Hα and Continuum Attenuation
The left panel of Figure 9 compares the nebular attenuation measured by Paβ/Hα with the continuum V -band attenuation A V . We also show Calzetti et al. (2000), Milky Way (Fitzpatrick 1999), and SMC (Gordon et al. 2003) attenuation curves for a stellar to nebular attenuation ratio of 0.44, as well as the Calzetti et al. (2000) curve for equal stellar and nebular attenuation. All models assume an intrinsic Paβ/Hα = 1/17.6. Most (8/11) of the galaxies have observed Paβ/Hα ratios that are < 3σ consistent with the dotted line, albeit with significant excess scatter that suggests a large diversity of stellar to nebular attenuation ratios. Three of the galaxies have significantly larger Paβ/Hα ratios than expected for their A V , and these are the same three galaxies Figure 7. Paβ and Hα fluxes for the 11 galaxies in our sample with matching TKRS optical spectroscopy, color-coded by stellar mass (left) and Balmer decrement (right). The dashed gray line indicates Paβ/Hα = 1/17.6, appropriate for Case B recombination with T = 10 4 K and ne = 10 4 cm −3 (Osterbrock 1989). Open circles show uncorrected fluxes and filled circles are dust-corrected fluxes, calculated using the observed Balmer decrement and a Calzetti et al. (2000) attenuation curve. Larger symbols represent objects in our primary sample with Paβ SNR>3, and smaller symbols represent objects in our secondary sample with reliable redshifts and marginal Paβ 1<SNR<3. About half of the objects have attenuation-corrected ratios of Paβ/Hα which are significantly larger than the expected ratio, over a wide range of stellar mass and Balmer decrement. This picture is supported by a scenario where the small aperture size of the TKRS measurements and the assumptions of invariant equivalent widths leads to potentially under-measured Balmer line fluxes. Some objects have large uncertainties in attenuation-corrected Hα flux due to highly uncertain Hβ flux measurements. with larger Paβ/Hα than expected for their Balmer decrements. These galaxies are likely to have measured A V values which underestimate the attenuation, and/or have enshrouded star-forming regions that are optically thick to optical line and continuum emission but apparent in the near-IR Paβ line.

Issues with Balmer-line measurements from TKRS slit spectroscopy
The small size of our sample, with only 11 galaxies that have both rest-optical spectroscopy for Hβ and Hα along with rest-IR spectroscopy for Paβ, makes it unclear if our observations are representative of the nebular attenuation properties in the broader population of galaxies. Our Paβ-selected sample is also likely to over-represent high ratios of Paβ/Hα: starting from a Hα-or Hβ-selected sample would likely result in a different distribution of Paβ/Hα. However, we find evidence that at least some galaxies (spanning the range of stellar mass and observed Hα/Hβ) have Balmer decrements that underestimate the attenuation and miss star-formation that is otherwise revealed by Paβ.
Another potential issue with the direct comparison of Paβ emission from CLEAR and Balmer emission from TKRS is the potential for slit losses in the Keck observations. The TKRS disk-integrated Balmer line fluxes are estimated from TKRS spectra in a way that attempts to account for slit losses by assuming that the emission-line equivalent widths are invariant across the stellar disk; see Section 2 of Weiner et al. (2007).
This assumption likely fails for extended objects, evident from the objects with the largest Paβ/Hα discrepancies from the expected models in the right panel of Figure 9 also being some of the most extended objects in the sample (see Figure  3). This is important as the slit-width used by TKRS is 1 arcsecond, while the galaxy isophotes extend over many arcseconds. The assumption of invariant emission-line equivalent widths across extended galaxies with large color gradients can lead to significantly underestimated Balmer emission for galaxies with central dust lanes and/or higher emission-line equivalent widths in their outer regions.
In this section we have discussed the advantages of using comparisons of Paβ to emission-line tracers and continuum indicators to show dust attenuation missed by these other methods. We showed that the ratio Paβ/Hα is a valuable Larger symbols represent objects in our primary sample with Paβ SNR>3, and smaller symbols represent objects in our secondary sample with reliable redshifts and marginal Paβ 1<SNR<3. The Paβ dust corrections are generally a factor of two or less, while the Hα dust corrections are often a factor of several or more, which explains the observed ratios exceeding the attenuation-corrected ratios. Even after dust correcting according to the Balmer decrement, many of the galaxies have dust-corrected Paβ fluxes that are significantly greater than the expected ratio of Paβ/Hα = 1/17.6. This indicates that grism-based Paβ picks up star formation missed by slit-based optical emission-line SFR indicators. We fit the attenuation-corrected points and observe an excess of Paβ star-formation consistent with a constant offset of ∼ 0.31 dex. attenuation indicator in moderate to dusty galaxies. A more nuanced analysis of these issues and the benefits of a three emission-line attenuation model using the CLEAR sample is discussed in a companion work Prescott et al. (2022). This equation includes an attenuation correction for the Paβ emission, A Paβ , calculated from the measured continuum attenuation A V with a Calzetti et al. (2000) attenuation curve and a stellar-to-nebular attenuation ratio of 0.44 (Calzetti et al. 1994). In the following subsections, we compare Paβ SFRs measured from these equations with continuum and Hα SFR estimates.

Paβ and Continuum SFR Indicators
Here we compare Paβ SFRs with attenuation-corrected UV continuum SFRs and UV + IR 'ladder' SFRs to investigate hidden SFR and star-formation histories (SFHs). The near-IR Paβ line is much less affected by dust attenuation than the UV continuum, and so Paβ can reveal starforming regions that are otherwise obscured by dust in light of shorter wavelengths (Kennicutt & Evans 2012). In addition, Paβ (and other hydrogen emission lines) probes recent (<10 Myr) star-formation, while the UV continuum probes star-formation over longer (100-500 Myr) timescales (Kennicutt & Evans 2012; Reddy et al. 2012). The comparison of Paβ and UV continuum SFRs yields an indicator of SFH stochasticity ("burstiness", e.g., Guo et al. 2012;Broussard et al. 2019). Figure 10 shows the relation between the attenuationcorrected Paβ SFRs and the UV+IR ladder SFRs from Barro et al. (2019). We note that for 18 of the 29 objects in our sample, the Paβ SFR is >3σ higher than the UV+IR ladder SFRs. The object with the highest attenuation-corrected Paβ SFR is GN 34456, which has a significant dust lane (see Figure 3).  2000), Gordon et al. (2003), and Fitzpatrick (1999) attenuation curves with a stellar to nebular attenuation ratio of 0.44, along with another Calzetti et al. (2000) attenuation curve with a stellar to nebular attenuation ratio of 1. Larger symbols represent objects in our primary sample with Paβ SNR>3, and smaller symbols represent objects in our secondary sample with reliable redshifts and marginal Paβ 1<SNR<3. Right: The relation between Paβ/Hα and the Balmer decrement for the same 11 objects with Balmer line fluxes. The blue, green, and orange lines indicate the expected ratios using intrinsic Case B ratios of Hα/Hβ = 2.86 and Paβ/Hα = 1/17.6, and Calzetti et al. (2000) Gordon et al. (2003), andFitzpatrick (1999) attenuation models. Eight of the 11 points have line ratios within 3σ consistent with the expectation. Our sample includes at least one highly dusty galaxy (in the upper right) for which Hα/Hβ cannot reliably measure dust attenuation due to high optical depth to Balmer emission. Figure 11 shows the relation between Paβ luminosity and the attenuation-corrected UV SFR from the CAN-DELS/SHARDS multiwavelength catalog (Barro et al. 2019), with three panels color-coded by stellar mass and Vband attenuation. Twelve of the galaxies in our sample have 24 µm detections, and their UV + IR ladder SFRs are shown as open diamonds (Wuyts et al. 2011). Both SFRs correlate with stellar mass, as expected given the well-known starformation mass sequence of galaxies (Noeske et al. 2007;Whitaker et al. 2012), but there is more apparent scatter between the two SFR indicators in low-mass galaxies. The Paβ luminosity also tends to be greater than expected from the UV SFR in galaxies with steep UV slopes and high attenuation. Figure 12 shows the "excess" Paβ SFR compared to the attenuation-corrected UV SFR, quantified as log SFR Paβ SFR corr UV , with stellar mass and V -band attenuation A V . Our Paβ-selected sample is generally only sensitive to galaxies with Paβ SFR similar to or greater than the UV SFR. We fit the attenuationcorrected detections (colored points) in each panel using lin-ear regression, as implemented by the linmix (Kelly 2007) Python package. We find no significant correlations (with a slope >3σ different from zero) between the log SFR Paβ SFR corr UV ratio and mass or A V . The Paβ SFR "excess" is instead consistent with a constant offset of ∼ 0.3 dex over the attenuationcorrected UV SFR.
Some galaxies in our sample tend to have Paβ SFRs that are ∼1-2 orders of magnitude greater than the attenuationcorrected UV SFRs (see Figure 12). These galaxies may have dust-enshrouded star formation that is not seen in UV light (for example, behind high optical depths) and is not accounted for by the attenuation correction of the UV SFR. Paβ emission can escape these dusty star-forming regions that are optically thick to UV light, revealing star-formation that is hidden at shorter wavelengths.
There is a broad range of Paβ SFR excess in Figure 12 for galaxies of varying stellar masses and A V . However, our measured Paβ/UV SFRs may be influenced by the Paβ detection limits, which are mostly closer to a ratio of Paβ to UV SFRs of unity at low stellar mass. Our Paβ-selected study is gen- Figure 10. The relation between attenuation-corrected Paβ SFR and UV+IR 'ladder' star-formation rates for galaxies in our sample. Larger symbols represent objects in our primary sample with Paβ SNR>3, and smaller symbols represent objects in our secondary sample with reliable redshifts and marginal Paβ 1<SNR<3. The solid black line indicates the one-to-one relation. We note that the attenuation-corrected Paβ SFRs are at lease 3 σ greater than the ladder SFRs for 18 of the 29 objects in our sample. The object with the highest attenuation-corrected Paβ SFR is GN 34456, which has a significant dust lane (see Figure 3). erally only sensitive to bursty star-formation that occurred within the last 10 Myr, detectable as Paβ emission (with SFR Paβ > SFR UV ). In contrast, a burst of star-formation occurring 10-100 Myr ago would lead to SFR Paβ < SFR UV that is generally not detectable in our sample (as shown by the black upward facing triangles that indicate the Paβ flux limit for each galaxy in Figure 12). Thus we interpret the scatter of Paβ/UV SFRs in Figure 12 as consistent with bursty starformation in some galaxies. Furthermore,the Paβ detection limits are consistent with an undetected population of lowmass galaxies with SFR Paβ < SFR UV that represent bursty star-formation occurring >10 Myr ago, analogous to detected SFR Paβ > SFR UV galaxies representing bursts within the last 10 Myr. We discuss this further in Section 4.2. Figure 13 shows the Paβ luminosity and attenuationcorrected UV SFR color-coded by galaxy size, Sérsic index (measured by van der Wel et al. 2012), and central density Σ 1kpc calculated from these values and their stellar masses. Figure 14 shows the log ratio of the Paβ to UV SFR versus the same morphology quantities. There are no significant correlations between the ratio of the two SFRs with galaxy size, Sérsic index, or Σ 1kpc . Instead, we observe the same constant vertical offset of ∼ 0.3 dex. We conclude that galaxy morphology does not play a dominant role in the ratio between the Paβ and UV SFR of a galaxy. We discuss this further in the next subsection.
We note the average (median) attenuation-corrected Paβ "excess" (log SFR Paβ SFR corr UV ) ratio for our sample of 29 galaxies is 0.3. We show this by the horizontal black dash-dot line in Figures 12 and 14. We have also performed a stack of 44 objects in CLEAR with 0.2 < z < 0.3 (both detections and non-detections) and found a stacked Paβ luminosity of 2.5 ± 2.0 × 10 39 erg/s, which corresponds to a SFR Paβ ∼ 0.21 ± 0.16 M /yr. For this sample of 44 objects, the median attenuation-corrected UV continuum SFR is 0.44 M /yr, which is comparable to that of Paβ to within about 1σ. This stacking analysis is consistent with the conclusions of the following subsection: namely that the Pab and attenuation-corrected UV SFRs are broadly consistent (especially when considering non-detections), but with large scatter in the ratio of SFRs among individual galaxies.
Our sample of Paβ detections might be biased toward high Paβ/UV ratios, since fainter Paβ emission would be undetected and excluded from the sample. In the following section we perform a survival analysis of the Paβ/UV ratio with a sample of 152 Paβ non-detections and our 29 Paβ detections in the CLEAR survey to remedy these sample-selection biases.
Future studies with deeper grism detections of Paschenlines (see discussion of JWST surveys in the summery) will also add valuable information on analyses of this kind.

Burstier Star-Formation at Low Stellar Mass
The different star-formation timescales probed by hydrogen recombination lines (∼5 Myr) and near-UV continuum emission (∼100 Myr) means that their comparison can be used to indicate the burstiness of star-formation. Observations have long shown that the average ratio of Balmer-line to UV SFRs decreases at lower stellar mass (Sullivan et al. 2000;Boselli et al. 2009;Lee et al. 2009;Guo et al. 2016). The observations are best explained by an increasing importance of bursty star-formation occurring on timescales of tens of Myr, in which star-formation that occurred 10-100 Myr ago is detected in UV emission but not in hydrogen emission lines (Weisz et al. 2012;Hopkins et al. 2014;Shen et al. 2014;Sparre et al. 2017). But Broussard et al. (2019) notes that comparisons between Balmer-line and UV SFRs can be biased by uncertainties in the amount of dust attenuation (especially the UV shape of the attenuation law and the ratio of nebular to continuum attenuation).
Because the 29 Paβ-detected galaxies shown in Figure 12 are selected by their Paβ fluxes, the sample might be biased to high Paβ/UV ratios. To test this bias, we identify the sample of Paβ non-detections in CLEAR in the same redshift range and perform a survival analysis of the Paβ/UV ratio with this total sample of detections plus non-detections. We perform a survival analysis instead of a stacking analysis to Figure 11. The relation between Paβ luminosity and continuum star-formation rates for galaxies in our sample, color coded by mass (left), and continuum Av (right). Filled circles correspond to the attenuation-corrected UV SFRs and open diamonds correspond to UV + IR ladder SFRs for the twelve galaxies with well-detected IR emission (Barro et al. 2019). Larger symbols represent objects in our primary sample with Paβ SNR>3, and smaller symbols represent objects in our secondary sample with reliable redshifts and marginal Paβ 1<SNR<3. The solid black line indicates the relation between Paβ luminosity and SFR calculated using Equation 3 (Kennicutt & Evans 2012). The dashed gray lines represent the 1-sigma detection limits at redshifts z = 0.28 and z = 0.1. The SFR measured from the Paβ luminosity is higher than the attenuation-corrected UV SFR in galaxies with higher attenuation, and there is more apparent scatter between the two SFRs in galaxies with low stellar mass.
preserve the scatter and eliminate loss of information from binning (Feigelson & Babu 2012).
To test the possibility that low-mass galaxies are just as likely to have suppressed Paβ/UV as enhanced Paβ/UV, we perform a survival analysis of the relation of the Paβ to attenuation-corrected UV SFR ratio to stellar mass in Figure 15. We include all objects with Paβ SNR>1, and generate survival analysis expectations for the 152 Paβ non-detections in CLEAR. The inclusion of the non-detections will populate the low Paβ/UV regime, where our Paβ-selected sample is biased against. The survival analysis reveals the large undetected population of low-Paβ/UV galaxies with a distribution (especially for galaxies with less attenuation) that is broadly consistent with previous work (cf. Figure 4 of Weisz et al. 2012).
For each non-detection we randomly generate expectation values from a half-normal distribution with standard deviation equal to the 1σ upper bounds of each of the nondetections (open diamonds). The survival analysis does not have a significant dependence on the assumed prior distribution of the resampled points. The gray lines show rolling me-dian absolute deviations of the sample and survival analysis points. We perform the following analyses both for the full sample of detections and randomly sampled non-detections (left panel of Figure 15) and again for detections and nondetections for relatively unattenuated galaxies with A V < 1 (right panel of Figure 15).
We fit the scatter of the Paβ to UV SFR ratio in 6 bins of stellar mass for the Paβ detections and the survival analysis draws of the non-detections. We resample the measurements and the survival analysis draws 200 times and calculate the standard deviation of the distribution of SFR ratios measured from the resampled points in each bin. We then fit the standard deviations of SFR ratios across the six bins and find a best-fit slope of m = −0.063 ± 0.039 (left panel) and m = −0.151 ± 0.048 (right panel), indicating higher scatter at lower stellar mass with ≈ 1.6σ and ≈ 3.1σ significance, respectively.
The higher significance of the negative correlation between the scatter of the Paβ/UV ratio and stellar mass in the right panel of Figure 15 indicates that the scatter cannot be due solely to differences in attenuation among these galaxies. In- Figure 12. The log ratio of the Paβ and attenuation-corrected UV SFRs with stellar mass (left) and continuum Av (right). Solid circles represent attenuation-corrected Paβ/UV SFR ratios, and hollow circles correspond to the same objects without attenuation-corrected Paβ. Larger symbols represent objects in our primary sample with Paβ SNR>3, and smaller symbols represent objects in our secondary sample with reliable redshifts and marginal Paβ 1<SNR<3. We perform a linear regression fit to the attenuation-corrected Paβ/UV SFR ratios in each panel, finding no significant correlations between the Paβ/UV SFR ratios and stellar mass or AV . The black dash-dot line shows the median attenuation-corrected Paβ/UV ratio of 0.3. stead, we interpret this as evidence of stochasticity of the starformation histories. This result of burstier SFHs at lower stellar mass is consistent with previous works using Balmer/UV ratios (Weisz et al. 2012;Hopkins et al. 2014;Shen et al. 2014;Sparre et al. 2017).
The Paβ line should offer a clearer approach to measuring the burstiness of star-formation. Because Paβ is much less affected by attenuation than the Balmer lines, we will be able to measure star-formation histories without dealing with the potential biases due to dust, a significant problem for optical emission-line to UV ratio star-formation history studies in the literature (Broussard et al. 2019).

SUMMARY AND CONCLUSIONS
We have analyzed a sample of 29 low-redshift (z < 0.287) Paβ emitting galaxies. The galaxies of our sample have been divided into two samples: a primary sample selected such that the HST G141 grism detects Paβ with an observed frame wavelength of λ ≤ 16500Å with a signal-to-noise ratio of σ ≥ 3 (20 galaxies), and a secondary sample of 9 galaxies with reliable spectroscopic redshifts from other lines (from either ground-based optical spectroscopy or from Hα emission in the G102 spectrum) that have Paβ detected with SNR > 1. We also required that the objects in our sample have minimally contaminated 2-D spectra by visual inspection.
We show that Paβ is a valuable indicator of star-formation rates and star-formation histories when compared to more widely used SFR indicators such as Hα flux and continuum emission, especially in moderately-dusty to highlydusty galaxies. Paβ as an indicator of SFR serves as a solution to the issues with dust attenuation experienced by SFRs based on optical emission line tracers such as Hα.
Our study of these Paβ emitting galaxies provides two primary findings: • SFRs calculated from Paβ probe a shorter timescale than probed by continuum emission, so we can draw conclusions about the star-formation histories of galaxies by comparing the two. We consider that the scatter in the Paβ/UV SFR ratio is due to the variability of dust attenuation, and rule this out based on observations of the full optical-to-mid-IR spectral energy distributions. Rather, we argue that the scatter in the Paβ/UV SFR ratio indicates increased stochasticity of the star-formation histories, which increases with decreasing stellar mass. This is substantiated by our sur- Figure 13. The relation between Paβ luminosity and attenuation-corrected UV star-formation rates for galaxies in our sample, color coded by effective radius (left), Sérsic index (center) and Σ 1kpc (right). Filled circles correspond to the attenuation-corrected UV SFRs and open diamonds correspond to UV + IR ladder SFRs for the twelve galaxies with well-detected IR emission (Barro et al. 2019). Larger symbols represent objects in our primary sample with Paβ SNR>3, and smaller symbols represent objects in our secondary sample with reliable redshifts and marginal Paβ 1<SNR<3. As in Figure 11, the solid black line is the relationship between Paβ luminosity and SFR following Equation 3 (Kennicutt & Evans 2012), and the dashed gray lines represent the 1-sigma detection limits at redshifts z = 0.28 and z = 0.1. None of these quantities have an apparent correlation with the ratio of Paβ to attenuation-corrected UV star-formation rates, shown quantitatively in Figure 14. Figure 14. The log ratio of the Paβ and attenuation-corrected UV SFRs with galaxy effective radius (left) Sérsic index (center), and central density log Σ1kpc (right). Solid circles represent attenuation-corrected Paβ/UV ratios, and open circles represent the same objects without attenuation-corrected Paβ. Larger symbols represent objects in our primary sample with Paβ SNR>3, and smaller symbols represent objects in our secondary sample with reliable redshifts and marginal Paβ 1<SNR<3. We find no significant correlations between the attenuation-corrected Paβ/UV ratio and effective radius, Sérsic index, or central density. The black dash-dot line shows the median attenuation-corrected Paβ/UV ratio of 0.3. vival analysis analysis, and agrees with other, independent observations in previous work.
• Paβ/Hα ratios serve as a valuable indicator of dust attenuation when we compare these to the Balmer decrement and continuum attenuation estimates, notably in moderately to severely dusty galaxies. Paβ/Hα has the same insensitivity to nuisance parameters such as metallicity, temperature, and density as the Balmer decrement, but does not risk miscalculating attenuation for ISM regions optically thick to Hβ. The Balmer decrement also develops large uncertainties in only moderately dusty galaxies due to poorly constrained Hβ emission.
Our results motivate future IR observations of Paschen series lines for measuring star-formation rate. The James Webb Space Telescope (JWST) will reach a flux limit that is an order of magnitude fainter than our CLEAR data for similar exposure times, enabling detection of fainter Paschen-line emission in low stellar mass galaxies where our work has to rely on survival analysis. In addition, the broad 1-5 µm spectroscopic coverage of JWST includes the Paα line, which is twice as bright than Paβ, for galaxies over z < 1.65. Future JWST observations of Paschen-line emission in galaxies are likely to reveal a much more complete picture of starformation and bursty formation histories, especially in galaxies with significant dust attenuation.
NIRCam grism slitless spectroscopy of the Balmer lines in conjunction with Paβ or Paα can offer more accurate Paschen-to-Balmer dust attenuation measurements, without the need for relying on ground-based optical measurements. This would offer a more accurate census of dust attenuation in highly dusty galaxies without making assumptions about slit losses which may fail for certain morphologies.
Software: linmix (Kelly 2007), grizli pipeline (Brammer et al. 2008), FAST (Kriek et al. 2009), EAZY (Brammer et al. 2008;Wuyts et al. 2011), GALFIT (Peng et al. 2010 Figure  12). The solid circles are all of the objects in our sample Paβ detections of SNR>1. Open diamonds are randomly generated data from these non-detections using a half-normal distribution with standard deviation equal to the 1σ Paβ upper limits of the non-detections. Both detections and the randomly generated survival analysis points are color-coded by AV . The gray lines show rolling median absolute deviations of the sample and survival analysis points. The right panel shows the same analysis for the subsample of galaxies with AV < 1. The survival analysis particularly populates the lower-left of the plot where we cannot detect Paβ-emitting galaxies. In order to quantify the scatter, we resample each point 200 times and take the standard deviation for in each of six stellar mass bins. We then fit the standard deviations for each resampling and find a median slope of m = −0.063 ± 0.039 (left panel) and m = −0.151 ± 0.048 (right panel), indicating higher scatter at lower stellar mass with ≈ 1.6σ and ≈ 3.1σ significance, respectively. The higher significance in the right panel, where moderate and highly dusty objects are eliminated, indicates that the scatter in the Paβ/UV ratio is not due to variable attenuation. We interpret this higher scatter at lower stellar mass to be evidence for burstier SFHs at lower stellar mass, consistent with previous work (Broussard et al. 2019;Guo et al. 2016;Weisz et al. 2014 Objects with spectroscopic redshifts are quoted without redshift uncertainty. Objects only with grism redshifts are quoted with uncertainties. ii From the 3D-HST catalog Skelton et al. (2014) iii From the CANDELS/SHARDS catalog Barro et al. (2019) iv From the GALFIT catalog van der Wel et al. (2012)