Dense Gas, Dynamical Equilibrium Pressure, and Star Formation in Nearby Star-forming Galaxies

As part of ALMA's Cycle 2 campaign, we observed HCN (1–0), HCO⁺ (1–0), CS (2–1), ¹³CO (1–0), and C¹⁸O (1–0) in four galaxies using ALMA's main array of 12 m antennas. For the remainder of the paper, we will refer to these lines as HCN, HCO⁺, CS, ¹³CO, and C¹⁸O. HCN, HCO⁺, and CS all have n_eff ∼ 10⁴–10⁵ cm⁻³ (see Table 1) and so are expected to trace mainly dense gas (though in the absence of such gas they can still emit; e.g., Shirley 2015; Leroy et al. 2017b). The CO isotopologues, ¹³CO and C¹⁸O, trace lower-density gas, n_eff ∼ 10³ cm⁻³. The contrast between the optically thin isotopologues and the optically thick ¹²CO emission constrains the optical depth and physical conditions in the bulk of the molecular gas. We make limited use of ¹³CO and C¹⁸O in this paper. These data are analyzed in detail by Jiménez-Donaire et al. (2017a). Fainter lines in the bandpass were analyzed by Jiménez-Donaire et al. (2017b).

Table 2 gives our adopted position, morphology, orientation, distance, beam size, and field of view for each target. We observed seven fields in a hexagonally packed mosaic pointed toward the center of each galaxy. The mosaic pattern used the default Nyquist spacing set by the ALMA observing tool.

We observed each galaxy with two spectral setups. The first covered lines tracing the dense gas: HCN, HCO⁺, and CS. The four spectral windows covered 85.4–87.2 GHz, 87.2–89.0 GHz, 97.2–99.0 GHz, and 99.0–100.8 GHz. The second spectral setup covered lines tracing the overall distribution of molecular gas: ¹³CO and C¹⁸O. Those four spectral windows covered 98.2–100.0 GHz, 96.6–98.4 GHz, 108.5–110.3 GHz, and 110.3–112.1 GHz. For both setups, we observed using a channel width of 976.6 kHz (∼3 km s⁻¹ at ν = 100 GHz) and bandwidth of 1.875 GHz, sufficient to cover the full velocity extent of each line in question.

We observed in a compact configuration in order to emphasize flux recovery and surface brightness sensitivity, reflecting the faint nature of the dense gas tracers. After calibration, the HCN observations had 703 (NGC 3351), 630 (NGC 3627), 561 (NGC 4254), and 561 (NGC 4321) baselines, with minimum and maximum unprojected lengths of 15 and 348 m, respectively, median baseline length of 90–100 m, and 20%–90% range of typically 50–195 m. For reference at the ν ∼ 89.5 GHz of the HCN and HCO⁺ lines, 50, 100, and 200 m correspond to ∼138, 69, and 35.

We processed the data using the the CASA software package (McMullin et al. 2007) and the observatory-provided calibration scripts. Most of the calibration occurred in CASA version 4.2.2, with one data set calibrated in CASA version 4.3.1. The calibration scripts were a mixture of the observatory-produced CASA scripts and calls to the formal ALMA pipeline. In all cases they are available, along with the data, from the ALMA archive. After inspecting the pipeline calibrated data, we imaged each line separately. For the final version of the imaging, we used CASA version 4.6.0.

We first subtracted continuum emission using the CASA task uvcontsub and avoiding the frequencies of known bright lines. We then imaged each cube using natural weighting, averaging in frequency to produce 10 km s⁻¹ wide channels, and applying a small u − v taper (2''–3'' depending on the line and target). The taper further emphasizes surface brightness sensitivity. The small loss of resolution is irrelevant to the science in this paper because the comparison to tracers of recent star formation already limits our work to ≳5'' resolution.

With the taper, but before any other processing, the synthesized beams in the deconvolved HCN images were 42 × 36 (NGC 3351), 44 × 37 (NGC 3627), 60 × 35 (NGC 4254), and 46 × 32 (NGC 4321). The pixel size adopted during imaging was always chosen to heavily oversample the beam. After imaging, we convolved each cube using the CASA task imsmooth to have a round 8'' × 8'' Gaussian beam. This allowed us to beam-match the poorer-resolution CO and infrared data that are crucial to the analysis. The final images used in this analysis all have 8'' (FWHM) beams, 10 km s⁻¹ channel width, and 05 pixels that heavily oversample the beam.

Given the u − v coverage mentioned above, we expect structures larger than ∼14'' in a single channel to suffer from spatial filtering in the ALMA main array data. To account for this, we combined our 12 m HCN, HCO⁺, ¹³CO, and C¹⁸O data with single-dish maps obtained as part of the IRAM EMPIRE Survey (Bigiel et al. 2016; Cormier et al. 2018, M. Jimenez Donaire et al. 2018, in preparation). To do this, we aligned the IRAM maps to the grid of the ALMA data, applied the primary beam response of the ALMA images to the IRAM data, and converted the IRAM data to have units of Jy beam⁻¹. Then we combined the two data sets using the CASA task feather. After the combination, we verified that the feathered cube indeed matched the spectral profile of the IRAM 30 m cube when both were convolved to a common 30'' resolution. Short-spacing data were not available for the CS, so those data are from ALMA's main 12 m array only in this paper.

Table 3 reports the total flux recovered for each line from each galaxy both with and without the addition of the IRAM 30 m data. We calculate the total flux by summing the pixels in the original data cubes. The table shows that ALMA recovers ≳95% of the flux for all lines in NGC 3351, which is dominated by a bright, compact nuclear source. In fact, for the two faintest lines, ALMA finds slightly more flux than the IRAM 30 m, indicating modest (∼10%) discrepancies in the flux calibration scale. ALMA recovers a lower fraction of the flux for galaxies with more extended, low signal-to-noise ratio (S/N) emission. On average, ALMA recovers ∼60%–80% of the flux found by the IRAM 30 m. Our worst case is NGC 4254, where the faint line emission appears extended and ALMA recovers only ∼30%–50% of the flux seen by the IRAM 30 m. This likely results from two factors: the large extent of the CO emission in the galaxy and the lack of a compact, bright nuclear source (in contrast to the other three targets). This variable, sometimes poor flux recovery emphasizes the importance of including short-spacing data.

Throughout the paper, we work with intensity in units of brightness temperature, T_B. We convert our final data cubes from their native units of Jy beam^–1 to T_B via the standard Rayleigh–Jeans formula:

$$\begin{eqnarray}&&{T}_{B}[{\rm{K}}]=\displaystyle \frac{{c}^{2}}{2\times {10}^{23}\,{k}_{B}\,{\nu }^{2}}\,{I}_{\nu }\,[\mathrm{Jy}/\mathrm{beam}],\end{eqnarray} \tag{ 1 }$$

where ν is the frequency of the line, k_B is Boltzmann's constant, and I_ν is the original intensity in Jy beam⁻¹.

Finally, we measure the rms noise from the signal-free region of each cube. At 8'' resolution, the statistical noise in each 10 km s⁻¹ channel is ∼5–10 mK for HCN, HCO⁺, and CS and ∼10 mK for ¹³CO and C¹⁸O. These vary somewhat from cube to cube, and we use the correct local value to construct uncertainty maps for the integrated intensity. The nominal flux calibration accuracy of ALMA in Band 3 during Cycle 2 was 5% according to the ALMA Cycle 2 Technical Handbook, though this may be somewhat optimistic. The IRAM 30 m intensity calibration for EMPIRE observations is internally consistent to ≈5% from night to night, but the absolute calibrations scale is uncertain at the 10%–15% level (M. Jimenez Donaire et al. 2018, in preparation).

We use ¹²CO (1–0) emission, hereafter referred to as CO, to trace the overall molecular gas reservoir. We draw these maps from literature data. In each case, our CO data include both interferometric and single-dish data, allowing us to reach our working ∼8'' resolution, but also recover zero- and short-spacing information. All of the CO data have a larger field of view than our ALMA maps, bandwidth that covers the entire CO line for the galaxy, and pixels that oversample the beam by more than the Nyquist rate.

For NGC 3351 and NGC 3627, we use CO observations from BIMA SONG (Helfer et al. 2003). These cubes include data from both the BIMA interferometer and the NRAO 12 m single-dish telescope on Kitt Peak. The combined BIMA+12 m cubes have native resolutions of 74 × 52 × 10 km s⁻¹ (NGC 3351) and 73 × 58 × 10 km s⁻¹ (NGC 3627). Helfer et al. (2003) quote a flux calibration uncertainty of 15%.

For NGC 4254, we use interferometric CO observations from CARMA STING (Rahman et al. 2011) and single-dish data from the CO extension to the IRAM EMPIRE survey (Cormier et al. 2018, M. Jiménez-Donaire et al. in preparation). Before convolution, the native resolution of the CARMA data is 33 × 27 × 5 km s⁻¹. Rahman et al. (2011) do not quote an amplitude calibration uncertainty but discuss ∼10% as a typical value. We combine the CARMA and IRAM data using the CASA task feather, which carries out a Fourier plane combination of the two cubes.

For NGC 4321, we use CO data provided as part of the ALMA science verification program. These include both main 12 m array and Atacama Compact Array short-spacing and total power data. We use the feathered "reference image" provided as part of the science verification release, which has resolution 39 × 25 × 5 km s⁻¹. As above, the nominal amplitude calibration accuracy for ALMA at Band 3 is ±5%, though in practice this seems optimistic. We convert all CO data to units of brightness temperature following Equation (1), convolve them to our working 8'' resolution, and align them to the astrometric and velocity grid of the ALMA dense gas data.

Whenever possible, we report our results in terms of observable quantities, e.g., I_CO, I_HCN, etc. We also express our results in terms of physical quantities, e.g., molecular gas mass and dense gas mass. The translation from observed to physical quantities carries substantial uncertainty and remains a topic of active research (e.g., Sandstrom et al. 2013; Usero et al. 2015; Leroy et al. 2017b). However, these physical quantities—e.g., Σ_H2 and Σ_dense—are of considerable interest and are central to our science goals. Therefore, following Usero et al. (2015), we will also express our results as best-estimate physical terms. Given a choice, we prefer the simplest possible translation from observed to physical quantities and will plot both axes whenever we can.

By default, we quote total molecular gas mass assuming a CO (1–0) to molecular gas mass conversion factor of ${\alpha }_{\mathrm{CO}}\approx 4.3\,{M}_{\odot }\,{\mathrm{pc}}^{-2}\,{({\rm{K}}\mathrm{km}{{\rm{s}}}^{-1})}^{-1}$. In order to derive a fiducial dense gas mass, we convert from HCN intensity to dense gas mass surface density assuming ${\alpha }_{\mathrm{HCN}}\approx 10\,{M}_{\odot }\,{\mathrm{pc}}^{-2}\,{({\rm{K}}\mathrm{km}{{\rm{s}}}^{-1})}^{-1}$. Both factors include helium.

Both of these conversion factors carry substantial uncertainty, though for different reasons. Our ${\alpha }_{\mathrm{CO}}\,\approx 4.3\,{M}_{\odot }\,{\mathrm{pc}}^{-2}\,{({\rm{K}}\mathrm{km}{{\rm{s}}}^{-1})}^{-1}$ is often taken as the Milky Way conversion factor and applied as a default to solar-metallicity massive galaxies (see Bolatto et al. 2013). However, dust-based studies that include our current targets suggest that the gas-rich regions in the central parts of some galaxies have lower α_CO (Sandstrom et al. 2013). This presumably reflects dynamical broadening of the CO line width, leading to more CO emission per unit mass.

Meanwhile, we consider α_HCN uncertain because the abundance and opacity of HCN as a function of density and environment remain poorly known. Gao & Solomon (2004) argue for ${\alpha }_{\mathrm{HCN}}\approx 10\,{M}_{\odot }\,{\mathrm{pc}}^{-2}\,{({\rm{K}}\mathrm{km}{{\rm{s}}}^{-1})}^{-1}$ to convert from HCN intensity to surface density of gas above ${n}_{{\rm{H}}2}\sim 3\times {10}^{4}$ cm⁻³ based on large velocity gradient modeling and invoking the virial theorem. However, if the inputs to these calculations, e.g., the abundance of HCN, the opacity of HCN, or the dynamical state of dense gas, vary from their assumed values, then the absolute value of α_HCN also changes (see, e.g., Leroy et al. 2017b). Constraints on these quantities in other galaxies remain very weak (see Martín et al. 2006; Jiménez-Donaire et al. 2017a), but observations of the Milky Way do demonstrate important variations (e.g., Pety et al. 2017).

In this paper, we focus on the dense gas fraction. Many environmental factors that affect CO might also affect HCN emission. Because we lack a useful prescription for α_HCN, it is unclear how we should implement variations in α_CO (for an in-depth discussion of how these quantities may relate, see Usero et al. 2015). Moreover, we aim to clearly present our results in terms of observed line ratios, e.g., ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$. This requires a simple α_CO and α_HCN.

We do consider how our multiline observations support the idea that HCN traces the dense gas. For a complete discussion of the variation in α_CO, α_HCN, and their effects on f_dense and SFE_dense, we refer the reader to Usero et al. (2015).

For each galaxy, we create a position–position–velocity mask using the CO data. We first identify pixels with S/N > 5 in at least two adjacent velocity channels. We then remove contiguous S/N > 5 regions that are small compared to the size of the beam. Next, we grow these regions to include adjoining pixels where S/N > 2. The resulting mask does a good job of capturing the region of bright CO emission that one would identify by eye.

CO emission tends to be brighter and easier to excite than emission from rarer isotopologues or dense gas tracers. Therefore, we take this region of bright CO emission to also represent the position–position–velocity region where we might find these fainter lines. We verify by eye that no clear dense gas tracer emission extends beyond this mask. Thus, regions outside this mask are taken to be signal free and used to estimate the rms noise in each cube.

We sum the masked line data cubes along each line of sight to produce maps of integrated intensity, in K km s⁻¹. In this sum, masked pixels have a value of zero, and lines of sight with no identified signal have an integrated intensity of zero. Thus, we also produce a 2D mask indicating which lines of sight include any bright signal in the CO cube.

We calculate the statistical uncertainty in the integrated intensity from the rms noise in an individual channel map (in K) times the channel width (in km s⁻¹) times the square root of unmasked voxels along a line of sight. Typical rms uncertainties due to statistical errors in the integrated lines are ∼0.2 K km s⁻¹ for the CO isotopologues and ∼0.1 K km s⁻¹ for the dense gas tracers. The CO data tend to have poorer sensitivity, with typical uncertainty ∼2 K km s⁻¹ in the integrated intensity per line of sight. However, because the CO is brighter, the CO data still have higher S/N.

The deep multiwavelength data available for our sample allow the prospect to estimate the surface density of recent star formation, Σ_SFR, in several ways. In our table of radial profiles, we provide the measurements necessary to reconstruct most popular monochromatic or "hybrid" tracers (see, e.g., review in Kennicutt & Evans 2012).

In the plots accompanying the main text, we present Σ_SFR estimated from a combination of Hα and 24 μm emission following Calzetti et al. (2007) and Leroy et al. (2012). We adopt this choice because for nearby nonstarburst galaxies it has emerged as a widely accepted estimator of Σ_SFR (see Kennicutt & Evans 2012).

Although Hα+24 μm is widely accepted, its calibration remains fundamentally empirical, with limited fundamental work on resolved targets (Calzetti et al. 2007; Kennicutt et al. 2007; Murphy et al. 2011). Much of the literature surrounding dense gas remains focused on infrared luminosity, considering the estimated bolometric luminosity (reprocessed by dust) of the young stellar population a more physical tracer of recent star formation activity. Therefore, as a complement to the Hα+24 μm estimates shown in the main plots, appendix investigates estimates of the total infrared (TIR) surface brightness. We convolve multiband Herschel and Spitzer imaging to a common, coarse resolution and apply the prescriptions of Galametz et al. (2013) for each galaxy. Based on this, we calibrate galaxy-by-galaxy conversions from the higher-resolution 70 μm data to Σ_TIR. These data are provided in the table of radial profiles and can be used to carry out a parallel analysis to the one presented in the main text.

The numerical results from these two main approaches to Σ_SFR differ at the ∼30% level. However, the choice of SFR tracer does not affect our qualitative conclusions. For either case, dense gas fraction appears to be a poor predictor of the bulk star formation efficiency, and the SFR per unit dense gas drops in the high surface density, inner parts of several of our targets.

2.5.1. Hα+24 $\mu m$

We calculate the surface density of recent star formation using the prescription of Calzetti et al. (2007) recast into surface brightness units. Their prescription (their Equation (7)) in luminosity units is

$$\begin{eqnarray}&&\mathrm{SFR}\,[{M}_{\odot }\,{\mathrm{yr}}^{-1}]=5.3\times {10}^{-42}({L}_{{\rm{H}}\alpha }+0.031\,{L}_{24}).\end{eqnarray} \tag{ 2 }$$

Here L_Hα is the line-integrated Hα luminosity, with units of erg s⁻¹ and no correction for internal extinction. ${L}_{24}\,\equiv \nu {L}_{\nu }(24\,\mu {\rm{m}})$ is the 24 μm luminosity, also in units of erg s⁻¹. This luminosity comes from multiplying the specific luminosity, L_ν, by the frequency at λ = 24 μm.

For a measured line-integrated intensity, I_Hα, in erg s⁻¹ cm⁻² sr⁻¹, across a beam with angular area Ω at distance d, the corresponding SFR surface density, Σ_SFR, is

$$\begin{eqnarray}&&{{\rm{\Sigma }}}_{\mathrm{SFR}}=\displaystyle \frac{\mathrm{SFR}}{{\rm{\Omega }}\,{d}^{2}},\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&=\,\displaystyle \frac{5.3\times {10}^{-42}{I}_{{\rm{H}}\alpha }\,{\rm{\Omega }}\,4\,\pi \,{d}^{2}}{{\rm{\Omega }}\,{d}^{2}}.\end{eqnarray} \tag{ 4 }$$

Adopting the standard units for Σ_SFR of ${M}_{\odot }\,{\mathrm{yr}}^{-1}\,{\mathrm{kpc}}^{-2}$, the d in the numerator is in cm, whereas the d in the denominator is in kpc. We correct this unit inconsistency by multiplying the numerator by 9.52 × 10⁴² cm² kpc⁻². This yields ${{\rm{\Sigma }}}_{\mathrm{SFR}}\,=634\,{I}_{{\rm{H}}\alpha }$ with Σ_SFR in ${M}_{\odot }\,{\mathrm{yr}}^{-1}\,{\mathrm{kpc}}^{-2}$ and I_Hα in erg s⁻¹ cm⁻² sr⁻¹.

To place the 24 μm term in surface brightness units, we combine the 634 prefactor from the previous calculation, the 0.031 relative weighting of luminosities from Equation (2), the conversion of $1\,\mathrm{MJy}\,{\mathrm{sr}}^{-1}={10}^{-17}$ erg s⁻¹ cm⁻² sr⁻¹ Hz⁻¹, and $\nu (24\,\mu {\rm{m}})\approx 1.25\times {10}^{13}$ Hz. Then

$$\begin{eqnarray}&&{{\rm{\Sigma }}}_{\mathrm{SFR}}=634\,{I}_{{\rm{H}}\alpha }+0.0025\,{I}_{24\mu {\rm{m}}},\end{eqnarray} \tag{ 5 }$$

again with Σ_SFR in ${M}_{\odot }\,{\mathrm{yr}}^{-1}\,{\mathrm{kpc}}^{-2}$, I_Hα in erg s⁻¹ cm⁻² sr⁻¹, and I_{24 μm} in MJy sr⁻¹.

We measure 24 μm intensity from the SINGS maps (FWHM ∼ 64; Kennicutt et al. 2003). To estimate I_Hα, we use the collection of literature Hα maps compiled and processed by Leroy et al. (2012). These draw heavily on LVL and SINGS (Kennicutt et al. 2003; Lee et al. 2009). This estimate assumes a Kroupa initial mass function (IMF). The correction to a Chabrier IMF is modest.

We convolve both the Hα and 24 μm maps to the 8'' resolution of the millimeter-wave data. To convolve the 24 μm, we use the kernels of Aniano et al. (2011). Because of the extended structure of the 24 μm point-spread function (PSF), translating from the native PSF to an 8'' Gaussian requires the use of an "aggressive" kernel. This may be another reason to prefer the 70 μm oriented alternative described below.

In appendix, we show that in our sources and target regions the 24 μm term dominates Σ_SFR. Hα itself typically contributes ∼30% to the overall estimate.

2.5.2. Indirect Total-infrared-based SFR Estimates

As an alternative to the Hα+24 μm approach, appendix works out an SFR estimate based on the TIR emission in our targets. We match the resolution of the Spitzer and Herschel maps at 24, 70, 160, and 250 μm at 30''.¹⁸ Then, we use the prescriptions of Galametz et al. (2013) to calculate Σ_TIR, the TIR luminosity surface density. This corresponds to the expected integral under the whole IR spectral energy distribution from 8 to 1000 μm. If the light in the region under consideration is dominated by an embedded young stellar population, then this offers an alternative tracer of recent star formation.

The full suite of IR bands is only available at poor (22'') resolution. To estimate the small-scale behavior of Σ_TIR, we use the high-resolution 70 μm maps from Herschel (Kennicutt et al. 2011). In appendix, we show that 70 μm emission correlates better with Σ_TIR than the 24 μm does at low resolution. We calibrate a simple conversion from 70 μm to TIR for each galaxy in our sample. This conversion is not perfect, but assuming that it holds subresolution allows us to estimate Σ_TIR at the 8'' common resolution of our millimeter-wave data.

We present these "TIR"-based estimates in the table or radial profiles, adopting the conversion from Σ_TIR to Σ_SFR calculated by Murphy et al. (2011):

$$\begin{eqnarray}&&\left(\displaystyle \frac{{{\rm{\Sigma }}}_{\mathrm{SFR},\mathrm{TIR}}}{{M}_{\odot }\,{\mathrm{yr}}^{-1}\,{\mathrm{kpc}}^{-2}}\right)=1.48\times {10}^{-10}\left(\displaystyle \frac{{{\rm{\Sigma }}}_{\mathrm{TIR}}}{{L}_{\odot }\,{\mathrm{kpc}}^{-2}}\right).\end{eqnarray} \tag{ 6 }$$

Murphy et al. (2011) derive this conversion assuming that the entire Balmer continuum is absorbed and re-radiated by dust and assuming a continuous star formation history with no contamination due to an old background stellar population. Based on comparison with UV emission, all of our targets appear heavily extinguished across the region of interest (so that adding a GALEX-based term would not influence ${{\rm{\Sigma }}}_{\mathrm{SFR},\mathrm{TIR}}$ significantly). However, variations in the star formation history or heating by older stellar populations remain a concern (though these are among the most gas-rich, IR-bright, actively star-forming regions in the KINGFISH Herschel survey).

This TIR estimate follows a less conventional approach than the Hα+24 μm, but it represents our preferred SFR estimate for comparison to the literature. We recommend adopting these values from the table of radial profiles when carrying out such a comparison. However, in the main text, we will plot results deriving SFR from Hα+24 μm.

2.5.3. Literature Star Formation Rate Estimates

Usero et al. (2015) and Bigiel et al. (2016) highlighted strong correlations with stellar structure, a possible driver of gas pressure. With this in mind, we compare our observations to 3.6 μm maps, which are dominated by light from old stars. We use maps obtained by Spitzer, mostly as part of SINGS (Kennicutt et al. 2003), and processed as part of the S⁴G survey (Sheth et al. 2010). Because Spitzer's 3.6 μm band can include emission associated with recent star formation, primarily dust emission, we use the contaminant-corrected maps of Querejeta et al. (2015), which are based on an algorithm developed by Meidt et al. (2012). Globally, the correction of Querejeta et al. (2015) and Meidt et al. (2012) attributes 10%–30% of the 3.6 μm emission to a component associated with star formation and removes it from the maps. The PSF of the IRAC 3.6 μm maps has FWHM ∼ 19 before any convolution.

We calculate the stellar surface density from the contaminant-corrected 3.6 μm data (Querejeta et al. 2015) assuming a mass-to-light ratio of 352 M_⊙ pc⁻² (MJy sr⁻¹) (Meidt et al. 2014). This roughly corresponds to 0.5 M_⊙ per L_⊙. Note that Usero et al. (2015) and Bigiel et al. (2016) used median-filtered versions of the 3.6 μm maps to correct for contamination, rather than the ICA-corrected versions from Querejeta et al. (2015). Usero et al. (2015) also used a different mass-to-light ratio. After accounting for the different mass-to-light ratio, the contaminant-corrected 3.6 μm map that we use here agrees reasonably with the median-filtered Usero et al. (2015) and Bigiel et al. (2016) data.

The distribution of atomic gas (H i) plays a role in calculating dynamical equilibrium pressure. We compare to the H i maps from THINGS (Walter et al. 2008) for NGC 3351, NGC 3627, and NGC 5194 and archival VLA maps for NGC 4254 and NGC 4321 (from Schruba et al. 2011; Leroy et al. 2013b). We begin with the 21 cm integrated intensity ("moment 0") maps and convert these to hydrogen column density assuming optically thin emission (see Walter et al. 2008).

The H i maps have coarser resolution than our other data. The naturally weighted THINGS maps for NGC 3351 and NGC 3627 have synthesized beams of 99 × 71 and 100 × 89. The archival maps of NGC 4254 and NGC 4321 have synthesized beams of 169 × 162 and 147 × 141. In all cases, the 21 cm map covers the full spectral extent of the galaxy, and the pixels heavily oversample the synthesized beam. However, we note that the VLA maps do not include zero-spacing data and can be expected to filter out emission with angular extent ≳30' in an individual channel. Despite this, Walter et al. (2008) show that for the THINGS maps the agreement between the flux measured in the interferometer-only cube and single-dish measurements is quite good. In practice, we suggest to consider the zero level of these profiles to be perhaps biased slightly low, by ≲2 M_⊙ pc⁻².

In our analysis, we will assume the H i to be smooth subresolution. This agrees with analysis showing the H i to be only weakly clumped (Leroy et al. 2013a) and with the relatively flat shape of the radial profiles (Schruba et al. 2011; see below). The main application of the H i profiles will be to calculate total gas surface density for use in pressure estimates, so the analysis is not terribly sensitive to this assumption.

Usero et al. (2015) and Bigiel et al. (2016), following Helfer & Blitz (1997), suggested that ISM pressure may play a central role in determining the gas density distribution and the role of gas at any particular density in star formation. We estimate the approximate midplane gas pressure needed to maintain a vertical dynamical equilibrium following Elmegreen (1989), Blitz & Rosolowsky (2006), Wong & Blitz (2002), Leroy et al. (2008), and Ostriker et al. (2010). Specifically, we take

$$\begin{eqnarray}&&{P}_{\mathrm{DE}}\approx \displaystyle \frac{\pi G{{\rm{\Sigma }}}_{\mathrm{gas}}^{2}}{2}+{{\rm{\Sigma }}}_{\mathrm{gas}}{(2G{\rho }_{\star })}^{1/2}\,{\sigma }_{\mathrm{gas},z}\,.\end{eqnarray} \tag{ 7 }$$

Here Σ_gas is the total atomic plus molecular gas surface density, ρ_⋆ is the volumetric mass density of stars and dark matter at the midplane, and ${\sigma }_{\mathrm{gas},z}$ is the vertical velocity dispersion of the gas.

P_DE expresses the pressure needed to balance the vertical gravity on the gas in the galaxy disk. The first term in Equation (7) reflects the gas self-gravity, and the second term reflects the weight of the gas in the potential well of the stars. For our target regions, the stellar potential will exceed the gas self-gravity, and we expect the second term to be dominant.

Only the stars and dark matter within the gas layer are relevant to the potential. We expect the density of dark matter within the gas layer to be small compared to that of the stars, and so we neglect that contribution to ρ_⋆. Thus, ρ_⋆ rather than Σ_⋆ enters Equation (7). Estimating the 3D structure of the disk represents the main challenge to gauging P_DE from observations. We need estimates of the local stellar scale height to translate our inferred Σ_⋆ into ρ_⋆ via

$$\begin{eqnarray}&&{\rho }_{\star }\approx \displaystyle \frac{{{\rm{\Sigma }}}_{\star }}{2\,{h}_{\star }}\end{eqnarray} \tag{ 8 }$$

following van der Kruit (1988).

A full treatment of the height of stellar disks as a function of morphology, mass, and galactocentric radius is beyond the scope of this paper. Observations still do not provide a perfect prescription, and more recent observations show that the topic may be more complex than suggested by earlier work (e.g., Comerón et al. 2012, 2014). Lacking a better prescription, we follow Leroy et al. (2008) and adopt the typical flattening ratio found by Kregel et al. (2002), ${l}_{* }/{h}_{* }=7.3\pm 2.2$, to relate the stellar scale height and the observed scale length for the disk component.

We adopt ${\sigma }_{z,\mathrm{gas}}\approx 15$ km s⁻¹. This is toward the high end of the observed gas velocity dispersions in galaxy disks, appropriate for large scales and the gas-rich, high surface density inner regions of galaxy disks (e.g., Tamburro et al. 2009; Caldú-Primo et al. 2013; Mogotsi et al. 2016).

For these assumptions, we calculate P_DE in each ring and report this in the online table described in Appendix C.

Finally, we provide a word on the interpretation of P_DE. P_DE is the time-averaged pressure needed to support the gas disk against its weight due to both self-gravity and the gravity of the stellar potential well. We include all gas in this calculation and expect Equation (7) to hold averaged over time and space (several gas scale heights), with all gas participating in this equilibrium. Thus, we expect P_DE to express the average state of a region of the ISM. Equation (7) does not distinguish between a gravitationally bound component and a diffuse component of the gas. This is consistent with recent simulations (e.g., Kim & Ostriker 2017) and observations that indicate that self-gravitating structures are quite transient (Kawamura et al. 2009; Schruba et al. 2010). Such transient structures may be thought of as fully participating in the ISM's large-scale dynamics when averaged over several scale heights and vertical crossing times. Similarly, for a turbulence-dominated system we do not think of P_DE as an external pressure that acts on cloud "surfaces"; rather, it is the mean pressure averaging over both high- and low-density regions.

We convolved all data to a working resolution of 8'', using the kernels provided by Aniano et al. (2011) when needed. After convolution, we aligned all of the broadband multiwavelength data to the common astrometric grid shared by the HCN, HCO⁺, CS, ¹³CO, and C¹⁸O data. As mentioned above, the H i data could not be beam-matched to an 8'' beam, and we use these at their native resolution, assuming a smooth distribution within the beam.

We use a combination of the convolved Spitzer 24 μm data and Hα data as our primary SFR tracer. Despite the modest ∼6'' core of the 24 μm PSF, the Spitzer PSF has significant extended structure, and an "aggressive" kernel is required to convolve a 24 μm image to have an 8'' Gaussian PSF (Aniano et al. 2011). As a result, mild artifacts are visible surrounding the nucleus of NGC 3351 in the convolved 24 μm image. The alternative 70 μm based tracer has a cleaner PSF. Here again, comparing the two approaches offers a useful way to assess systematic effects.

Note that we use this 8'' working resolution for NGC 3351, NGC 3627, NGC 4254, and NGC 4321. We use a lower 30'' working resolution for NGC 5194, set by the resolution of the IRAM 30 m data. Otherwise, the procedure is the same.

The most striking contrast in our data is between the bright central regions and the surrounding disk, which includes fainter emission from gas in the spiral arms and bars. Even ALMA struggles to detect emission from the extended disk at high significance. To improve our S/N, we conduct much of our analysis using radial profiles of line emission and supporting data. Although azimuthal (e.g., arm–interarm) contrast will also be interesting for future studies, our present data lack the sensitivity to recover faint interarm emission.

We set the width of each radial bin to be half the FWHM of the working beam (15'' for NGC 5194 and 4'' for all other galaxies). This should critically sample the radial structure of the galaxy. It also yields profiles in which adjacent rings are correlated, so that we have twice as many data in each plot than independent radial measurements.

In each ring, we calculate two values: (1) the mean intensity of all pixels within the 2D version of the CO mask, ignoring all other values, and (2) the mean intensity of all pixels in the ring. For profile 2, we treat pixels outside the map as having intensity 0 for the molecular lines, so that the ring-averaged value represents the values in the mask diluted by the empty space outside the mask. For the ancillary data (3.6 μm, Hα, 24 μm, etc.), we calculate the mean value only within the 2D CO-based mask for profile 1 and the mean value over the whole ring for profile 2. We report both values in the table described in Appendix C. Note that for these other data, we do not set any pixels to 0 when taking the mean over the whole ring, but instead use all of the data.

We calculate the uncertainty for these mean intensities by propagating the error from individual intensity measurements. In this case the uncertainty is

$$\begin{eqnarray}&&{{\rm{\Sigma }}}_{{ii}}=\displaystyle \frac{\sqrt{{\sigma }_{1}^{2}+{\sigma }_{2}^{2}\,+\ldots +{\sigma }_{N}^{2}}}{N}\times \sqrt{O}\,.\end{eqnarray} \tag{ 9 }$$

Here Σ_ii is the uncertainty in the integrated intensity averaged over the ring, and N is the number of pixels in the ring. O is a factor (the "oversampling factor") designed to account for the nonindependence of the pixels. If all data are independent with the same uncertainty, ${\sigma }_{{II}}={\sigma }_{1}/\sqrt{N}$ as expected for the uncertainty in the mean of N data. Generally O should be the number of pixels per beam. Here, given that our rings are thinner than a beam, we take O equal to the ring thickness in pixels times the FWHM beam size in pixels. We never allow N/O < 1, so that σ_II is at most the uncertainty in an individual point. When calculating this mean intensity for the line data, pixels outside of the mask are ignored even for profile 2. That is, our uncertainties do not reflect uncertainty in mask construction.

The intensity profiles of each line, along with profiles of the supporting multiwavelength data, are a main observational result of the paper. These are available in full as online material. All surface brightness and surface densities quoted in the paper have been corrected for the effects of inclination.

Even with ALMA, emission from dense gas tracers remains hard to detect. To overcome this, we carry out most of our analysis using azimuthally averaged radial profiles (see Section 2.10). Note that the masked regions are zeros in this averaging. When we report ratios, these are the ratios among the azimuthal profiles taken after radial profile construction.

The middle panels in Figures 1–5 show these profiles for our target lines (left) and key measures of the galaxy and ISM structure (right). Even with azimuthal averaging, we still do not detect the fainter lines in some rings. In these cases we show the 2σ value as an open symbol to indicate an upper limit.

The absolute intensities trace the distribution of material. To capture the changing physical state of the gas, the bottom rows of Figures 1–5 plot ratios among radial profiles. For the lines, we show the intensity of each relative to CO. For the galaxy structure, we show the fractional content of gas (as opposed to stars), the rate of star formation per unit gas, and the balance between atomic and molecular gas. Note that we multiply Σ_SFR by a factor of 10¹⁰ to show it on the same scale as the other quantities.

3.1.1. Galaxy and ISM Structure

3.1.2. Line Ratio Profiles

3.2.1. Dense Gas–Star Formation Correlation

We obtained these observations to test the idea that the dense gas fraction determines the star formation efficiency of molecular gas, ${\mathrm{SFE}}_{\mathrm{mol}}\equiv \mathrm{SFR}/{M}_{\mathrm{mol}}$. In a simple threshold model, the mass of dense gas determines the SFR. Such a relationship would explain the striking correlation between TIR luminosity (commonly used to trace star formation in the HCN literature) and HCN luminosity (tracing dense gas). This correlation extends from Galactic cores all the way to starburst galaxies (e.g., Gao & Solomon 2004; Brouillet et al. 2005; Wu et al. 2010; Rosolowsky et al. 2011; García-Burillo et al. 2012; Kepley et al. 2014; Bigiel et al. 2015, 2016; Usero et al. 2015; see Figure 6).

Zoom In Zoom Out Reset image size

To place our rings on this correlation, we translate our SFRs calculated from Hα+24 μm to the corresponding TIR luminosity expected from an embedded population with this SFR (using Equation (6)). Though indirect, this approach is at least internally consistent with the fiducial SFR estimates in the main text. We use this conversion only for Figure 6. An alternative approach, using our best estimates of the local TIR luminosity based on 70 μm emission (Section 2.5.2), yields similar results and can be carried out with the data provided in the table or profiles.

Figure 6 shows that our rings do populate the HCN–IR luminosity correlation seen at smaller and larger scales (Gao & Solomon 2004; Lada et al. 2012). Our new data occupy the intermediate regime between whole galaxies and Milky Way cores or individual GMCs. In this, they resemble results for the disk pointings of Usero et al. (2015), the ∼30'' NGC 5194 data of Bigiel et al. (2016) and Chen et al. (2015), and recent studies of M82 and the Antennae galaxies (Kepley et al. 2014; Bigiel et al. 2015). Our spatial resolution is finer than that of Bigiel et al. (2016), but our use of rings raises the luminosity in any individual data point.

Taking all of the literature data together, Figure 6 shows that that SFR correlates with HCN luminosity, with about a factor of two scatter across almost 10 orders of magnitude. Even if the luminosity of a specific ring in a galaxy lacks physical meaning, the narrow spread in the TIR-to-HCN ratio is significant. The right panel shows that over a large range, most data in the literature exhibit TIR-to-HCN ratios within ∼±0.3 dex of one another. Although systematic differences are already evident both within and among studies, there is clearly a strong relationship between HCN luminosity and recent star formation.

3.2.2. Does the Dense Gas Fraction Predict the Star Formation Efficiency of Molecular Gas?

One popular approach to thinking about star formation, following, e.g., Gao & Solomon (2004) and Lada et al. (2012), is to imagine a gas density threshold, above which the SFR scales with the mass of dense gas. In such a model, the star formation per unit mass of total molecular gas can and does vary, but star formation will correlate tightly with the mass of dense gas traced, e.g., by HCN. In such a model, any variation in SFR_mol ≡ Σ_SFR/Σ_mol should be due to variations in f_dense (traced here with ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$).

We explore this prediction in Figure 7. We show Σ_SFR/Σ_mol, or the overall efficiency of star formation, as a function of Σ_dense/Σ_mol, or f_dense. We plot results for our ALMA targets, NGC 5194, the Usero et al. (2015) pointings, and the Gao & Solomon (2004) survey.

Zoom In Zoom Out Reset image size

These data all show a similar picture: a correlation exists between ${{\rm{\Sigma }}}_{\mathrm{SFR}}/{{\rm{\Sigma }}}_{\mathrm{mol}}$ and Σ_dense/Σ_mol with the expected sense, that a higher dense gas fraction does imply a higher star formation efficiency. However, the scatter about this relation is large compared to the dynamic range in Σ_dense/Σ_mol or Σ_SFR/Σ_mol.

The weak correlation does not appear to be an artifact of limited dynamic range. We picked our targets to span a large range of Σ_SFR/Σ_mol, with the goal of testing the dense gas threshold hypothesis. Indeed, Figure 7 shows that our galaxies exhibit an order-of-magnitude spread in Σ_dense/Σ_mol. This is a large fraction of the total variation seen in the local universe (e.g., Saintonge et al. 2011; Leroy et al. 2013b, 2015). We also see roughly an order-of-magnitude variation in Σ_dense/Σ_mol across our sample. If a sharp correlation between Σ_SFR/Σ_mol and Σ_dense/Σ_mol was present, we would expect to see it in Figure 7.

The black line in the left panel of Figure 7 shows the expectation for a fixed SFR per unit dense gas. We set the slope to the median Σ_dense/Σ_mol value of our sample and illustrate the ±1σ range as a shaded gray region. This median ratio, Σ_SFR/Σ_dense = 1.10 × 10⁻⁸ with 0.25 dex scatter, is significantly higher than that found by Gao & Solomon (2004), 9.33 × 10⁻¹⁰, but comparable to those found by Usero et al. (2015) and García-Burillo et al. (2012), 1.12 × 10⁻⁸ and 1.45 × 10⁻⁸, respectively. Recall that each of these works calculates SFR differently and that the Gao & Solomon (2004) and García-Burillo et al. (2012) samples focus on starburst systems, whereas our sample consists of more normal star-forming galaxies.

The red dashed line in Figure 7 shows an ordinary least-squares bisector fit to our data. The best-fit power law, weighting all significant detections equally, is

$$\begin{eqnarray}&&{\mathrm{log}}_{10}\displaystyle \frac{{{\rm{\Sigma }}}_{\mathrm{SFR}}}{{{\rm{\Sigma }}}_{\mathrm{mol}}}\approx -8.05+0.92{\mathrm{log}}_{10}\displaystyle \frac{{{\rm{\Sigma }}}_{\mathrm{dense}}}{{{\rm{\Sigma }}}_{\mathrm{mol}}}.\end{eqnarray} \tag{ 10 }$$

This relation should be taken as indicative of the scaling between Σ_SFR/Σ_mol and Σ_dense/Σ_mol in our data. As is visible from the plot, the choice of which data to fit (starbursts, whole galaxies, individual regions, etc.), the weighting applied to different data (e.g., by area, luminosity, or galaxy), and the model used (e.g., including astrophysical scatter or not) all have the potential to alter the fit substantially. This would not be true if Σ_dense/Σ_mol perfectly predicted Σ_SFR/Σ_mol in a universal way across all systems.

The green, blue, and purple lines show similar fits applied to only the Usero et al. (2015), Gao & Solomon (2004), and García-Burillo et al. (2012) data. Considering only the Usero et al. (2015) and García-Burillo et al. (2012) data, we fit slightly steeper slopes than for our data, but we find a consistent normalization near log₁₀f_dense ∼ −1.5. The fit to the Gao & Solomon (2004) fit has a steeper slope and lower intercept than our data, consistent with the offset in measured Σ_SFR/Σ_dense mentioned above. The fits for both Usero et al. (2015) and García-Burillo et al. (2012) are consistent with the fit to our data within the uncertainties.

Thus, our best-fit slope of ∼1 agrees with the schematic picture that Σ_dense/Σ_mol as traced by Σ_HCN/Σ_CO should predict Σ_SFR/Σ_mol in a simple way. However, our results also show that Σ_dense/Σ_mol is not a highly accurate predictor of Σ_SFR/Σ_mol. Across our whole sample, the rank correlation between Σ_SFR/Σ_mol and Σ_dense/Σ_mol is only 0.13. Meanwhile, the overall scatter in Σ_SFR/Σ_mol is 0.33 dex.

Note that individual galaxies do show tighter or looser individual trends in Figure 7, and we report the rank correlation coefficients and logarithmic fits for all of the scaling relations considered in this work in Tables 5–8. ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ can be a good predictor of Σ_SFR/Σ_mol for some galaxies. But overall, ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ remains only a moderately good predictor of Σ_SFR/Σ_mol. The relationship appears particularly poor for our radial profiles, which emphasize the contrast between nuclear regions and the surrounding disks.

Table 5.  Σ_SFR/Σ_mol vs. Σ_dense/Σ_mol

Data Set	Rank Corr.	log10 ΣSFR/Σdense
Radial profiles
NGC 3351	0.77 (0.052)	−7.94 (±0.19)
NGC 3627	−0.01 (0.488)	−7.73 (±0.34)
NGC 4254	0.10 (0.342)	−7.89 (±0.26)
NGC 4321	0.21 (0.202)	−8.26 (±0.23)
NGC 5194a	−0.26 (0.236)	−8.12 (±0.28)
All profiles	0.13 (0.003)	−7.89 (±0.33)
Usero et al. (2015)	0.42 (0.005)	−7.95 (±0.25)
Gao & Solomon (2004)	0.65 (0.000)	−9.03 (±0.24)
García-Burillo et al. (2012)	0.35 (0.075)	−7.84 (±0.20)
All data	0.30 (0.000)	−7.96 (±0.31)

Note. Rank correlation quotes with p value in parentheses; log₁₀Σ_SFR/Σ_dense in units of ${\mathrm{log}}_{10}\,{\mathrm{yr}}^{-1}$. We quote the mean of the logarithm of the ratio and the ±1σ scatter in the log of the ratio. "All data" treats all data points with equal weight, regardless of the spatial scale sampled.

^aBigiel et al. (2016).

Download table as:  ASCIITypeset image

Table 6.  Σ_SFR/Σ_dense vs. Other

Data Set	Rank Corr. (versus Σ*)	Rank Corr. (versus Σmol)
Radial profiles
NGC 3351	−1.0 (0.002)	−1.0 (0.001)
NGC 3627	−0.74 (0.000)	−0.76 (0.000)
NGC 4254	−0.69 (0.001)	−0.69 (0.001)
NGC 4321	−0.90 (0.000)	−0.74 (0.000)
NGC 5194a	−1.0 (0.000)	−1.0 (0.000)
All profiles	−0.71 (0.000)	−0.73 (0.000)
Usero et al. (2015)	−0.60 (0.000)	−0.55 (0.000)
All data	−0.66 (0.000)	−0.72 (0.000)

Note. Rank correlation quotes with p value in parentheses. We quote the mean of the logarithm of the ratio and the ±1σ scatter in the log of the ratio. "All data" treats all data points with equal weight, regardless of the spatial scale sampled.

^aBigiel et al. (2016).

Download table as:  ASCIITypeset image

Table 7.  ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ vs. Other

Data Set	Rank Corr. (versus Σ*)	Rank Corr. (versus ICO)
Radial profiles
NGC 3351	1.0 (0.001)	1.0 (0.001)
NGC 3627	0.81 (0.000)	0.86 (0.000)
NGC 4254	0.43 (0.037)	0.43 (0.037)
NGC 4321	0.71 (0.001)	0.72 (0.001)
NGC 5194a	1.0 (0.000)	1.0 (0.000)
All profiles	0.67 (0.000)	0.56 (0.000)
Usero et al. (2015)	0.78 (0.000)	0.72 (0.000)
All data	0.67 (0.000)	0.56 (0.000)

Note. Rank correlation quotes with p value in parentheses. We quote the mean of the logarithm of the ratio and the ±1σ scatter in the log of the ratio. "All data" treats all data points with equal weight, regardless of the spatial scale sampled.

^aBigiel et al. (2016).

Download table as:  ASCIITypeset image

Table 8.  Other vs. Pressure

Data Set	Rank Corr.	Rank Corr.
	(ΣSFR/Σdense versus)	(Σdense/Σmol versus)
Radial profiles
NGC 3351	−1.0 (0.001)	1.0 (0.001)
NGC 3627	−0.65 (0.000)	0.79 (0.000)
NGC 4254	−0.69 (0.001)	0.43 (0.037)
NGC 4321	−0.70 (0.001)	0.67 (0.002)
NGC 5194a	−1.0 (0.000)	1.0 (0.000)
All profiles	−0.71 (0.000)	0.57 (0.000)
Usero et al. (2015)	−0.47 (0.001)	0.67 (0.000)
All data	−0.54 (0.000)	0.36 (0.000)

Note. Rank correlation coefficients with p value in parentheses. We quote the mean of the logarithm of the ratio and the ±1σ scatter in the log of the ratio. "All data" treats all data points with equal weight, regardless of the spatial scale sampled.

^aBigiel et al. (2016).

Download table as:  ASCIITypeset image

Murphy et al. (2015) find similar results for NGC 3627. They find that the SFE_mol decreases as a function of f_dense for the nucleus and two extranuclear star-forming regions in the disk of NGC 3627. Furthermore, they find that the velocity dispersion of the dense gas decreases with increasing SFE_mol. They conclude that the dynamical state of the dense gas, not its abundance, is what sets the SFR.

Certainly, stars form in the densest parts of local molecular clouds. Why is ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ only a weak predictor of the star formation efficiency of molecular gas? In the next two sections, we consider whether ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ indeed traces the dense gas mass fraction (or at least indicates the mean ISM density) and then present evidence for a context-dependent role for dense gas in star formation (building on García-Burillo et al. 2012; Usero et al. 2015; Bigiel et al. 2016).

Does the surprising lack of correlation between ${{\rm{\Sigma }}}_{\mathrm{TIR}}/{I}_{\mathrm{CO}}$ and ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ reflect a problem with our interpretation of ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$? We take this ratio to trace the fraction of dense gas, f_dense. In this picture, CO traces the bulk of the molecular gas, while HCN traces material with ${n}_{{\rm{H}}2}\gtrsim {10}^{4}$ cm⁻³ (e.g., Gao & Solomon 2004; Krumholz & Thompson 2007). This is the standard interpretation for this ratio, because the most effective density for HCN emission, n_eff, is much higher than both the n_eff for CO and the mean density of a molecular cloud in the solar neighborhood. Meanwhile, CO traces the overall molecular gas mass, though not without caveats (see Bolatto et al. 2013).

Reality may be more complicated. We expect a wide distribution of densities in each beam, and the shape of that distribution can affect our interpretation of the line ratio. Variations in the optical depth, temperature, and chemical abundance of one or both species can also produce line ratio variations that mimic those expected from changing density. Although we lack a "smoking gun," our observations provide several indirect ways to test whether ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ indeed reflects the density distribution in our data.

First, in the left panel of Figure 8 we compare ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ to the average surface density on much larger scales, traced by the average surface brightness of CO in a ring. This large-scale surface density convolves the distribution of clouds and the internal density of each cloud, but it still provides a coarse measure of how dense the molecular ISM is in a given region of our Galaxy. The figure shows that ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ indeed correlates with I_CO, suggesting that high dense gas fractions occur in regions of high mean density.

Zoom In Zoom Out Reset image size

This correlation appears stronger within individual galaxies than taking all of our targets together (see Table 7). Conversion factor variations (Sandstrom et al. 2013) and differences in the sub-beam ISM structure (e.g., Leroy et al. 2016) likely contribute to these galaxy-to-galaxy variations. In that sense, our 8'' ∼ 500 pc resolution still limits this test. A stronger version of this test will be to correlate ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ with cloud-scale (tens of parsec) surface and volume density estimates based on high-resolution CO imaging. We will present such a comparison in M. Gallagher et al. (2018, in preparation).

Our line suite allows us to construct several ratios that constrain the shape of the density distribution and the potential impact of abundance variations. We show these in Figure 9. There we plot ${I}_{\mathrm{HCO}+}/{I}_{\mathrm{CO}}$ and ${I}_{13\mathrm{CO}}/{I}_{\mathrm{CO}}$, each as a function of ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$. These ratios offer a check that other dense gas tracers show the same results as HCN and offer some constraint on the sub-beam density distribution.

Zoom In Zoom Out Reset image size

What do we see? ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ and ${I}_{\mathrm{HCO}+}/{I}_{\mathrm{CO}}$ correlate throughout our data, with a median ratio of ${I}_{\mathrm{HCO}+}/{I}_{\mathrm{HCN}}\,=0.76$ and a scatter of 0.24. The solid black line illustrates a constant ${I}_{\mathrm{HCO}+}/{I}_{\mathrm{HCN}}$ equal to this median ratio, with the black dashed lines showing the scatter in our data. We would expect abundance variations unrelated to density to induce scatter in this relation. The observed correlation suggests that our results would be qualitatively similar if we used ${I}_{\mathrm{HCO}+}/{I}_{\mathrm{CO}}$ to trace dense gas instead of ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$.

Correlations between the ${I}_{\mathrm{HCO}+}/{I}_{\mathrm{HCN}}$ ratio and environment have been studied in great detail in luminous infrared galaxies (LIRGs; e.g., Graciá-Carpio et al. 2006, 2008; Imanishi et al. 2006; Privon et al. 2015). These studies reveal some source-to-source variations. But consistent with our results here, starburst-dominated systems show a mean ratio ${I}_{\mathrm{HCO}+}/{I}_{\mathrm{HCN}}\approx 0.77\pm 0.24$ (Privon et al. 2015) and an overall good correlation between ${I}_{\mathrm{HCO}+}$ and I_HCN. The gray line and region in Figure 9 show this ratio and scatter from Privon et al. (2015).

${I}_{\mathrm{HCO}+}/{I}_{\mathrm{CO}}$ variations appear modestly weaker than ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ variations. That is, the slope in the left panel of Figure 9 appears slightly sublinear. Meanwhile, the right panel shows that the magnitude of ${I}_{13\mathrm{CO}}/{I}_{\mathrm{CO}}$ variations is much smaller than that in ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$. As discussed in Leroy et al. (2017b), this is the expected behavior for a distribution of gas densities similar to a lognormal. The n_eff of optically thin ¹³CO is higher than that of CO, which is usually significantly optically thick. But both n_eff are lower than the mean density of the medium. Thus, the fractional change in the gas above these two densities will vary only weakly as the density distribution shifts. Meanwhile, HCO⁺ has a slightly lower n_eff than HCN. For a curving distribution like a lognormal, we expect variations in HCN emission to be stronger than those in HCO⁺ as the mean of the density distribution shifts (Leroy et al. 2017b). Conversely, the modest changes in ${I}_{13\mathrm{CO}}/{I}_{\mathrm{CO}}$ indicate that a power law or other self-similar distribution would be a poor description of the entire density distribution.

Thus, our data offer support for the idea that ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ traces the dense gas fraction: high ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ occurs in regions of high surface brightness, ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ tracks with ${I}_{\mathrm{HCO}+}/{I}_{\mathrm{CO}}$, and the sense of line ratio variations is that expected from a lognormal-like distribution subresolution. However, variations in the quantitative translation of ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ to f_dense remain likely. Temperature, abundance, and optical depth variations of both lines remain plausible and hard to quantify (e.g., Jiménez-Donaire et al. 2017a). These could affect both α_HCN and α_CO. If both lines are optically thick, then temperature variations might be expected to cancel out to some degree, but this is far from certain. We proceed quoting results in terms of the observable ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ and noting the corresponding f_dense for our fiducial α_HCN and α_CO, but we note the latter as uncertain.

Wide-field mapping of local molecular clouds does indicate α_HCN variations within, e.g., Orion (see Kauffmann et al. 2017; Pety et al. 2017). Usero et al. (2015) considered a wide range of possible α_HCN. They show that such variations are highly unlikely to conspire in a way that returns a fixed star formation efficiency of dense gas. Still, the quantitative translation of line ratios into density variations remains uncertain. More multiline wide-area maps of the Milky Way will be key in this regard. Meanwhile, synthetic multiline modeling of the full line suite (e.g., building on Leroy et al. 2017b), comparison to small-scale ISM structure measurements (M. Gallagher et al. 2018, in preparation), and improved constraints on optical depth (Jiménez-Donaire et al. 2017a) are next goals of this project and the EMPIRE survey.

Our data follow the broad correlation between SFR and HCN luminosity, but we see clear evidence of variations in Σ_SFR/Σ_dense, the apparent efficiency with which dense gas forms stars.

3.4.1. Implied Depletion Time and Efficiency per Freefall Time

3.4.2. Variation of Dense Gas Star Formation Efficiency with Environment

Usero et al. (2015) and Bigiel et al. (2016) showed that stellar surface density (Σ_*) and the molecular-to-atomic-gas ratio act as "third parameters" in the relationship between f_dense and SFE_dense. In their studies, regions with high stellar surface density or high molecular gas fraction also have high dense gas fractions. However, in these same high stellar surface density regions dense gas appears less efficient at forming stars.

Such behavior might be expected if clouds with a high mean density (n₀) still only form stars in local overdensities within the cloud. This would be true if these clouds have a typical, roughly virialized, dynamical state at large scales despite their high n₀. Then, where n₀ is high, the immediately star-forming gas engaged in direct collapse and star formation may lie at still higher densities. Similar arguments have been proposed to explain the behavior of our own Galactic center (e.g., Kauffmann et al. 2013; Kruijssen et al. 2014; Rathborne et al. 2015).

In this scenario, SFE_dense for "dense" gas defined by some fixed density threshold will anticorrelate with n₀. Thus, e.g., we would expect Σ_SFR/I_HCN to go down as the mean density of the ISM goes up and HCN traces more "normal" and less directly star-forming gas. In this case SFE_dense should anticorrelate with environmental quantities related to the mean density of the ISM.

To explore this, we plot Σ_SFR/Σ_dense, where Σ_dense is a scaled I_HCN, or SFE_HCN, as a function of the surface densities of stars and I_CO, tracing molecular gas surface density. Figure 10 shows that Σ_SFR/Σ_dense indeed decreases with increasing Σ_mol. A reasonable fit to the relationship is

$$\begin{eqnarray}&&{\mathrm{log}}_{10}\displaystyle \frac{{{\rm{\Sigma }}}_{\mathrm{SFR}}}{{{\rm{\Sigma }}}_{\mathrm{dense}}}=-6.72-0.72{\mathrm{log}}_{10}{{\rm{\Sigma }}}_{\mathrm{mol}}\,.\end{eqnarray} \tag{ 11 }$$

Here Σ_mol in our sample is the azimuthally averaged intensity of CO emission multiplied by α_CO. This traces molecular gas surface density on large scales (recall that our resolution is already hundreds of parsecs). A more rigorous test of this idea using higher-resolution imaging will be presented in M. Gallagher et al. (2018, in preparation).

Zoom In Zoom Out Reset image size

Stellar surface densities exceed the gas surface density over the whole area of our survey, and the stars play a key role in setting the gravitational potential within which the gas exists. If the gas is in some semblance of vertical dynamic equilibrium, then one expects high Σ_* to also lead to higher average midplane pressure and so a higher average gas density. The right panel of Figure 10 shows Σ_SFR/Σ_dense as a function of Σ_*. As in Usero et al. (2015) and Bigiel et al. (2016), we observe a clear anticorrelation between Σ_SFR/Σ_dense and Σ_*, so that in regions with high stellar surface density dense gas appears worse at forming stars.

Σ_SFR/Σ_dense and Σ_* have a rank correlation coefficient of −0.71, indicating a stronger relationship than between Σ_SFR/Σ_dense and Σ_dense/Σ_mol. Using the bisector method, we fit a power law of

$$\begin{eqnarray}&&{\mathrm{log}}_{10}\displaystyle \frac{{{\rm{\Sigma }}}_{\mathrm{SFR}}}{{{\rm{\Sigma }}}_{\mathrm{dense}}}=-6.08-0.72{\mathrm{log}}_{10}{{\rm{\Sigma }}}_{* }\,.\end{eqnarray} \tag{ 12 }$$

Both Equations (11) and (12) should be taken as indicative. We use all of the data available to us, but note that our sample is in no way statistically representative of the whole galaxy population.

Both panels in Figure 10 support the hypothesis that in a high mean density ISM, the gas that emits HCN is less effective at forming stars than in regions of low mean density. Caveats related to the ability of HCN to trace dense gas apply; we cannot rule out that instead this gas may simply be much better at emitting HCN. However, as discussed in Section 3.3, indirect evidence in our data set, as well as the discussion in Usero et al. (2015), gives us some confidence that we do capture density effects. And as discussed by Usero et al. (2015), HCN emission has so far been widely used to trace dense gas in the extragalactic literature (Gao & Solomon 2004). These results complicate any threshold-style interpretation based on those observations.

We focus on I_CO ∼ Σ_mol and Σ_* because they relate to the mean density of the gas. Assuming that the weight of the gas is balanced by its random motions, this relationship may be captured by considering the dynamical equilibrium pressure, P_tot ∼ P_DE. P_DE will depend on both Σ_* and Σ_mol and may in turn relate to the mean density of a cloud, n₀. We return to this point in Section 3.6.

We also observe ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$, our observational tracer of f_dense, to vary significantly across our sample. Adopting our fiducial conversion factors, our observed ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ implies f_dense between ∼0.5% and ∼25% with a median f_dense ∼ 4%.

This median value is lower than the median value of f_dense ∼ 8% from Usero et al. (2015) and the median value of f_dense ∼ 12% from Gao & Solomon (2004). This deviation could be expected based on the construction of samples, as f_dense is not expected to remain constant. Gao & Solomon (2004) study gas-rich starburst galaxies. Meanwhile, the pointings of Usero et al. (2015) were chosen to be bright in CO and hence rich in gas. We target a wider area with less of a selection bias and thus recover lower dense gas fractions. Still wider area maps covering whole galaxy disks may return even lower median f_dense (a main goal of EMPIRE; see M. Jimenez Donaire et al. 2018, in preparation).

Our median f_dense ∼ 4% is lower than the typical f_dense ∼ 10% found for local clouds by Lada et al. (2012). This may reflect different ways of tracing dense gas, as they use extinction rather than HCN emission and adopt a different nominal density. It may also reflect more sensitivity of our observations to an extended molecular component.

Following the logic laid out in the previous section, we would expect the ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$, tracing f_dense, to correlate with both Σ_* and I_CO, tracing Σ_mol. In Figure 8, we plot ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ as a function of both quantities. The figure shows that ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ indeed increases with increasing I_CO. Using the bisector method, we find a power-law fit for ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ versus I_CO of

$$\begin{eqnarray}&&{\mathrm{log}}_{10}\displaystyle \frac{\mathrm{HCN}}{\mathrm{CO}}=-2.42+0.72{\mathrm{log}}_{10}{I}_{\mathrm{CO}}\,.\end{eqnarray} \tag{ 13 }$$

Again, we emphasize that Equation (13) is specific to our data and that I_CO traces the molecular gas surface density averaged on moderately large (few hundred parsec to roughly kiloparsec) scales.

We also expect to find a higher mean density and hence perhaps more dense gas where the gravitational potential well is deeper. We compare ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ to Σ_* in the right panel of Figure 8. We observe a clear correlation between ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ and Σ_*. This has the sense that in regions with high stellar surface density there is more gas dense enough to produce HCN emission. ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ and Σ_* exhibit a rank correlation coefficient of 0.74, indicating a stronger relationship between f_dense and Σ_* than between ${{\rm{\Sigma }}}_{\mathrm{TIR}}/{I}_{\mathrm{CO}}$ and ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$.

A reasonable fit to our data, including NGC 5194 and the data of Usero et al. (2015), is

$$\begin{eqnarray}&&{\mathrm{log}}_{10}\displaystyle \frac{\mathrm{HCN}}{\mathrm{CO}}=-3.28+0.62{\mathrm{log}}_{10}{{\rm{\Sigma }}}_{* }\,.\end{eqnarray} \tag{ 14 }$$

This is moderately steeper than the relation found by Usero et al. (2015):

$$\begin{eqnarray}&&{\mathrm{log}}_{10}{\displaystyle \frac{\mathrm{HCN}}{\mathrm{CO}}}_{\mathrm{Usero}}=-2.58+0.42{\mathrm{log}}_{10}{{\rm{\Sigma }}}_{* }\,.\end{eqnarray} \tag{ 15 }$$

Thus, following Usero et al. (2015) and Bigiel et al. (2016), but now for maps of the disk and central regions of five galaxies, our observations reveal strong correlations between disk structure and ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$. These correlations extend over roughly an order of magnitude in ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ and span two orders of magnitude in stellar surface density and I_CO. We interpret this as evidence for a strong link between the dense gas fraction and local disk structure.

The large-scale gas surface density, traced by I_CO, and the stellar surface density, Σ_*, appear to affect both SFE_dense and f_dense. Above, we hypothesized that these relationships emerge because the mean density of molecular clouds reflects their large-scale environment. Specifically, following Elmegreen (1989), Ostriker et al. (2010), Kim et al. (2011), and others, we noted that the state of the ISM should self-regulate so that the pressure averaged over some relatively large time and length scales will approach the "dynamical equilibrium" value. In this case, molecular clouds will represent overpressured regions against a background pressure that reflects the environment. Helfer & Blitz (1997) argued for a version of this idea to explain the variable dense gas fractions observed in nearby galaxies, and both Usero et al. (2015) and Bigiel et al. (2016) noted this as one plausible family of explanations for their observations. Galactic center studies have considered similar scenarios.

To explore the role of pressure directly, we estimate the dynamical equilibrium pressure, P_DE, for each ring in our sample (see Section 2.8). In Figure 11, we plot ${{\rm{\Sigma }}}_{\mathrm{dense}}/{{\rm{\Sigma }}}_{\mathrm{mol}}(\propto {I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$) and ${{\rm{\Sigma }}}_{\mathrm{SFR}}/{{\rm{\Sigma }}}_{\mathrm{dense}}(\propto {{\rm{\Sigma }}}_{\mathrm{SFR}}/{I}_{\mathrm{HCN}}$) as a function of P_DE.

Zoom In Zoom Out Reset image size

A reasonable fit relating P_DE and Σ_dense/Σ_mol is

$$\begin{eqnarray}&&{\mathrm{log}}_{10}\displaystyle \frac{{{\rm{\Sigma }}}_{\mathrm{dense}}}{{{\rm{\Sigma }}}_{\mathrm{mol}}}=-4.01+0.48{\mathrm{log}}_{10}{P}_{\mathrm{DE}}\,.\end{eqnarray} \tag{ 16 }$$

P_DE and ${{\rm{\Sigma }}}_{\mathrm{dense}}/{{\rm{\Sigma }}}_{\mathrm{mol}}$ relate via a rank correlation coefficient of 0.60. The sense of this correlation, as expected, is that regions with high P_DE exhibit high apparent dense gas fractions. This is the same as the sense of the correlations between I_CO and ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ and between Σ_* and ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$. The strength of the correlations does not recommend one of these as more fundamental, but we note that the correlation with P_DE has a physical explanation that would be expected to yield the other two.

High P_DE corresponds to a higher dense gas fraction, but also a lower star formation efficiency in the dense gas. In our data, a bisector fit between P_DE and Σ_SFR/Σ_dense yields

$$\begin{eqnarray}&&{\mathrm{log}}_{10}\displaystyle \frac{{{\rm{\Sigma }}}_{\mathrm{SFR}}}{{{\rm{\Sigma }}}_{\mathrm{dense}}}=-5.15-0.50{\mathrm{log}}_{10}{P}_{\mathrm{DE}}\,.\end{eqnarray} \tag{ 17 }$$

The two quantities show a rank correlation coefficient of −0.71, indicating a reasonably strong anticorrelation. Again, this goes in the same direction as the anticorrelation between Σ_* and Σ_SFR/Σ_dense and between I_CO and Σ_SFR/Σ_dense. In detail, the slopes differ, with somewhat shallower slopes for the P_DE relations than those treating Σ_* or I_CO as the independent variable.

Should P_DE affect the mean cloud density? Almost all regions that we consider have ${P}_{\mathrm{DE}}/{k}_{B}\gtrsim {10}^{5}$ cm⁻³ K. We expect that the average pressure of the ISM must approach these values. Compare this to the typical internal pressure for a local molecular cloud with Σ_mol ∼ 100 M_⊙ pc⁻². Taking ${P}_{\mathrm{int}}\sim G\pi {{\rm{\Sigma }}}_{\mathrm{mol}}^{2}/2$, we see ${P}_{\mathrm{in}}/{k}_{B}\sim 3\times {10}^{5}$ cm⁻³ K. In the dynamical equilibrium picture, we would expect the internal pressure of molecular clouds to be somewhat higher than P_DE, as molecular clouds represent overpressured regions in a dynamical ISM. In this sense, almost all of our regions have P_DE high enough that P_int cannot remain fixed for clouds across our survey while still balancing the weight of the gas in the stellar potential.

Our result for the dense gas fraction agrees well with the prediction by Meidt (2016), shown by the black line in Figure 11. Meidt considers isothermal clouds bounded by interstellar pressure with ρ ∼ r⁻² density profiles. We show her predicted ${f}_{\mathrm{dense}}={{\rm{\Sigma }}}_{t,\mathrm{HCN}}/{{\rm{\Sigma }}}_{t,\mathrm{mol}}$. Here ${{\rm{\Sigma }}}_{t,\mathrm{HNC}}$ is the column density above the critical density of HCN and ${{\rm{\Sigma }}}_{t,\mathrm{mol}}$ is the column density of the molecular material at the interface with the ambient ISM and is characterized by the ISM pressure. In this view, our calculated dynamical equilibrium pressure is interpreted as an ambient interstellar pressure; this assumption should work as long as the cloud cores modeled by Meidt (2016) do not make up too much of the mass. Figure 11 shows good agreement with most of our data and with the data of Usero et al. (2015). The corresponding prediction for SFE_dense also captures our measured slope and intercept reasonably well.

Is there a threshold pressure for environment to affect dense gas? Conversely, our data do not yield much insight into the question whether the apparent dense gas fraction and SFE_dense continue to correlate with P_DE at the lower values of P_DE expected over most of the area in galaxy disks. One scenario would be that when P_DE falls far below P_int of a typical molecular cloud and H i dominates the ISM, the clouds decouple and converge to a fairly universal population. In this case, the threshold picture of Lada et al. (2012) might still hold for analogs of the solar neighborhood.

Our study, which is optimized to study the contrast between inner molecular disks and the bright central regions of our targets, does not test this scenario. The full EMPIRE survey (M. Jimenez Donaire et al. 2018, in preparation) and full-disk, sensitive ALMA follow-up of our targets should offer more insight.

Our main observational result is that Σ_dense/Σ_mol and Σ_SFR/Σ_dense show correlations with environment that are stronger than the correlation between Σ_SFR/Σ_mol and Σ_dense/Σ_mol. A natural interpretation for this is that the density distribution is a strong function of environment while the role of any particular density in star formation is context dependent.

Environment and Density: We find that the fraction of dense gas depends on environment. In this paper, f_dense is traced mainly by ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$, but HCO⁺ and CS show similar results. Here "environment" means that f_dense appears to increase as Σ_* and I_mol ∼ Σ_mol increase.

Stellar and gas surface densities both relate to the mean pressure, P_DE, needed for the ISM to support its own weight in the (mostly stellar) potential well. Inferring P_DE from observables, we show that f_dense does correlate with P_DE. As P_DE increases, f_dense appears to increase. This follows closely on the results of Helfer & Blitz (1997) that gas becomes denser in a deeper gravitational potential well where the ISM pressure is higher.

We suggest the following sketch to explain this correlation. In order for approximate equilibrium to hold, a disk must self-regulate to have $\sim {P}_{\mathrm{DE}}$ over a moderate length or timescale (Ostriker et al. 2010; Kim et al. 2011, 2013; Ostriker & Shetty 2011; Kim & Ostriker 2015). This dynamical equilibrium captures an average behavior. Molecular clouds will represent overdensities relative to this average. P_DE over large parts of our sample already corresponds to the typical internal pressure of a standard local cloud, ${P}_{\mathrm{int}}/{k}_{B}\sim \pi G{{\rm{\Sigma }}}_{\mathrm{gas}}^{2}/2{k}_{B}\sim 3\,\times {10}^{5}$ cm⁻³ K. For molecular clouds to represent overpressured regions (relative to the time and space average) in our data, their internal pressure must vary, increasing with P_DE.

Our view is that the mean internal pressure of molecular clouds tracks P_DE, and that as P_DE shifts so does the density distribution. That is, in high-pressure regions, the mean density is higher and the overall density distribution also shifts to higher values. We would also expect these regions to have higher line widths, and likely higher Mach numbers, which could also yield a larger dense gas fraction (see, e.g., Padoan & Nordlund 2002; Krumholz & Thompson 2007).

Our own observations, along with those of Usero et al. (2015) and Bigiel et al. (2016), support a scenario where the mean density of molecular clouds varies with the surface density of the disk. Viewing such a scenario through the lens of self-regulation to some mean P_DE based on disk structure and gas content seems reasonable. This explanation is not unique, but it does offer a logical scenario to explain what we see.

Mean Density, Dense Gas Fraction, and Star Formation: While ${I}_{\mathrm{HCN}}/{I}_{\mathrm{CO}}$ increases with Σ_*, ${I}_{\mathrm{CO}}\sim {{\rm{\Sigma }}}_{\mathrm{mol}}$, and P_DE, Σ_SFR/I_HCN decreases. That is, as the dense gas fraction, and so likely the mean density, increases, the apparent rate of star formation per unit dense gas decreases. Taken at face value, this implies that any particular density—here specific n_eff of HCN—has a context-dependent role in star formation.

A context-dependent role for density in star formation could be expected if clouds are always marginally bound or virialized but vary in their mean density. Following, e.g., Krumholz & Kruijssen (2015) and many others, overdensities will collapse out of the turbulent cascade and form stars. This immediately star-forming material seems likely to exhibit a power-law distribution of densities (see, e.g., Kainulainen et al. 2009; Federrath & Klessen 2013). In clouds with high mean density, the immediate condition for these star-forming structures may be shifted to higher densities (e.g., Kruijssen et al. 2014). That is, stars should form out of overdensities in turbulent molecular clouds.

If HCN always traces the same gas density, n_HCN, and if that density does not always capture the collapsing, power-law tail of the density distribution, then we could expect that the ratio of star formation per unit HCN will change as the density distribution changes. In regions of low P_DE and low mean cloud density, HCN traces high-density gas more immediately related to star formation. In high-P_DE regions with high mean density, HCN may trace a large fraction of the bulk mass in the cloud.

As above, this is not a unique interpretation. The mapping between HCN emission and dense gas may change with environment. And this interpretation requires that HCN not always sample a self-similar power law tail of the density distribution that behaves the same in any environment (see, e.g., Leroy et al. 2017b). There is support for a scenario like what we describe from the Milky Way center (Kepley et al. 2014; Rathborne et al. 2015), but much more work remains to be done.

Turbulence or Pressure? We have emphasized mean density in this section, but note that the Mach number also plays a key role in the density distribution for a turbulent cloud (Padoan & Nordlund 2002; Krumholz & McKee 2005). Usero et al. (2015) showed that, using the model of Krumholz & Thompson (2007), plausible variations in the Mach number could explain many of their observations, along with those of García-Burillo et al. (2012) and Gao & Solomon (2004). Observations suggest that as the mean density changes, the turbulent velocity dispersion will also change (see, e.g., Leroy et al. 2016), so many of the same considerations described above will apply. We intend to return to this topic using high-resolution CO observations of our targets.

Resolved maps of density-sensitive line ratios across galaxies are new. We have focused on scaling relations that relate closely to observables. Moreover, we placed particular emphasis on HCN, the brightest high critical density millimeter-wave line, and we use our other lines to investigate the robustness of our physical interpretation. We also target a limited area and set of targets designed to highlight the contrast between galaxy disks and central starburst regions. In this section, we highlight some of the caveats and logical next steps related to these choices.

First, our survey does not yet offer a representative sample. It was designed to capture variations in Σ_SFR/Σ_mol rather than to span whole galaxy disks or the local galaxy population. Thus, our scaling relations should be taken as indicative rather than authoritative. In particular, we do not yet know how they would extend to the low surface density, outer parts of galactic disks. We also lack sizable samples of barred and unbarred galaxies to quantitatively estimate the impact of a strong bar. We expect the full EMPIRE sample to help with this, as well as future deeper, wide-area mapping by the GBT and ALMA.

Second, we adopt a simple translation from HCN and CO to dense and total molecular gas mass. Given the novelty of wide-area HCN mapping, we emphasize reporting the basic scaling relations. We also checked that our other dense gas tracers yield qualitatively similar results. Moving forward, a next direction for this field is clearly to combine a large suite of lines to model the density distribution; see Leroy et al. (2017b) for a first step in this direction, following Krumholz & Thompson (2007) and Narayanan et al. (2008). Just as important, higher J transitions are needed to understand the impact of the excitation effect. Optically isotopologue observations are needed to constrain the optical depth of the dense gas tracers, which directly impacts their n_eff (see Jiménez-Donaire et al. 2017a). Extragalactic observations lack an external constraint on the dense gas mass (such as dust observations can provide in local clouds). As a result, wide-area Milky Way mapping of diverse environments (Pety et al. 2017) has a crucial role to play in informing the translation from line ratios to the density distribution.

We also adopt a simplified treatment of α_CO in order to remain close to observables. Based on Sandstrom et al. (2013), who studied some of our targets, we know that this represents an oversimplification. α_CO often drops toward the bright centers of barred galaxies. However, as discussed by Usero et al. (2015), implementing α_CO corrections with no knowledge of α_HCN variations may obscure the basic scaling relations. Still, our targets have also been mapped by ALMA with the explicit goal of deriving α_CO from resolved multiline modeling and dust comparison. We expect to revisit this topic in detail in future work.

Other physical parameter estimates may also induce uncertainty. We adopt a combination of Hα and 24 μm to trace star formation. Much of the starburst literature (Gao & Solomon 2004; García-Burillo et al. 2012) has focused on bolometric infrared emission as a star formation tracer, motivating us to explore the impact of adopting other SFR tracers (e.g., IR tracer) in appendix. While there are subtle differences in using TIR, UV and 24 μm, or Hα and 24 μm, none of them appear to change our core results.

Finally, we make frequent reference to the structure of turbulent clouds, discussing the relationship between mean density, dynamical state (i.e., virial parameter), and the density distribution probed by our spectroscopy. Our galaxies are also targets of high-resolution (θ ∼ 1'') ALMA CO (2–1) imaging ("PHANGS-ALMA"; A. K. Leroy et al. 2018, in preparation). A next main direction for these kinds of observations will be to directly compare the turbulent structure of clouds at ∼50–100 pc scales from high-resolution imaging to the density distribution inferred from spectroscopy. The same data set will allow detailed comparison of the structure of the nuclear region and gas flow along the bars to the density distribution.