# Spatially Resolved Dust, Gas, and Star Formation in the Dwarf Magellanic Irregular NGC 4449

D. Calzetti, G. W. Wilson, B. T. Draine, H. Roussel, K. E. Johnson, M. H.

Heyer, W. F. Wall, K. Grasha, A. Battisti, J. E. Andrews, et al.

► To cite this version:

D. Calzetti, G. W. Wilson, B. T. Draine, H. Roussel, K. E. Johnson, et al. Spatially Resolved Dust, Gas, and Star Formation in the Dwarf Magellanic Irregular NGC 4449. The Astrophysical Journal, 2018, 852, 10.3847/1538-4357/aaa1e2. insu-03747718

HAL Id: insu-03747718

https://insu.hal.science/insu-03747718v1

Submitted on 19 Aug 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Spatially Resolved Dust, Gas, and Star Formation in the Dwarf Magellanic Irregular NGC 4449*

D. Calzetti<sup>1</sup> <del>1</del>, G. W. Wilson<sup>1</sup> <del>1</del>, B. T. Draine<sup>2</sup> <del>1</del>, H. Roussel<sup>3</sup> <del>1</del>, K. E. Johnson<sup>4</sup> <del>1</del>, M. H. Heyer<sup>1</sup>, W. F. Wall<sup>5</sup> <del>1</del>, K. Grasha<sup>1</sup> <del>1</del>, A. Battisti<sup>1</sup> <del>1</del>, J. E. Andrews<sup>6</sup>, A. Kirkpatrick<sup>7</sup> <del>1</del>, D. Rosa González<sup>5</sup> <del>1</del>, O. Vega<sup>5</sup> <del>1</del>, J. Puschnig<sup>8</sup> <del>1</del>, M. Yun<sup>1</sup> <del>1</del>, G. Östlin<sup>8</sup> <del>1</del>, A. S. Evans<sup>4,9</sup> <del>1</del>, Y. Tang<sup>1</sup> <del>1</del>, J. Lowenthal<sup>10</sup> <del>1</del>, and D. Sánchez-Arguelles<sup>5</sup> <del>1</del>

Received 2017 May 17; revised 2017 November 28; accepted 2017 December 12; published 2018 January 12

## Abstract

We investigate the relation between gas and star formation in subgalactic regions, ~360 pc to ~1.5 kpc in size, within the nearby starburst dwarf NGC 4449, in order to separate the underlying relation from the effects of sampling at varying spatial scales. Dust and gas mass surface densities are derived by combining new observations at 1.1 mm, obtained with the AzTEC instrument on the Large Millimeter Telescope, with archival infrared images in the range 8–500 µm from the *Spitzer Space Telescope* and the *Herschel Space Observatory*. We extend the dynamic range of our millimeter (and dust) maps at the faint end, using a correlation between the far-infrared/millimeter colors $F(70)/F(1100)$ (and $F(160)/F(1100)$) and the mid-infrared color $F(8)/F(24)$ that we establish for the first time for this and other galaxies. Supplementing our data with maps of the extinction-corrected star formation rate (SFR) surface density, we measure both the SFR-molecular gas and the SFR-total gas relations in NGC 4449. We find that the SFR-molecular gas relation is described by a power law with an exponent that decreases from ~1.5 to ~1.2 for increasing region size, while the exponent of the SFR-total gas relation remains constant with a value of ~1.5 independent of region size. We attribute the molecular law behavior to the increasingly better sampling of the molecular cloud mass function at larger region sizes; conversely, the total gas law behavior likely results from the balance between the atomic and molecular gas phases achieved in regions of active star formation. Our results indicate a nonlinear relation between SFR and gas surface density in NGC 4449, similar to what is observed for galaxy samples.

**Key words:** dust, extinction – galaxies: dwarf – galaxies: general – galaxies: individual (NGC 4449) – galaxies: ISM – galaxies: star formation

## 1. Introduction

The evolution of the baryonic component of galaxies across cosmic times is inextricably linked to the relation between star formation and the gas reservoir. Large observational and theoretical efforts have therefore been expended on both measuring this relation and understanding its physical underpinning.

Star formation and gas reservoirs are related by a tight scaling (the Schmidt-Kennicutt law, or SK law; Kennicutt 1998), which holds at both low and high redshifts, when these quantities are averaged over entire galaxies (e.g., Kennicutt & Evans 2012; Daddi et al. 2010; Genzel et al. 2010, 2015). The same relation breaks down when zooming into smaller regions, and recently formed stars appear spatially uncorrelated with molecular clouds over scales smaller than ~100–200 pc (e.g., Momose et al. 2010; Onodera et al. 2010; Schruba et al. 2010). The breakdown may be due to the onset of two effects: (1) the increasing scatter in tracers of both star formation rate (SFR)

and gas clouds, due to small number statistics in the random sampling (Calzetti et al. 2012; Kruijssen & Longmore 2014); and (2) the finite timescale of the physical association between young stars and their natal cloud (e.g., Kawamura et al. 2009; Whitmore et al. 2014). Within molecular clouds, the dust-enshrouded star formation is most closely associated with the densest gas components (e.g., Gao & Solomon 2004; Lada et al. 2010, 2012; Evans et al. 2014), while it appears only indirectly related to the molecular gas reservoir (e.g., Heiderman et al. 2010; Guttermuth et al. 2011; Faesi et al. 2014; Hony et al. 2015). At the kiloparsec or larger scale, where multiple clouds are averaged together, the relation becomes better defined, but its exact functional form, which carries physical information, remains controversial (Kennicutt et al. 2007; Bigiel et al. 2008; Rahman et al. 2012; Lada et al. 2013; Leroy et al. 2013; Shetty et al. 2014a). The ill-quantified CO-dark gas in galaxies (Pineda et al. 2013) and the presence of non-star-forming, diffuse molecular gas (Shetty et al. 2014a; Mogotsi et al. 2016) contribute to complicating an already-complex picture.

The interlocking of these still-disjointed pieces of evidence requires the investigation of the relation between star formation

<footnote n="1"><sup>1</sup> Department of Astronomy, University of Massachusetts-Amherst, Amherst, MA 01003, USA; calzetti@astro.umass.edu</footnote>

<footnote n="2"><sup>2</sup> Department of Astrophysical Science, Princeton University, Princeton, NJ 08544, USA</footnote>

<footnote n="3"><sup>3</sup> Institut d'Astrophysique de Paris, 75014 Paris, France</footnote>

<footnote n="4"><sup>4</sup> Department of Astronomy, University of Virginia, Charlottesville, VA 22904, USA</footnote>

<footnote n="5"><sup>5</sup> Instituto Nacional de Astrofísica, Optica, y Electrónica (INAOE), Tonantzintla, 72840 Puebla, Mexico</footnote>

<footnote n="6"><sup>6</sup> Department of Astronomy, University of Arizona, Tucson, AZ 85719, USA</footnote>

<footnote n="7"><sup>7</sup> Department of Astronomy, Yale University, New Haven, CT 06511, USA</footnote>

<footnote n="8"><sup>8</sup> Department of Astronomy, Oskar Klein Centre, Stockholm University, SE-106 91 Stockholm, Sweden</footnote>

<footnote n="9"><sup>9</sup> National Radio Astronomy Observatory, Charlottesville, VA 22903, USA</footnote>

<footnote n="10"><sup>10</sup> Department of Astronomy, Smith College, Northampton, MA 01063, USA</footnote>

<footnote>* Based on observations obtained with the Large Millimeter Telescope Alfonso Serrano—a binational collaboration between INAOE (Mexico) and the University of Massachusetts-Amherst (USA).</footnote>
and gas reservoir across a range of spatial scales, probing a variety of galactic environments and star formation activity, from the molecular cloud size regions all the way to whole galaxies. Liu et al. (2011) derive the SK law as a function of spatial scale for two late-type star-forming spirals: NGC 5194 (M51) and NGC 3521, covering the range from ~250 pc to 1.3 kpc. When expressing the SK law as $\Sigma_{SFR} \propto \Sigma_{H_2}^\gamma$, with $\Sigma$ indicating the surface density and $H_2$ indicating that only the molecular phase is included, Liu et al. find that the exponent $\gamma_{H2}$ is a decreasing function of increasing spatial scale, converging toward unity beyond ~1–1.3 kpc. They also find that the scatter of the data about the best-fitting relation decreases for increasing spatial scale. Wall et al. (2016) measure $\gamma_{H2}$ in both NGC 5194 and NGC 5236 and find values that are in agreement with those of Liu et al. (2011), once the differences in the fitting routines are factored in (see Calzetti et al. 2012, for a quantification of the impact of fitting routine choices).

Calzetti et al. (2012) model these results as the effect of stochastic sampling of the molecular clouds’ mass function and of star formation within galaxy regions. At region sizes larger than ~1 kpc, those authors find that the cloud mass function becomes fully sampled, and the relation between $\Sigma_{SFR}$ and $\Sigma_{H2}$ converges to an exponent $\gamma_{H2} = 1$ irrespective of the underlying, intrinsic relation between SFR and gas reservoir. The exact size of the galactic region where the convergence to $\gamma_{H2} = 1$ occurs depends on the maximum size of the molecular clouds formed in the galaxy, with region sizes in the range of ~1–2 kpc for $M_{cloud}(max) \sim 10^6–10^{7.5} M_{\odot}$. Similar results, at least qualitatively, are obtained by Kruijssen & Longmore (2014), where they model their galaxy regions by centering on either a star-forming region or a molecular cloud; they are, however, quantitatively different from those of Calzetti et al. (2012), who sample their regions randomly within galaxies. Both studies agree that once the full sampling of the cloud mass function is achieved, the relation between $\Sigma_{SFR}$ and $\Sigma_{H2}$ is only “counting clouds”; the physical information is contained in the trends of $\gamma_{H2}$ and the scatter of the data as a function of the spatial scale of galaxy regions.

The literature on the analysis of the observed SFR–gas relation as a function of region size is currently sparse; the relevant studies have so far concentrated on large, star-forming spiral galaxies (see above) and, albeit not in a systematic manner, on the Magellanic Clouds (Hony et al. 2015; Jameson et al. 2016). In order to expand on those studies, we have recently mapped the nearby starburst dwarf NGC 4449 with the Large Millimeter Telescope, using the AzTEC instrument at 1.1 mm, which we combine with imaging data at shorter wavelengths to investigate the relation between SFR and gas densities at subkiloparsec and kiloparsec scales.

NGC 4449 is a well-studied local Magellanic irregular dwarf galaxy hosting a central starburst. Its basic parameters are listed in Table 1. The specific SFR (SFR/stellar mass) places NGC 4449 solidly, about a factor of 2.5–9, above the main sequence for local star-forming galaxies, with the exact offset depending on the main-sequence relation used (Whitaker et al. 2012; Cook et al. 2014). The starburst has been possibly triggered by a minor merger (Lelli et al. 2014) or by the interaction with another galaxy (Hunter et al. 1998). It has subsolar metallicity with a modest gradient (Table 1 of Pilyugin et al. 2015). The central oxygen abundance value determined by Pilyugin et al. (2015) is in agreement with that of Berg et al.

<table><caption>**Table 1**

Adopted Parameters for NGC 4449</caption><thead><tr><td>Parameter (Units)</td><td>Value</td><td>Reference<sup>a</sup></td></tr></thead><tbody><tr><td>Distance (Mpc)</td><td>4.2</td><td>1, 2</td></tr><tr><td>Inclination (deg)</td><td>68</td><td>3, 4</td></tr><tr><td>Stellar mass ($M_{\odot}$)</td><td>1 × 10<sup>9</sup></td><td>5</td></tr><tr><td>H I mass ($M_{\odot}$)<sup>b</sup></td><td>2 × 10<sup>9</sup></td><td>6</td></tr><tr><td>SFR ($M_{\odot}$ yr<sup>-1</sup>)<sup>c</sup></td><td>0.5</td><td>7</td></tr><tr><td>12 + log(O/H)<sup>d</sup></td><td>8.26</td><td>8</td></tr><tr><td>Metall. gradient (dex kpc<sup>-1</sup>)<sup>d</sup></td><td>-0.055</td><td>8</td></tr></tbody></table>

**Notes.**

<sup>a</sup> References: (1) Karachentsev et al. 2003 (from TRGB); (2) Tully et al. 2013 (from TRGB); (3) Hunter et al. 2002; (4) Hunter et al. 2005; (5) Calzetti et al. 2015; (6) Huchtmeier & Richter 1989; (7) Lee et al. 2009; (8) Pilyugin et al. 2015.

<sup>b</sup> Atomic hydrogen mass associated with the galaxy.

<sup>c</sup> Star formation rate from the extinction-corrected ultraviolet luminosity.

<sup>d</sup> Central oxygen abundance and abundance gradient, respectively. The central oxygen abundance has an uncertainty of ±0.01 (Pilyugin et al. 2015).

(2012), 12 + log(O/H) = 8.26 ± 0.09. Our analysis concentrates on the inner ~5 kpc, for which, when necessary, we adopt a mean value of the oxygen abundance 12 + log(O/H) ~ 8.2, or 30% of the solar value.<fnref n="11" />

The low metallicity of NGC 4449 is consistent with its relatively low dust content. The IR/UV ratio indicates that only about 40% of the light from young stars is absorbed by dust in this galaxy (Hao et al. 2011; Grasha et al. 2013). Reines et al. (2008) find that, with one exception, even embedded sources have $A_V < 1.5$ mag. This characteristic enables us to effectively decouple tracers of SFR from those of gas content. In dusty galaxies, the use of the blue and red sides of the IR spectral energy distribution (SED) to derive SFR and dust (gas) masses, respectively, can potentially introduce a covariance between the two physical quantities (Rojas Bolivar et al. 2016). The low extinction of NGC 4449 ensures that most of the SFR will be traced at UV–optical (e.g., $H\alpha$) wavelengths; this decouples the SFR measurements from those of the dust mass, the latter used to derive gas masses. Despite being relatively metal deficient, NGC 4449 is still sufficiently metal-rich that measurements of dust masses yield reliable measurements of gas masses (Draine et al. 2007; Rémy-Ruyer et al. 2014; Accurso et al. 2017).

The properties listed above, together with the relative proximity of the galaxy (10" subtends a linear scale of 204 pc), make NGC 4449 an ideal target for investigating the relation between SFR and gas surface densities as a function of spatial scale.

## 2. Observations and Archival Data

### 2.1. LMT Data

#### 2.1.1. Observations

Observations of NGC 4449 (R.A.(2000) = 12h28m13s60, decl.(2000) = +44°05′35″14) were carried out with the AzTEC instrument on the 32 m Large Millimeter Telescope near the city of Puebla (Mexico) between 2015 February 25 and March 9, in excellent weather conditions with a median

<footnote n="11">¹¹ We adopt a solar oxygen abundance of 12 + log(O/H) = 8.69 (Asplund et al. 2009).</footnote>
<table><caption>**Table 2**
Infrared and Millimeter Flux Densities for NGC 4449</caption><thead><tr><th>Wavelength (μm)</th><th>Facility/<sup>a</sup> Instrument</th><th>PSF<sup>b</sup> (arcsec)</th><th>Galaxy<sup>c</sup> Flux Density (Jy)</th><th>Central<sup>d</sup> Flux Density (Jy)</th></tr></thead><tbody><tr><td>8</td><td>SST/IRAC</td><td>2.8</td><td>1.41 ± 0.17</td><td>1.35 ± 0.16</td></tr><tr><td>24</td><td>SST/MIPS</td><td>6.4</td><td>3.04 ± 0.34</td><td>2.93 ± 0.33</td></tr><tr><td>70</td><td>HSO/ PACS</td><td>5.7</td><td>51.2 ± 6.8</td><td>47.3 ± 6.2</td></tr><tr><td>100</td><td>HSO/PACS</td><td>7.0</td><td>74.4 ± 9.1</td><td>68.5 ± 8.6</td></tr><tr><td>160</td><td>HSO/PACS</td><td>11.2</td><td>74.2 ± 8.9</td><td>61.3 ± 7.4</td></tr><tr><td>250</td><td>HSO/SPIRE</td><td>18.2</td><td>34.3 ± 4.1</td><td>25.6 ± 3.1</td></tr><tr><td>350</td><td>HSO/SPIRE</td><td>24.9</td><td>15.5 ± 2.0</td><td>11.3 ± 1.5</td></tr><tr><td>350</td><td>Planck</td><td>278.</td><td>15.5 ± 0.4</td><td>...</td></tr><tr><td>500</td><td>HSO/SPIRE</td><td>36.1</td><td>6.0 ± 0.7</td><td>4.2 ± 0.5</td></tr><tr><td>550</td><td>Planck</td><td>290.</td><td>4.6 ± 0.1</td><td>...</td></tr><tr><td>850</td><td>Planck</td><td>292.</td><td>1.24 ± 0.09</td><td>...</td></tr><tr><td>1100</td><td>LMT/<br>AzTEC<sup>e</sup></td><td>8.5</td><td>0.38 ± 0.16</td><td>0.25 ± 0.09</td></tr><tr><td>1380</td><td>Planck</td><td>301.</td><td>0.25 ± 0.04</td><td>...</td></tr></tbody></table>

**Notes.**

zenith opacity, $\tau_{225}$ GHz, of 0.055. AzTEC is a 144-element bolometer array configured to observe in the 1.1 mm atmospheric window (225 GHz; Wilson et al. 2008), with a half-power beam width of 8"5 (corresponding to a beam area of 164 arcsec²; Table 2). The total observing time was 9.03 hr, divided into 23 individual maps. Each map consisted of a combination of a fast motion of the telescope boresight in a Lissajous pattern superimposed on a slow raster scan of the center of the Lissajous pattern over the field. Both motions were performed in the azimuth-elevation coordinate system. This observational technique, known colloquially as a “rastajous scan,” allows for deep imaging of moderate-sized fields with high scan speeds and good cross-linking of pixels but without the inefficiencies of lost observing time due to rapid turnarounds of the telescope during fast raster scan.

Corrections to the LMT’s native pointing model were made by interspersing the maps with pointed observations of the quasar 1203+480, which has an unresolved flux of 190 mJy. Pointing corrections were typically ±3" in azimuth and -4" to 11" in elevation. These corrections are typical for AzTEC observations owing to a known boresight mismatch between the AzTEC and the LMT’s pointing model. Reconstructed

pointing errors following these corrections are typically of order 1", which is much smaller than the beam size (8"5).

Data were calibrated based on beam map observations of the asteroids Ceres and Pallas taken during each evening of observations along with observations of the pre-planetary nebula CRL 618. All calibrations are self-consistent, and we estimate the calibration error to be approximately 10%.

### 2.1.2. Data Processing and Final Maps

The data were processed using the recently developed AzTEC C++ pipeline, which is now the standard pipeline for AzTEC data on the LMT. This new pipeline offers the same reduction processes reported in previous AzTEC studies (e.g., Perera et al. 2008; Scott et al. 2008, 2012; Austermann et al. 2010; Liu et al. 2010), including removal of atmospheric contributions (the dominant noise source in the AzTEC maps), correction of pointing errors, data calibration, and co-addition of the ensemble of maps, but with improvements to achieve better performance and speed. The updated pipeline also implements a new approach to estimate the contribution of the atmosphere into the observed raw time streams, aiming at improving the recovery of extended emission in the AzTEC data.

The new approach is based on the Cottingham method (Cottingham 1987), where an unbiased estimator of the atmosphere signal is constructed using a set of B-spline basis vectors and the pointing information from individual detectors. Hincks et al. (2010) demonstrate that this method is a maximum likelihood estimator for both the atmosphere signal and the astronomical surface brightness distribution. Our implementation has been adapted to maximize the compatibility with the existing pipeline code and other atmosphere subtraction techniques; therefore, it only estimates and subtracts the atmosphere component using Equation (7) from Hincks et al. (2010). Details on the technique are given in the Appendix.

In addition to a surface brightness map, a noise map is constructed from jackknifed noise realizations of each AzTEC map in which the time stream is randomly multiplied by ±1. The effect of the jackknife step is to remove sources while retaining noise properties. Each jackknifed time stream is converted into a map and co-added in weighted quadrature to produce a final image of rms noise values, $\sigma$.

The primary data products from the AzTEC pipeline are an image of surface brightness at 1.1 mm, in units of Jy beam⁻¹, and a corresponding weight image (1/σ²) that reflects both instrumental and atmospheric contributions to the noise budget. The final images are constructed in equatorial coordinates (J2000). As in all such maps, the outer edges have larger errors owing to fewer numbers of collected samples at these positions, and therefore less accumulated time, relative to the central regions. The map has a diameter of about 600" (~12 kpc at the distance of NGC 4449), with a lower noise region of about 450" diameter (Figure 1). Our galaxy fits well within this central region: NGC 4449 has a 500 µm diameter <450", as measured from the Herschel SPIRE map (see next section), and an optical major axis of ~370".<fnref n="12" /> The central 450" diameter region has a median sensitivity of 1 mJy beam⁻¹, with 1σ variations of less than 10%.

<footnote><sup>a</sup> Facility/instrument combination: SST = Spitzer Space Telescope, IRAC and MIPS Cameras; HSO = Herschel Space Observatory, PACS and SPIRE instruments; Planck = Planck Observatory; LMT = Large Millimeter Telescope with the AzTEC instrument.</footnote>

<footnote><sup>b</sup> Point-spread function (PSF) FWHM. The SST and HSO PSFs are from Aniano et al. (2011); the Planck PSFs (effective FWHMs) are from the Planck Collaboration, as listed at https://wiki.cosmos.esa.int/planckpla/index.php/Effective_Beams.</footnote>

<footnote><sup>c</sup> Flux density of the entire galaxy. Except for the Planck photometry, all measurements are performed on images matched in angular resolution to the HSO/SPIRE 500 µm PSF. The Planck photometry is retrieved from the Planck Legacy Archive at https://www.cosmos.esa.int/web/planck/pla.</footnote>

<footnote><sup>d</sup> Flux density of a central, rectangular region with size ~4.7 × 4.7 kpc², corresponding to a deprojected size of ~4.7 × 12.5 kpc² (Figure 2). All measurements are performed on images matched in angular resolution to the HSO/SPIRE 500 µm PSF.</footnote>

<footnote><sup>e</sup> The LMT/AzTEC flux density obtained by including only regions with S/N ≥ 3.5 in the map is 0.078 ± 0.004 Jy. The sum area of these regions is ~10% of the area of the central region, or 1.46 × 1.46 kpc² (corresponding to a deprojected size of 1.46 × 3.90 kpc²).</footnote>

<footnote n="12">¹² From NED, the NASA Extragalactic Database.</footnote>
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Side-by-side comparison between a three-color composite (left) and the AzTEC 1.1 mm emission (right) of NGC 4449. The three colors in the left panel are GALEX far-ultraviolet (blue), SST/MIPS 24 µm (green), and HSO/PACS 70 µm (red) images, all at the native resolution of ~5"-6".5. The size of each panel is about 480", corresponding to about 9.8 kpc at the distance of the galaxy. The cyan circle has a diameter of 450", i.e., the size of the photometric aperture we use to measure the total flux from the galaxy. The yellow square marks the location and size of the central region analyzed in detail in this work. The dynamic range of the AzTEC image shown in this figure is [-0.5, +0.5] MJy sr⁻¹, and the size of the AzTEC beam is shown about 30% larger than actual, for display purposes. North is up, east is left.</figcaption></figure>

The measured fluxes for NGC 4449 at 1.1 mm (=1100 µm) are listed in Table 2. The AzTEC bandpass excludes the CO (2–1) line, so no contamination is expected from this transition. Estimates of the free–free and synchrotron emission by extrapolating the measurements at centimeter wavelengths for this galaxy<fnref n="13" /> yield a ≲8% contribution of these processes to the light in the AzTEC band. Our estimate uses spectral indices of −0.5 for the synchrotron emission and −0.1 for the free–free emission, from Srivastava et al. (2014). Thus, the measured 1100 µm flux can be attributed almost entirely to dust emission.

## 2.2. Archival Images

For our analysis we make use of several images and maps from archives, covering a wide range of wavelengths: ground-based Hα and H I, Spitzer Space Telescope (SST) 3.6–24 µm, and Herschel Space Observatory (HSO) 70–500 µm.

The ground-based Hα image was obtained at the Bok Telescope as part of the LVL (Local Volume Legacy) project, a Spitzer Legacy program that has observed almost 300 galaxies in the local 11 Mpc (Dale et al. 2009; Lee et al. 2011). The image is available through NED, already continuum-subtracted and calibrated. We correct the fluxes for the effects of foreground Milky Way (MW) extinction (E(B − V) = 0.017, corresponding to A_Hα = 0.043 mag), and of [N II] contamination ([N II]/Hα = 0.23; Kennicutt et al. 2008). The final image has a point-spread function (PSF) FWHM ~ 1"5.

The SST IRAC 3.6 µm, 4.5 µm, 5.8 µm, and 8.0 µm and MIPS 24 µm, 70 µm, and 160 µm images were also obtained as part of the LVL project (Dale et al. 2009) and are available already fully processed and calibrated through IRSA and NED.

All Spitzer images are in units of MJy sr⁻¹. For our pixel-resolved analysis, we make use of the highest angular resolution images from Spitzer, i.e., up to and including the 24 µm image. For the purpose of fitting the dust SED, we remove the stellar contribution from both the 8 and 24 µm images, using the formulae of Helou et al. (2004) and Calzetti et al. (2007; see also Calapa et al. 2014): $f_{\nu,D}(8) = f_{\nu}(8) - 0.25 f_{\nu}(3.6)$, and $f_{\nu,D}(24) = f_{\nu}(24) - 0.035 f_{\nu}(3.6)$, where the flux densities are in units of Jy, and the subscript “D” indicates the dust-only emission component. Using these formulae, the typical stellar continuum contribution to the fluxes at 8 and 24 µm of individual pixels is about 6%–7% and <1%, with maximum values of about 20% and 2.4% at 8 and 24 µm, respectively. The contamination of the 3.6 µm image by the 3.3 µm polycyclic aromatic hydrocarbon (PAH) emission feature is minimal (~5%–15%; Meidt et al. 2012), and the IRAC 3.6 can be used as a reasonable stellar continuum tracer for the purpose of subtraction from the 8 and 24 µm fluxes. When required by the analysis (see later sections), we will instead use the uncontaminated 4.5 µm IRAC image as a reliable tracer of the stellar continuum distribution.

The HSO PACS 70 µm, 100 µm, and 160 µm and SPIRE 250 µm, 350 µm, and 500 µm images were retrieved from the Herschel archive and reprocessed by us using Scanamophos v24.0 (Roussel 2013), which dramatically improves the final maps by removing the low-frequency noise associated with them. The PACS images are in units of Jy pixel⁻¹, and the SPIRE images are in units of Jy beam⁻¹. For photometry, we use the following SPIRE beam areas: 469 arcsec², 831 arcsec², and 1804 arcsec² for 250 µm, 350 µm, and 500 µm, respectively.<fnref n="14" />

<footnote n="13">¹³ Photometry at centimeter wavelengths is from NED, the NASA Extragalactic Database. See references therein.</footnote>

<footnote n="14">¹⁴ See http://herschel.esac.esa.int/twiki/bin/view/Public/SpirePhotometerBeamProfile2.</footnote>
Color corrections to the *SST* and *HSO* photometry for the typical SEDs of star-forming galaxies are at the level of 10%–12% for the *HSO/PACS70* band and less than 5% for all other *HSO* and *SST* bands. We apply color corrections only to the *HSO/PACS70* photometry, as the corrections for the other bands are small relative to the measurement uncertainties. In addition, we implement aperture corrections of about 8% to the *SST/IRAC* 8 µm, as appropriate for partially extended sources, to compensate for the scattered light into this band.<fnref n="15" />

HI maps of NGC 4449 were obtained as part of the THINGS (The HI Nearby Galaxy Survey; Walter et al. 2008) project, which observed about three dozen nearby galaxies with the NRAO Very Large Array.<fnref n="16" /> The highest angular resolution velocity-integrated HI maps for NGC 4449, obtained with the “robust weighted data” (Walter et al. 2008), have a beam of 13′′74 × 12′′50 and are calibrated in units of Jy beam⁻¹ m s⁻¹.

## 2.3. Additional Processing and Photometry

After verifying the alignment in sky coordinates, we convolve all the images to common resolution, using the kernels of Aniano et al. (2011). We then resample all images to common pixel sizes that we refer to as “spaxels,” since a crude spectrum made of multiwavelength photometry is associated with each resampled pixel. The spaxels are used for the spatially resolved analysis. We select two resolutions for our photometry: that of the SPIRE 500 µm images for large-scale analysis, and that of the PACS 160 µm images for the spatially resolved analysis (Table 2). The latter choice is dictated by a compromise between exploiting the excellent angular resolution of the AzTEC images and preserving the largest possible number of photometric data points for obtaining reliable fits of the dust SEDs.

The two resolutions dictate the choice of photometric measurements included in each dust SED fit: for the large-scale fits we include all dust measurements between 8 and 1100 µm, while for the spaxel-based fits we include only measurements at 8 µm, 24 µm (from the *SST*), 70 µm, 100 µm, 160 µm (from the *HSO*), and 1100 µm (from the LMT). Table 2 lists, for each of the dust emission bands used in this work, the central wavelength, the facility/instrument combination that has provided the measurement, the PSF’s FWHM, the whole-galaxy photometry,<fnref n="17" /> and the photometry of a central rectangular region with projected size 4.7 × 4.7 kpc² and deprojected size 4.7 × 12.5 kpc² for a 68° inclination (Hunter et al. 2002, 2005; Table 1 and Figure 2).

The central rectangular region is chosen to include all AzTEC pixels with S/N > 3.5 that also have PACS 160 µm values with S/N > 10, when the data are binned in spaxels of 11" squares. Not coincidentally, this is also very close to the region imaged by Kohle et al. (1998) and Kohle (1999) in CO

(1–0) (see Figure 1 in Böttner et al. 2003). In Section 6.2, we use this CO measurement to derive an *M(H₂)–L(CO)* conversion factor for the central region.

We compare our photometry in the AzTEC band with similar measurements available either from archives or from the literature. For the whole galaxy, we use the *Planck* photometry at the position of NGC 4449 (Planck Collaboration 2016a), as retrieved from the *Planck* Legacy Archive. We employ the archival photometry extracted using an aperture with Gaussian weights and with axes of about 5', close to the FWHM of the effective beam in the four bands at wavelengths between 350 and 1380 µm. The *Planck* photometry is listed in Table 2. Interpolation between the two *Planck* bands at 850 and 1380 µm predicts a flux of 0.53 Jy at 1100 µm for this galaxy, implying that our recovered flux of ~0.4 Jy (Table 2) is a little over 70% of the total. Although the *Planck* data may suffer from a modest flux increase owing to potential foreground and background source contamination (Dale et al. 2017), our AzTEC observations are below the *Planck* interpolated flux, although not significantly when uncertainties are included. Some of the discrepancy may be due to the fact that we do not recover all the extended millimeter emission from the galaxy. This is expected, due to the large spatial scale filtering applied to the millimeter data in order to remove contributions from the atmosphere (see Section 2.1); this process also filters out extended, low surface brightness emission from the galaxy and is a common outcome of all ground-based millimeter and submillimeter camera systems. The size of the AzTEC map, which is only slightly larger than the galaxy itself, may also contribute to underestimating the total flux. We further verify this hypothesis in Section 5, where we model the low-S/N spaxels in the AzTEC image. For the central region, we compare our AzTEC photometry, which yields a flux density of 0.25 ± 0.09 Jy (Table 2), with the value of 0.26 ± 0.04 Jy obtained by Kohle (1999) at 1.2 mm for roughly the same region, after removal of the CO (2–1) line contamination (Böttner et al. 2003). The two values are identical within the uncertainties.

We note that while the *HSO/PACS100* and *HSO/PACS160* whole-galaxy photometries are virtually identical to the IRAS100 (Hunter et al. 1989) and *SST/MIPS160* photometries, respectively, the *HSO/PACS70* is ~30% higher than the *SST/MIPS70* measurement. The discrepancy between MIPS70 and PACS70 has already been reported for other galaxies (Aniano et al. 2012; Draine et al. 2014), with striking variations observed in M31, where the ratio between the two bands covers the range 0.5–2 across the body of the galaxy, with no obvious correlation with any specific property (Draine et al. 2014). Figure 3 shows the ratio between the *HSO/PACS70* and *SST/MIPS70* images for NGC 4449, in spaxels 20" in size (comparable to the FWHM of the MIPS70 PSF; Aniano et al. 2011). We observe spaxel-to-spaxel variations in the PACS70/MIPS70 ratio between ~0.9 and 1.8, about half as large as those observed in M31. We use the PACS70 data as our default in what follows, so as to leverage the higher angular resolution of the *HSO*, but bearing in mind this potential difficulty. In Section 6.1, we also assess the impact on the inferred dust masses of using the *SST/MIPS70* photometry instead of the *HSO/PACS70* one.

For the spaxel-based analysis, we make use of the HI map at the original resolution, since the FWHM of the HI beam is close to the FWHM of the Herschel PACS 160 µm PSF.

<footnote n="15">¹⁵ The sources in our field of view have a typical extent of 4"–8", and we apply the aperture corrections described in Section 4.11 of the IRAC Instrument Handbook (http://irsa.ipac.caltech.edu/data/SPITZER/docs/irac/iracinstrumenthandbook/1/).</footnote>

<footnote n="16">¹⁶ Fully processed maps are publicly available at http://www.mpi.ae/THINGS/Overview.html.</footnote>

<footnote n="17">¹⁷ For the whole-galaxy photometry at millimeter wavelengths, we measure the total flux within an aperture of diameter 360" on the AzTEC map; this is the diameter that encompasses all pixels detected with signal-to-noise ratio (S/N) ≥ 3 in the *HSO/SPIRE500* map of NGC 4449. For the other images, which are significantly larger than the extent of the galaxy, we use the curve-of-growth technique. In all cases, the image/sky background is subtracted.</footnote>
We divide the central region into 21 × 21 square spaxels with size 11" each. These are the primary resolution elements used in our spatially resolved analysis. Before binning the central region, we convolve all images to the HSO/P160 PSF (Aniano et al. 2011), with the exclusion of the SPIRE images and Planck data that are not included in the spatially resolved analysis. Out of 441 total spaxels, 429 have S/N ≥ 10 at 8, 24, 70, 100, and 160 µm. The flux distribution of the 429 spaxels in the AzTEC image is consistent with a Gaussian, as expected; 42 spaxels (~10%) have S/N ≥ 3.5, 131 have S/N ≥ 2, and 91 have negative or null 1100 µm flux. We also analyze larger spaxels with sides of 22", 33", and 44" (2, 3, and 4 times the size of the original spaxels). This binning provides information as a function of increasing spatial scale, in the range from ~360 pc to 1.5 kpc. At the largest binning size, there are 25 spaxels in the central region. We do not attempt to go beyond 1.5 kpc, as the number of spaxels becomes too small for statistical purposes.

In Section 6, we use the photometry of both the total galaxy and the central region (see above and Figures 1 and 2) from 8 to 1100 µm to perform SED fits and derive physical parameters, including dust masses. These are used to validate the procedure we utilize later to derive spaxel-based dust mass surface densities. In what follows, photometric values are expressed as $F(\lambda) = \nu F_{\nu}(\lambda)$, in units of erg s⁻¹ cm⁻², unless otherwise noted.

### 3. Star Formation Rates

We derive dust-corrected, spaxel-based SFRs using the hybrid indicator that combines Hα emission with 24 µm emission. We use the calibration of Calzetti et al. (2007), which is appropriate for H II regions and subkiloparsec star-forming regions:

$$ \text{SFR}(M_{\odot} \text{ yr}^{-1}) = 5.4 \times 10^{-42} [L(H\alpha) + 0.031 L(24)], \quad (1) $$

where both luminosities are in units of erg s⁻¹, and the SFR is calculated in each spaxel. The calibration factor in this equation assumes a Kroupa (2001) stellar initial mass function (IMF) in the stellar mass range 0.1–100 $M_{\odot}$. The SFRs are then converted to SFR surface densities, $\Sigma_{SFR}$, by dividing the SFR by the deprojected area of each spaxel. We use the inclination of 68° from Table 1. The advantage of using ionized gas and dust emission for deriving the SFR is that both ionized gas and neutral gas have the same inclination relative to the plane of the sky for this galaxy (Hunter et al. 2002). Stars have a different inclination in NGC 4449, and could be viewed face-on, if the galaxy is dominated by a strong bar (Hunter et al. 2005). This characteristic makes the use of hybrid SFR indicators that combine direct stellar light (e.g., UV) with dust emission (e.g., 24 µm emission) more complicated, and we avoid this.

At each position, the Hα light and 24 µm light may receive contributions from populations other than those associated with the local current star formation, and which may form a diffuse background to the actual light from star formation (e.g., Calzetti et al. 2005). In the case of Hα, ionizing photon leakage from massive stars at other locations (up to ≳1 kpc) can contribute to the observed ionized gas emission (commonly known as diffuse ionized gas, or DIG). Part of the 24 µm emission can be contributed by nonequilibrium heating of the dust by preexisting stellar populations, older

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Maps of the central region in NGC 4449 (the yellow rectangle in Figure 1) at 1.1 mm, from AzTEC (top row), at 24 µm, from SST/MIPS (middle left and bottom left), and at Hα, from ground-based observations (middle right and bottom right). The top and middle rows show the images after convolution to the HSO/PACS 160 µm PSF and resampling in spaxels of 11" × 11". The images in the bottom row are shown at native PSF and resolution to help orient the reader. The AzTEC map is in a linear scale, while the 24 µm and Hα maps are shown with a logarithmic stretch. The dynamic range of the AzTEC image is [−0.5, +0.9] MJy sr⁻¹. The low surface brightness areas of the 24 µm and Hα images have been depressed to highlight the structure of the bright star-forming regions. The top right panel, in heat scale, shows the location of the 11" × 11" spaxels in the AzTEC image with S/N ≥ 3.5. The central region, ~4.7 × 12.5 kpc² in deprojected size, is almost coincident with the area mapped in CO (1–0) by Kohle et al. (1998) and Kohle (1999). North is up, east is left.</figcaption></figure>

than the typical age of the current event of star formation (Calzetti et al. 2005; Draine et al. 2007; Calapa et al. 2014). The impact of these spurious contributions to $\Sigma_{SFR}$ is quantified in Calzetti et al. (2012). An estimate of the contribution of older stellar populations to the heating at 24 µm is obtained from the best-fit solution to the dust SED of the central region, presented in Section 6.2: about 30% of the emission at this wavelength is provided by the general, diffuse starlight (Section 4), which is unrelated to the current star formation.

We use GALFIT (Peng et al. 2002, 2010) to produce a smooth map of the galaxy from the 4.5 µm image; this step removes the low-level, small-scale contribution of recent star formation to the 4.5 µm emission and produces an image that traces the old, distributed stellar population (Kendall et al.
<figure><img src="image_4.png" /><figcaption>**Figure 3.** 5′3 × 5′3 region centered on NGC 4449, showing the full resolution HSO/PACS70 image in gray scale (left); the HSO/PACS70 image resampled to 20″ × 20″ pixels, i.e., roughly the SST/MIPS70 PSF size, in gray scale (middle); and the PACS70/MIPS70 ratio image, also in 20″ × 20″ pixels, in heat scale (right). The left and middle images are shown to orient the reader as to the location of the galaxy relative to the ratio image. Scale bars are shown below the middle image (in units of Jy) and the left image (adimensional). The PACS70/MIPS70 ratio is typically between 0.9 and 1.8, with a galaxy-wide ratio of ~1.3.</figcaption></figure>

2008). The GALFIT best-fit parameters include an exponen-
tial disk (Sérsic index = 1.2), with an effective radius of
~1.4 kpc and a position angle of 45°. The effective radius
(corresponding to a scale length of ~0.8 kpc) and the position
angle we recover are in agreement with the parameters of the
stellar disk obtained by Hunter et al. (1999, 2002). The
smooth stellar disk map produced by GALFIT is then
subtracted from both the 24 µm and Hα images, by removing
30% and 10% of the total emission from the central region at
24 µm and Hα, respectively (Calzetti et al. 2005, 2012;
Draine et al. 2007; Liu et al. 2011). The 30% contaminating
fraction to the 24 µm emission is discussed in the previous
paragraph. The 10% fraction of diffuse Hα emission is
consistent with the fraction of DIG present in starburst
galaxies (Oey et al. 2007). With these fractions, the number
of oversubtracted spaxels, i.e., containing zero or negative
SFR, is less than 5%. The negative/null spaxels are
distributed at the edges of the central region, in agreement
with most of the star formation being clustered toward the
center of the galaxy. After subtraction of the smooth
background emission, we combine the Hα and 24 µm
emission to obtain dust-free Σ<sub>SFR</sub> estimates at each spaxel.
We verify that in NGC 4449 most of the emission from star
formation emerges directly, instead of being reprocessed by
dust. Indeed, for 85% of the 11" × 11" spaxels, 60%–100%
of the star formation is recovered through the Hα emission
and less than 40% through the 24 µm emission. For the
remaining 15% of pixels, the light from star formation
emerges through each of the Hα emission and 24 µm
emission in roughly equal fractions.

4. Models and Fits of the Infrared
Spectral Energy Distribution

For our infrared SED fits, we employ the models of Draine
& Li (2007), as implemented in Draine et al. (2007).<fnref n="18" /> The
models consist of mixtures of dust grains that combine
carbonaceous grains (including PAHs) and amorphous

silicates, with size distributions that aim at reproducing the
MW, Large Magellanic Cloud (LMC), and Small Magellanic
Cloud (SMC) extinction curves, respectively. For each grain
distribution (extinction curve), a range of PAH dust mass
fractions is considered, between $q_{\text{PAH}} = 0.01\%$ and $4.6\%$, in
several discrete values. The highest $q_{\text{PAH}}$ value is consistent
with MW-type dust, while the SMC- and LMC-type dust has
lower $q_{\text{PAH}}$ values. In the models, the size distribution of PAHs
and the ratio of neutral to ionized species are held constant,
consistent with results that show little variation in the neutral/
ionized PAH ratio in galaxies that do not host AGNs (Smith
et al. 2007). Even in low-metallicity galaxies like the SMC,
which has ~1.5 times lower metal content than NGC 4449, the
PAH estimates based on the 8 µm photometry are in good
agreement with those based on *Spitzer*/IRS spectroscopy
(Sandstrom et al. 2010).

The dust mixture is heated by the combination of two
starlight intensity components: the diffuse starlight that
permeates the interstellar medium, described by the energy
density parameter $U_{\text{min}}$, and a range of regions with a power-
law distribution of intensities, $dM_{\text{dust}}/dU \propto U^{-\alpha}$, between
$U_{\text{min}}$ and $U_{\text{max}}$ and slope $\alpha = 2$. The two starlight intensity
components are added together in proportion to $(1 - \gamma)$ and
$\gamma$, where $0 \le \gamma \le 1$. The parameter $\gamma$ is related to the fraction
of starlight intensity due to current star formation or other
activity. A more intuitive parameter is $f_{\text{PDR}}$, the fraction of
total dust luminosity radiated from regions with $U > 10^2$, as
defined in Draine & Li (2007, Equation (29)); $f_{\text{PDR}}$ is a
measure of the dust emission from the photodissociation
regions (PDRs) of star-forming regions and is a combination
of $\gamma$, $U_{\text{min}}$, and $U_{\text{max}}$. We use $f_{\text{PDR}}$ instead of $\gamma$ in the analysis
that follows. The emission in each band will also be
proportional to the total amount of dust $M_{\text{dust}}$. Thus, for
each given extinction curve, there are a total of five
parameters for the models: $q_{\text{PAH}}$, $U_{\text{min}}$, $U_{\text{max}}$, $\gamma$ (or $f_{\text{PDR}}$), and
$M_{\text{dust}}$. The parameters increase to six if we include the
extinction curve as one of them.

For the large-scale SED fits (Section 6), we employ
photometric data in nine bands, implying three degrees of
freedom (or two if the extinction curve is included among the

<footnote n="18">¹⁸ The models are publicly available and can be retrieved from http://www.astro.princeton.edu/~draine/dust/irem.html.</footnote>
parameters to be fit). Thus, we have enough data points to avoid underconstrained fits. We implement $\chi^2$ minimization for the best fits (Draine & Li 2007) and derive the uncertainties on our best-fit quantities from the shape of the reduced $\chi^2$ probability distribution. For the $\chi^2$ minimization, we apply the standard technique of weighting each data point by the inverse of the square of the measurement uncertainty. For the spaxel-based SED fits (Section 7), we employ photometric data from six bands, implying that the fits are underconstrained (zero degrees of freedom). Thus, for the spaxels we fix several free parameters using the large-scale results for guidance.

Karczewski et al. (2013) and Rémy-Ruyer et al. (2015) derive $q_{\text{PAH}} \sim 2\%$ and $3\%$, respectively, from the global photometry of NGC 4449. These authors include SST spectroscopy in their fits, although they do not have a millimeter point. Our SED fits for the global photometry (Section 6.1) yield a best extinction curve/$q_{\text{PAH}}$ combination of MW/3.2%. Fits that use MW/2.5% and LMC/2.4% also yield acceptable results, with reduced $\chi^2$ values that are within a factor of 1.5–2 of those of the MW/3.2% parameter combination. All other combinations produce significantly larger (by factors of 3 or more) reduced $\chi^2$ values, and we consider these unacceptable. On account of the sparse sampling in both the extinction curve and $q_{\text{PAH}}$ in the models, it is difficult to assign an uncertainty to these parameters. Our best-fit values of $q_{\text{PAH}} \sim 2.4\%–3.2\%$ are in excellent agreement with those of Karczewski et al. (2013) and Rémy-Ruyer et al. (2015). The $q_{\text{PAH}}$ values from the best fits of Karczewski et al. (2013) and Rémy-Ruyer et al. (2015) are consistent with NGC 4449 being a low-metallicity galaxy, where $q_{\text{PAH}}$ is expected to be lower than the MW $q_{\text{PAH}} \sim 4.6\%$; however, NGC 4449 is still above the threshold of 12 + log(O/H) ~ 8.1, below which the PAH emission is observed to be significantly suppressed in galaxies (Draine et al. 2007; Smith et al. 2007). In summary, for the modeling of the infrared SEDs, we adopt the MW extinction curve with a fraction of 3.2% of PAHs.

# 5. Mid/Far-infrared Color Trends

## 5.1. Color Trends within NGC 4449

The dynamic range covered by the AzTEC image for the spaxels with S/N ≥ 3.5 is relatively small, about a factor of 2 (Figure 2, top right panel); furthermore, 90% of the spaxels in the central region have insufficient S/N to be considered reliable, i.e., they are below S/N = 3.5. In order to increase the number of spaxels used in our analysis, we leverage properties of the infrared to millimeter SEDs of regions in this galaxy. Figure 4 shows that the flux ratios $F(70)/F(1100)$ and $F(160)/F(1100)$ are anticorrelated with $F(8)/F(24)$ for all the high-S/N$_{1100}$ µm spaxels. This can be understood if regions of high SFR correspond to increased values of $F(24)$, $F(70)$, and $F(160)$ relative to $F(8)$ and $F(1100)$. $F(8)$ and $F(1100)$ have appreciable contributions from dust heated by the diffuse starlight, with $F(8)$ being produced by single-photon heating and $F(1100)$ coming from the larger grains with nearly steady temperatures maintained by the diffuse starlight. This is shown in Figure 5, where we derive the best SED fits for the two extreme (i.e., the lowest $F(8)/F(24)$ and the lowest $F(70)/F(1100)$) spaxels with S/N$_{1100}$ µm ≥ 3.5. The spaxel with the lowest $F(8)/F(24)$ ratio corresponds to a

region of high SFR surface density, which is contributing significantly to the heating of the dust ($f_{\text{PDR}} = 0.57$). The spaxel with the lowest $F(70)/F(1100)$ ratio corresponds to a region of very low SFR density, where the dust is entirely (within our uncertainties) heated by the general interstellar radiation field, i.e., $f_{\text{PDR}} = 0$. Figure 5 also identifies the location of both spaxels within the central region, confirming that they correspond to high and low local star formation, respectively. In summary, the trends marked by the high-S/N$_{1100}$ µm spaxels in Figure 4 identify a smooth progression from regions of high star formation, with low $F(8)/F(24)$ and high $F(70)/F(1100)$ (or $F(160)/F(1100)$) ratios, to regions of low to negligible star formation, with high $F(8)/F(24)$ and low $F(70)/F(1100)$ (or $F(160)/F(1100)$) ratios.

Comparison of the observed color–color distributions with the expectations from the Draine & Li (2007) models provides further support to the above picture (bottom two panels of Figure 4). Spaxels with low $F(8)/F(24)$ and high $F(70)/F(1100)$ (or $F(160)/F(1100)$) correspond to models with high $U_{\text{min}}$, high $U_{\text{max}}$, and high $\gamma$ values, i.e., with a high contribution to the dust heating by current star formation. Conversely, spaxels with high $F(8)/F(24)$ and low $F(70)/F(1100)$ (or $F(160)/F(1100)$) correspond to models with relatively low $U_{\text{min}}$, low $U_{\text{max}}$, and low $\gamma$ values, as expected from a larger contribution to the dust heating by the diffuse starlight. The location of the data relative to the models is also consistent with the fact that the spaxels in the central region of NGC 4449 typically require values of $U_{\text{min}} > 2$, as appropriate for an actively star-forming galaxy (see Section 7; Draine et al. 2007). The DL07 models do not produce values of $\log[F(8)/F(24)] > 0.5$, because of the assumed properties for the PAH population in the models, with single-photon heating contributing $F(24) \sim (1/3)F(8)$ (see Figure 15 of DL07). This assumption does not correspond to a physical limitation, and, indeed, several of the spaxels in NGC 4449 have $\log[F(8)/F(24)] > 0.5$.

When plotting the observed 8–24 µm luminosity ratio as a function of both the $H\alpha$ luminosity surface density and the $H\alpha$-to-4.5 µm luminosity ratio for all the 11" × 11" spaxels in the central region, irrespective of the S/N$_{1100}$ µm value, we find a general anticorrelation between the two quantities (Figure 6). The $H\alpha$ luminosity surface density is a proxy for the SFR surface density, while the ratio $H\alpha/4.5$ µm is a proxy for the spaxel-based specific SFR, since the 4.5 µm emission traces the old, diffuse stellar population; both quantities measure the local strength of the star formation.

The strength of the correlations is evaluated by applying the nonparametric Kendall $\tau$ test to the data. In this work, we prefer Kendall $\tau$ over Spearman $\rho$, as it is less sensitive to outliers (Feigelson & Jogesh Babu 2012; Ivezic et al. 2014); thus, Kendall $\tau$ works well also in the case of small data sets (e.g., Section 7). Figure 6 lists, in each panel, the probability $p_{\tau}$ that the data are uncorrelated: small values of $p_{\tau}$ correspond to high probability of correlation, with $p_{\tau} = 0.05$ indicating a 2$\sigma$ significance level for a two-tailed test. In our case, the $F(8)/F(24)$ luminosity ratio is highly anticorrelated with both $\Sigma(H\alpha)$ and the $F(H\alpha)/F(4.5)$ luminosity ratio, with a likelihood of $6.4 \times 10^{-18}$ and $4.4 \times 10^{-25}$, respectively, that the pairs of data are drawn from a random distribution. Thus, more strongly star-forming regions have lower $F(8)/F(24)$ ratios, in agreement with the results of Figures 4 and 5. As the plots show, this
<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Top: Infrared to millimeter color–color plots for the 11" × 11" spaxels (corresponding to deprojected ~223 pc × 595 pc) in the central region of NGC 4449. Red asterisks indicate the 131 spaxels with S/N$_{1100}$$_{µm}$ ≥ 2.0, and blue circles indicate the 42 spaxels with S/N$_{1100}$$_{µm}$ ≥ 3.5. Representative 1σ error bars are shown to the right in each panel. The data follow a well-defined trend in both color–color plots, through which a best-fit curve is drawn (black line). The expression for each curve is given in the text. The large black triangle with error bars is the mean value of the entire central region of Figure 2. Bottom: same as the top panels, but showing only the spaxels with S/N$_{1100}$$_{µm}$ ≥ 2.0 (red asterisks) with their 1σ error bar. The locus occupied by the Draine & Li (2007) models for the range U$_{min}$ = 0.1, ..., 25 (bottom to top), U$_{max}$ = 10$^{3}$ (light gray), 10$^{4}$ (dark gray), 10$^{5}$ (magenta), and 10$^{6}$ (cyan), and γ = 0.0, ..., 1.0 (right to left), at fixed q$_{PAH}$ = 3.2%, is shown by the colored tracks. The limit for log[νF$_{ν}$ (8)/νF$_{ν}$ (24)] ≤ 0.5 is built into the Draine & Li (2007) models.</figcaption></figure>

result holds for all spaxels in the central region, irrespective of their S/N$_{1100}$$_{µm}$. Regions of low 1100 µm flux and, in first approximation, low surface dust density are still well described by the behavior of the 8 µm/24 µm color ratio observed for higher-S/N$_{1100}$$_{µm}$ spaxels.

The smaller scatter of the F(70)/F(1100) versus F(8)/F(24) correlation relative to the F(160)/F(1100) versus F(8)/F(24) one in Figure 4 is also in line with the findings above, since F(160) receives a smaller contribution than F(70) from dust heated by recent star formation (e.g., Calzetti et al. 2010; Li et al. 2013). A minimum χ² fit through the highest

(S/N$_{1100}$$_{µm}$ ≥ 3.5) significance data of Figure 4 gives

$$
\log[F(70)/F(1100)] = 4.19 - \{(0.65 \pm 0.05)\log[F(8)/F(24)] + (0.75 \pm 0.04)\}^5 \quad (2)
$$

$$
\begin{equation}
\begin{aligned}
& \log[F(160)/F(1100)] = 3.67 \\
& \quad - \{(0.72 \pm 0.05)\log[F(8)/F(24)] + (0.62 \pm 0.04)\}^4.
\end{aligned}
\tag{3}
\end{equation}
$$
<figure><img src="image_7.png" /><figcaption>**Figure 5.** The top left panel of Figure 4 is now shown at the center of this figure. On the left and right sides are the data and the best-fit SED (top) and the identification of the location within the central region of the S/N<sub>1100 µm</sub> ≥ 3.5 spaxels with the most extreme color combinations in Figure 4 (bottom). The spaxel at the left of the infrared to millimeter color–color plot corresponds to a local peak of star formation, whose SED is well described by models that include 57% PDR contribution. The spaxel at the right of the infrared to millimeter color–color plot corresponds to a region of very low current star formation, and its SED is well described by a model that includes emission only from cold dust, i.e., *f*<sub>PDR</sub> = 0. The SFR surface density, Σ<sub>SFR</sub>, is indicated within the SED plots of each spaxel.</figcaption></figure>

When including also the spaxels with S/N<sub>1100 µm</sub> ≥ 2.0 in the 1100 µm emission, the fits still agree within the 1σ uncertainty. Thus, we consider these results sufficiently robust to be used to extrapolate F(1100), and therefore dust masses, also to spaxels with insufficient S/N (S/N<sub>1100 µm</sub> < 3.5) from the AzTEC observations.

## 5.2. Comparison with Other Galaxies

We compare the IR–millimeter relations of Equations (4) and (5) with the trends of an independent galaxy sample. We concentrate on the relation between log[F(8)/F(24)] and log[F(70)/F(1100)], which has the lowest dispersion of the two. Dale et al. (2017) collected and homogenized the UV, optical, infrared, submillimeter, and radio measurements for the galaxies in the KINGFISH sample (Kennicutt et al. 2011), augmented by the galaxies in the SINGS sample (Kennicutt et al. 2003) that are not already in KINGFISH. The KINGFISH sample and its parent SINGS sample are representative of star-forming galaxies within the local 30 Mpc (Kennicutt et al. 2003) and are therefore reasonable templates against which to compare our spatially resolved results. Measurements at 850 µm obtained with SCUBA are available for 27 of the 79

galaxies in Dale et al. (2017); an additional 14 galaxies have Planck photometry. Of these 41 galaxies, 38 were observed with both the SST and HSO, and the remaining 3 (MRK 33, NGC 5033, and NGC 7552) were observed with SST only.

For the galaxies with Planck photometry, we simply interpolate between the two closest frequency bands to the AzTEC frequency. For the galaxies with SCUBA data, we extrapolate the SCUBA 850 µm photometry to 1100 µm using the DL07 models appropriate for each galaxy,<fnref n="19" /> as described in Draine et al. (2007) and Dale et al. (2009). For the 70 µm emission of the galaxies, we use the HSO measurements when available; for the three galaxies with SST-only observations, we increase the MIPS70 values by 15% (0.06 in log scale), to account for the slight systematic discrepancy in the photometry of the two instruments (see photometric values listed in Dale et al. 2017).<fnref n="20" /> Photometry is corrected as discussed in Section 2.2.

<footnote n="19">¹⁹ Assumptions, found sometimes in the literature, that the dust emission follows a Rayleigh–Jeans tail approximation beyond 500 µm should be taken with caution. The deviation of the HSO SPIRE 500 µm emission from the R–J tail approximation ranges from a factor of 2.2 at 20 K to 68% at 30 K.</footnote>

<footnote n="20">²⁰ Our results do not change significantly whether this small shift in flux is included or not for all three galaxies: MRK 33, NGC 5033, and NGC 7552.</footnote>
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Ratio of 8 µm/24 µm luminosity as a function of both the Hα luminosity surface density (left) and the Hα/4.5 µm luminosity ratio (right), for all 11" × 11" spaxels in the central region (black), those with S/N<100 µm ≥ 2 (red), and those with S/N<1100 µm ≥ 3.5 (blue). The number of spaxels in each regime of S/N is listed in the figure. The anticorrelations observed in both panels suggest that regions of high SFR (traced by Hα) relative to both the subtended area and the luminosity of the old stellar populations (traced by the 4.5 µm emission) have proportionally higher 24 µm luminosity relative to the 8 µm one. The nonparametric Kendall probability pτ, shown in the panels, is the probability that the two quantities along the axes are uncorrelated.</figcaption></figure>

Figure 7 shows the results for the 41 KINGFISH+SINGS galaxies, in comparison with the S/N<1100 µm ≥ 2 spaxels from NGC 4449. The general locus occupied by the galaxies is consistent with that of the spaxels in NGC 4449, supporting the use of Equations (4) and (5) as tools to obtain dust mass measurements below the detection limits of the AzTEC observations.

The KINGFISH+SINGS galaxies, however, span a smaller range of F(8)/F(24) colors than the spaxels in NGC 4449, with the maximum value for the galaxies being log[F(8)/F(24)] ~ 0.53. This is consistent with the picture that the KINGFISH+SINGS galaxies are star-forming, while several spaxels in NGC 4449 are areas of little current star formation (and F(8)/F(24) > 0.53; see also Figure 5). The galaxy data also justify the assumption in the DL07 models to only include PAH populations that yield log[F(8)/F(24)] ≤ 0.5. In order to apply the DL07 models in a regime where they are valid, in the following sections we only use spaxels with log[F(8)/F(24)] ≤ 0.55.<fnref n="21" /> The right panel of Figure 7 shows the distribution of F(8)/F(24) values for all the spaxels in the central region of NGC 4449, irrespective of their S/N<1100 µm, and 91% of them (total of 393) are below log[F(8)/F(24)] = 0.55, suggesting that our upper limit to the F(8)/F(24) ratio is not expected to impact the generality of our results.

5.3. Extending the Millimeter Flux to Low-significance Regions

Figure 8 provides additional sanity checks for the model-
derived 1100 µm flux. The left panel compares the original
1100 µm measurements, irrespective of their S/N<1100 µm, with
the model-derived fluxes, from the application of Equation (2).
As expected, the spread between the two sets of data is large,

but the two sets show a general positive trend, i.e., larger
observed 1100 µm fluxes get assigned larger values also by the
model, as also highlighted by the one-to-one relation drawn on
the plot. A Kendall τ test returns a probability of 2 × 10⁻³⁰
that the two are drawn from a random distribution. This result,
which indicates a ~14.5σ correlation, lends additional cred-
ibility to the model of Equation (2).

There are highly discrepant points from the general relation of the left panel of Figure 8, which we investigate further. Highly deviant spaxels, i.e., more than a factor of 10 away from the one-to-one relation, have S/N<1100 µm < 1.8, and they are found at the edges of the starburst, in the NW and SE areas of the central region. Thus, the large deviations from our model flux predictions are attributable to the low significance of the measured fluxes. In all cases, the model of Equation (2) provides a far more physical 1100 µm flux than the (low-significance) measurements. Spaxels located in the lower right region of the left panel of Figure 8 have unphysically high 1100 µm observed fluxes, i.e., fluxes that are larger than any of the shorter-wavelength fluxes, while the spaxels located in the upper left region of the left panel of Figure 8 have too low 1100 µm observed fluxes, which would require unphysically low emissivity for the dust relative to that at shorter wavelengths. Based on all these considerations, we consider the model of Equation (2) a robust approach to derive 1100 µm flux in the low-significance spaxels of the central region of NGC 4449.

The right panel of Figure 8 compares the model 1100 µm fluxes with the 24 µm fluxes. This comparison is crucial for our analysis, since a linear correlation between the logarithm of the two quantities would undermine any conclusion relating the SFR to the gas density. The reason for this is because the derivation of the SFR includes use of the 24 µm flux and the surface density of gas is derived from the 1100 µm intensity (see below); however, F(24) and F(1100) are also linked by being measured from the dust emission SED and could

<footnote n="21">²¹ We include in our fits data with F(8)/F(24) ratios about 10% larger than the DL07 model limit of log[F(8)/F(24)] = 0.5, to make allowance for measurement uncertainties.</footnote>
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Left: same as the top left panel of Figure 4, but now shown in gray scale for the spaxels of NGC 4449. The black triangle with error bars is the mean value of the central region of NGC 4449, as shown in Figure 4, as well. The infrared to millimeter color–color integrated data of 41 nearby galaxies from Dale et al. (2017) are shown as blue and red symbols, together with $1\sigma$ error bars. The 1100 µm data of the galaxies are obtained from either the interpolation of the *Planck* data (red) or the extrapolation of the SCUBA 850 µm measurements (blue). We use the distribution of local galaxies to set an upper limit of $log[F(8)/F(24)] = 0.55$, indicated by the vertical cyan line, in order to include only spaxels to which the DL07 models can be applied. Right: histogram of the $F(8)/F(24)$ color distribution of the spaxels in the central region of NGC 4449, irrespective of $S/N_{1100\mu m}$. Most (91%) of the spaxels are below the value $log[F(8)/F(24)] = 0.55$ (vertical cyan line), i.e., have values consistent with those of the nearby galaxies.</figcaption></figure>

<figure><img src="image_10.png" /></figure>

therefore be covariant. The right panel of Figure 8 serves the purpose of verifying that the two quantities are not covariant in a simple way, i.e., through a linear correlation, although they are linked via a complex (polynomial) relation. From this, we conclude that our results, based on the use of $F(1100)$ for deriving dust (and gas) masses and of $F(24)$ as the nondominant contributor to the SFR, are not trivially covariant. This is further helped by the fact that for most spaxels the SFR is mainly traced by $H\alpha$.

In Section 7, we derive dust mass surface densities for all spaxels with $log[F(8)/F(24)] \le 0.55$, using directly the measured AzTEC 1100 µm flux for the 36 $11'' \times 11''$ spaxels with $S/N \ge 3.5$ and Equation (2) for the remaining 357 spaxels. Maps of the model 1100 µm flux for the central region at $11''$ and at $44''$ resolution are shown in Figure 9. The sum of the values in all the spaxels gives a flux density of 0.28 Jy, about 10% higher than the value directly measured from the original AzTEC map, but well within the uncertainties of the central region 1100 µm flux. If we extend the model of Equation (2) to all spaxels in the galaxy, we obtain a total 1100 µm flux of ~0.47 Jy. This is a little over 20% higher than the flux directly measured from the AzTEC image (Table 2), although consistent within the uncertainties. We caution the reader that this result is obtained by applying Equation (2) to spaxels that are outside the central region, with $S/N < 10$ at 8, 24, and 70 µm, i.e., outside the regime where Equation (2) has been derived. The trend for larger regions to yield larger discrepancies between the measured and modeled 1100 µm flux supports our earlier statement that some of the extended flux in the original image is lost owing to spatial filtering. The flux of 0.47 Jy is about 10% lower than the value of 0.53 Jy from the interpolation of the *Planck* photometry. As already stated in Section 2.3, this small discrepancy could be due to

contamination from foreground and background sources in the *Planck* data (Dale et al. 2017).

## 6. Global Dust Masses

We derive both the total and central region (Figures 1 and 2) dust masses for NGC 4449, in order to constrain some of the free parameters in the DL07 models, to be used in the spaxel-based analysis. For the large-scale fits, we increase the number of constraints in the fits by using three additional photometric data points (the three HSO/SPIRE bands) that cannot be used for the spatially resolved analysis owing to resolution limitations.

### 6.1. The Whole Galaxy

For the whole-galaxy photometry, the best-fit values of $U_{min}$, $U_{max}$, and $f_{PDR}$ (Figure 10, left panel, and Table 3) are consistent with those derived by Draine et al. (2007) for the starburst galaxies in the SINGS sample (Kennicutt et al. 2003). As with these authors, our results are not very sensitive to the specific value of $U_{max}$, although we continue to carry this as a free parameter. We do not include the *Planck* photometry in the fits, for uniformity of treatment of both the whole galaxy and central region. The left panel of Figure 10 shows the data points with the best fit from the models. Although not used for the fit, the *Planck* photometry is in general agreement with the model expectation: there is a small positive offset, by about 10%, for the data longward of 350 µm, which has been discussed in the previous section. The best-fit parameters for the entire galaxy remain virtually unchanged whether we include our millimeter data point or not, as the data from SPIRE provide sufficient constraints to the long-wavelength behavior of the SED. The total dust mass we derive (Table 3) is consistent with the value
<figure><img src="image_11.png" /><figcaption>**Figure 8.** Left: measured flux at 1100 µm compared with the flux, $F_{1100,mod}$, obtained from the model of Equation (2). Data with $S/N_{1100}$ µm ≥ 3.5 retain the original measurements in both axes (red circles), while the model is applied to all other data. The measurements are reported for all $S/N_{1100}$ µm values (black circles), including negative/null values, which have been assigned an arbitrarily low flux for plotting purposes (blue circles). There is a general positive trend between the two sets of measurements, as highlighted by the one-to-one solid black line. The dashed black lines mark the locus of factor of 10 deviation from the one-to-one relation. The relation between the two quantities is highly significant; the Kendall τ probability of a random distribution is indicated by the parameter $p_τ$ in the panel and corresponds to a 14.5σ deviation from the null hypothesis. The Kendall τ test is applied only to the data with $S/N_{1100}$ µm < 3.5 (i.e., excluding the red circles). This lends further support to the relation linking mid-IR to far-IR/millimeter fluxes. Highly deviant points from the relation are discussed in the text. Right: 1100 µm measurements from the model of Equation (2) as a function of the flux at 24 µm. Typical 1σ error bars are shown for both F(24) (top of the panel) and F(1100) (right of the panel); for the latter we include the dispersion in the relation of Equation (2) as inferred from Figure 4. Red circles show the location of the 1100 µm measurements with $S/N > 3.5$, to which no model is applied.</figcaption></figure>

of $3.5 \times 10^6 M_⊙$ obtained by Karczewski et al. (2013), once their value is rescaled to our preferred distance of 4.2 Mpc.

If we perform the fits using the SST/MIPS70, rather than the HSO/PACS70, i.e., include a 30% lower 70 µm photometric point, we obtain results that are consistent with those that use the PACS70 photometry, within the uncertainties. As an example, the best-fit dust mass now is $M_{dust} \sim 3.8 \times 10^6 M_⊙$, about 12% higher than the mass obtained when using the PACS70 point, but consistent with it within the 1σ error bar.

Recent results from Planck indicate that the dust masses derived using the Draine & Li (2007) models need to be reduced by some factor, in order to reconcile the implied dust opacities with the measurements for the MW (Fanciullo et al. 2015; Planck Collaboration 2016b). The reduction factor is a function of $U_{min}$, decreasing for increasing $U_{min}$ and with an average value of ~2. The analysis in both papers only includes values of $U_{min} \lesssim 1$, which are well below what we measure in NGC 4449. A simple linear extrapolation of the observed trends for $U_{min} > 1$ suggests that the reduction factor for our case is $\lesssim 1.5$. We carry forward our analysis into the next section adopting two cases: a no-reduction case and a case in which dust masses need to be reduced by a factor of ~1.5 (Table 3).

## 6.2. The Central Region

The best-fit model to the IR SED of the central region gives parameters that are similar to those of the whole galaxy, with the exception of the dust mass, which is about 2/3 in value (Figure 10, right panel). A result of the fit is also that at 24 µm about 1/3 of the emission is contributed by dust heated by the diffuse starlight, which is not directly related to the current star

formation event and needs to be removed when using the 24 µm emission as an SFR estimator.

In order to further verify that our dust-fitting parameters are reasonable, we derive a CO-to-H$_2$ conversion factor from the dust mass of the central region and compare it with other derivations in the literature. The dust mass is related to the total (atomic+molecular) hydrogen mass via the relation

$$M_{\text{dust}} = (\text{D/H})M_{\text{H}}, \qquad (4)$$

where $M_H = [M(H\,I) + M(H_2)]$ is the hydrogen mass and D/H is the dust-to-hydrogen ratio, i.e., it is the dust-to-gas mass ratio uncorrected for the helium factor. We will refer all our calculations in this section to the hydrogen mass, rather than the gas mass. Taking into account the subsolar metallicity of NGC 4449, we expect

$$\text{D/H} \sim 0.01 \frac{(\text{O/H})}{(\text{O/H})_{\text{MW}}} \sim 0.003 \qquad (5)$$

(Draine et al. 2007), where $(\text{O/H})/(\text{O/H})_{\text{MW}}$ is the ratio of the oxygen abundance in NGC 4449 to the MW oxygen abundance, which we adopt to be 0.3 on average (Table 1). Our assumed value for D/H is consistent with the preferred fit to the D/G versus metallicity data of Rémy-Ruyer et al. (2014), once we account for the factor of 1.36 for helium and metals included in the models by those authors. The inferred hydrogen mass is listed in Table 3; the uncertainties include the 20% variation within the central region due to the gradient in oxygen abundance. The H\,I mass (uncorrected for helium), as calculated from the THINGS map for the central region, represents about 40% (or 60% if we consider the Planck
<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /><figcaption>**Figure 9.** Maps at 1100 µm of the central region of NGC 4449 with 11" × 11" resolution (left) and 44" × 44" resolution (right). The 1100 µm fluxes are from the original measurements for spaxels with S/N$_{1100}$$_{µm}$ ≥ 3.5 and from the model of Equation (2) for the other spaxels; only spaxels with log[F(8)/F(24)] ≤ 0.55 are assigned valid F(1100) fluxes, for a total of 393 spaxels in the maps to the left. The units of the maps are shown by the scale to the right of each panel.</figcaption></figure>

correction factor for the dust mass) of the hydrogen mass. Subtracting the H I mass from the total hydrogen mass yields the molecular hydrogen masses listed in Table 3. The ratio M(H$_2$)/M(H I) ~ 1.5 and 0.7, for the uncorrected and Planck-corrected dust mass cases, respectively, is consistent with the range of molecular-to-atomic hydrogen ratio found in nearby galaxies (Leroy et al. 2005; Saintonge et al. 2011).

The relation between H$_2$ mass and CO luminosity is

$$ \frac{M(\text{H}_2)}{M_{\odot}} = 3.2 \left( \frac{X_{\text{CO}}}{X_{\text{CO,MW}}} \right) \frac{L(\text{CO})}{\text{K km s}^{-1} \text{ pc}^2} \quad (6) $$

(e.g., Bolatto et al. 2013), where X$_{CO}$ is the “X”-factor, with the MW value being X$_{CO,MW}$ = 2 × 10$^{20}$ cm$^{-2}$ (K km s$^{-1}$)$^{$_{−}$1}$; the factor 3.2 reflects our choice not to include helium corrections. The CO (1–0) luminosity of the central region is L(CO) = 1.02 × 10$^{7}$ K km s$^{-1}$ pc$^{2}$ (Böttner et al. 2003), which yields the X-factors listed in Table 3. Within the uncertainties, both values are consistent with the value (X$_{CO}$/X$_{CO,MW}$) ≈ 11 determined$^{22}$ by Böttner et al. (2003) and with the range 5–25 expected for galaxies with metallicity ~1/3 solar (Bolatto et al. 2013). However, the model by Accurso et al. (2017) predicts (X$_{CO}$/X$_{CO,MW}$) ~ 7 for a galaxy with the oxygen abundance and the main-sequence offset of NGC 4449; this value is marginally more consistent with the lower X$_{CO}$/X$_{CO,MW}$ derived with a 1.5 reduction factor on the dust mass (last column of Table 3). One additional consideration is that, based on the above, the reduction factor is unlikely to be as large as the factor of ~2 appropriate for the diffuse ISM (Planck Collaboration 2016b), which would yield a molecular gas mass consistent with zero within the uncertainties and inconsistent with the presence of strong star formation within this galaxy.

Outside of the central region, the total hydrogen mass is between (3.6 ± 1.5) × 10$^{8}$ M$_{⊙}$ and (2.5 ± 1.1) × 10$^{8}$ M$_{⊙}$, the latter after inclusion of the Planck correction factor. In the same outer region, the H I mass is (1.9 ± 0.1) × 10$^{8}$ M$_{⊙}$. Since the central region contains almost 90% of the star formation in the galaxy, we expect the external regions to contain little molecular gas, and the total hydrogen mass should be the same as the H I mass. In this case, the total hydrogen mass reduced by the 1.5 Planck factor is marginally more consistent with the H I mass than the unreduced case, leading to a preference for using the former value over the latter.

The results above, from the (X$_{CO}$/X$_{CO,MW}$) comparison within the central region and the M$_{dust}$/M(H I) comparison outside the central region, show better consistency when the dust masses derived from the Draine & Li (2007) models are decreased by a factor of 1.5. We adopt this correction in the rest of this work.

## 7. Spatially Resolved Dust and Gas Masses

Spaxel-based dust masses are derived by fitting the luminosity densities at 8, 24, 70, 100, 160, and 1100 µm in each spaxel. Since we are only fitting six photometric data points for each spaxel, we are only allowed to have four free parameters in order not to overfit the data. Our results from Section 6 enable us to constrain several model parameters, thus simplifying the fits. We fix the extinction curve/q$_{PAH}$ combination to be MW/3.2%, as in the whole galaxy and central region best fits (Section 6). As already mentioned above and in other papers (Draine et al. 2007), the derived dust masses are not sensitive to the choice of U$_{max}$, which we fix at a value of 10$^{5}$. We find that most of the regions with S/N$_{1100}$$_{µm}$ ≥ 3.5 are fit by U$_{min}$ in the range 2–8; we thus adopt U$_{min}$ = 5, as derived for the central region fit (Section 6.2), for all spaxels in the region. Variations in the derived dust masses for different choices of U$_{min}$ in the range 2–8 are about 30%, much less than our typical uncertainty in the dust mass value. In summary, we can fix four of the six free parameters in the models (extinction curve, q$_{PAH}$, U$_{min}$,

<footnote>$^{22}$ Böttner et al. (2003) do not provide an uncertainty for their X$_{CO}$ value, although they state that the dust mass from which it is derived has a large uncertainty, within a factor of 3.</footnote>
<figure><img src="image_14.png" /><figcaption>**Figure 10.** Dust SED fits of the whole galaxy (left) and of the central region (right; ~4.7 × 4.7 kpc² deprojected to ~4.7 × 12.5 kpc²). The black lines are the best fits through the data, the red triangles are the resulting broadband photometry from the best-fit models, and the blue circles are the photometric data used in the fits with the 1σ error bars indicated. For the whole galaxy we also report the Planck photometric measurements with their uncertainties as magenta circles; these data are not used in the fits, but they are used here to confirm the robustness of the model fits. The best-fit SEDs are constructed from a dust mixture heated by two starlight intensity components: the general interstellar medium (black dotted line), described by the energy density parameter $U_{min}$, and regions with a power-law distribution of intensities, $dM_{dust}/dU \propto U^{-2}$, between $U_{min}$ and $U_{max}$ (black dashed line). Recent results from the Planck Collaboration indicate that the dust masses from the Draine & Li (2007) models, and listed in the two panels, may need to be reduced by a factor of as much as ~1.5 (see Section 6 for details).</figcaption></figure>

<table><caption>**Table 3**

Derived Quantities for NGC 4449</caption><thead><tr><th>Parameter (Units)</th><th>Galaxy<sup>a</sup></th><th>Central Region<sup>a</sup></th></tr></thead><tbody><tr><td><i>U</i><sub>min</sub><sup>b</sup></td><td>4 ± 1</td><td>5 ± 1</td></tr><tr><td><i>U</i><sub>max</sub><sup>b</sup></td><td>~10<sup>5</sup></td><td>~10<sup>5</sup></td></tr><tr><td><i>f</i><sub>PDR</sub><sup>b</sup></td><td>0.12 ± 0.04</td><td>0.15 ± 0.04</td></tr><tr><td><i>M</i><sub>H 1</sub> (<i>M</i><sub>☉</sub>)<sup>c</sup></td><td>(5.0 ± 0.3) × 10<sup>8</sup></td><td>(3.1 ± 0.1) × 10<sup>8</sup></td></tr><tr><td></td><td>D&L07</td><td>D&L07+Planck</td><td>D&L07</td><td>D&L07+Planck</td></tr><tr><td><i>M</i><sub>dust</sub> (<i>M</i><sub>☉</sub>)<sup>d</sup></td><td>(3.4 ± 0.6) × 10<sup>6</sup></td><td>(2.3 ± 0.4) × 10<sup>6</sup></td><td>(2.3 ± 0.4) × 10<sup>6</sup></td><td>(1.5 ± 0.3) × 10<sup>6</sup></td></tr><tr><td><i>M</i><sub>H</sub> (<i>M</i><sub>☉</sub>)<sup>e</sup></td><td>(11.3 ± 2.9) × 10<sup>8</sup></td><td>(7.6 ± 1.9) × 10<sup>8</sup></td><td>(7.7 ± 2.0) × 10<sup>8</sup></td><td>(5.1 ± 1.3) × 10<sup>8</sup></td></tr><tr><td><i>M</i><sub>H2</sub> (<i>M</i><sub>☉</sub>)<sup>f</sup></td><td>(6.3 ± 3.0) × 10<sup>8</sup></td><td>(2.6 ± 2.0) × 10<sup>8</sup></td><td>(4.6 ± 2.0) × 10<sup>8</sup></td><td>(2.0 ± 1.3) × 10<sup>8</sup></td></tr><tr><td><i>X</i><sub>CO</sub>/<i>X</i><sub>CO,MW</sub><sup>g</sup></td><td>...</td><td>...</td><td>14.4 ± 6.1</td><td>6.1 ± 3.9</td></tr></tbody></table>

**Notes.**

<sup>a</sup> Derived quantities, from best-fit dust SED models (Draine & Li 2007), for the entire galaxy (second column) and for the central region as defined in the text and in Table 2 (third column).

<sup>b</sup> The minimum and maximum ($U_{min}$ and $U_{max}$) energy density parameters and the fraction ($f_{PDR}$) of dust luminosity due to current star formation and other activity, as defined in Draine & Li (2007). The SED fits are not strongly sensitive to the value of $U_{max}$.

<sup>c</sup> H I mass, in solar masses and uncorrected for helium, derived from the THINGS maps of Walter et al. (2008), using the relation $M(H I) = 2.356 \times 10^5 D^2 \Sigma_i S_i \delta v$, where D is the galaxy's distance in Mpc and $\Sigma_i S_i \delta v$ is the integral of the intensity along the line-of-sight velocity, in units of Jy km s⁻¹. The uncertainty in the H I mass estimate is ~6%, after including registration errors. For the “galaxy” measurement, the H I mass is determined in the same area used for the dust mass, and the resulting mass is about 1/4 of the total atomic hydrogen mass associated with the galaxy (Huchtmeier & Richter 1989).

<sup>d</sup> The mass in dust for the entire galaxy and for the central region, in solar masses, as derived from the models of Draine & Li (2007, hereafter D&L07) and from the same models reduced by a factor of 1.5, to account for recent findings by the Planck Collaboration (D&L07+Planck; Sections 6.1 and 6.2).

<sup>e</sup> Total gas mass, uncorrected for helium, from the dust mass using Equations (4) and (5).

<sup>f</sup> The molecular hydrogen mass, uncorrected for helium, derived from subtracting $M_{H 1}$ from $M_H$.

<sup>g</sup> The ratio of the NGC 4449 CO-to-H₂ conversion factor, or X-factor $X_{CO}$, to the MW's $X_{CO,MW}$, derived from Equation (6), using the H₂ molecular hydrogen mass in this table together with the CO luminosity measured by Böttner et al. (2003). We adopt an MW value $X_{CO,MW} = 2 \times 10^{20} \text{ cm}^{-2} (\text{K km s}^{-1})^{-1}$.

and $U_{max}$) without losing generality; the two remaining parameters, $f_{PDR}$ and $M_{dust}$, are left as free parameters in the fits. Figure 11 shows, for the spaxels with $S/N_{100} \mu m \ge 3.5$, the range of $f_{PDR}$ we fit as a function of the SFR surface density

(Section 3). Dust mass surface densities, $\Sigma_{dust}$, in units of $M_{\odot} \text{ pc}^{-2}$, are derived by correcting each area by the inclination of 68° (Table 1). We further divide the dust mass surface densities by a factor of 1.5, as discussed in Section 6. Figure 12 shows
<figure><img src="image_15.png" /><figcaption>Figure 11. Range of $f_{PDR}$ values, the fraction of the dust luminosity radiated from regions with $U > 100$, as a function of the SFR surface density, $\Sigma_{SFR}$, for the $36\ 11'' \times 11''$ spaxels in the central region with $S/N_{1100\mu m} \ge 3.5$ and $\log[F(8)/F(24)] \le 0.55$. The positive correlation between the two quantities is expected, since stronger star-forming regions will contribute a higher fraction of the dust luminosity.</figcaption></figure>

the maps of the fitted values of $f_{PDR}$ and $\Sigma_{dust}$ at $11''$ and $44''$ resolution, together with the maps of $\Sigma_{SFR}$. The sum of the dust masses derived individually for each spaxel is within 12% of the dust mass derived from the global fit of the central region. This discrepancy is smaller than the typical uncertainties on the derived dust masses, including the global ones, supporting the robustness of our approach (see discussion in Galametz et al. 2012).

The use of the $1100\ \mu m$ data point, or other submillimeter/millimeter dust-dominated emission, provides important constraints to the dust masses, as already found and discussed in Draine et al. (2007). Figure 13 shows the effect of including or excluding the $1100\ \mu m$ data in the SED fits of the spaxels with $S/N_{1100\mu m} \ge 3.5$. For this experiment, we have allowed $U_{min}$ to vary in the fitting. The fits yield higher dust surface densities when the $1100\ \mu m$ data are excluded from the fits, with a median offset of about 60%. The reason is that, in the absence of a submillimeter or millimeter data point, the best fits tend to yield lower $U_{min}$ values, and therefore larger dust masses, relative to the cases that include the long-wavelength data (Draine et al. 2007). If we exclude the data at $1100\ \mu m$, the best-fit $U_{min}$ is in the range 0.5–8 for the $S/N_{1100\mu m} \ge 3.5$ spaxels, as opposed to 2–8, as reported above. Figure 13 also shows a mild trend for larger offsets in dust surface densities for regions with larger $\Sigma_{SFR}$, which can directly impact any conclusion that relates SFRs to gas (dust) masses. Specifically, for $\Sigma_{SFR} > 0.015\ M_{\odot}\ yr^{-1}\ kpc^{-2}$, the dust surface densities are overestimated by about a factor of 2, in the absence of the $1100\ \mu m$ data point.

Gas mass surface densities, $\Sigma_{H1+H2}$, are derived from the dust mass surface densities using Equation (4), after implementation of a metallicity-dependent gas-to-dust ratio (Equation (5)), using the metallicity gradient derived by

Pilyugin et al. (2015). Variations in the D/H ratio can be as large as a factor of ~3 within galaxies (Roman-Duval et al. 2014, 2017), which is larger than the systematic change from the metallicity gradient ($\Delta(O/H)/(O/H) < 0.15$) across the central region; however, they are smaller than the dynamic range of our gas mass data. We multiply all gas mass densities by 1.36 to include contributions from helium and metals, for easier comparison with results from the literature. The gas surface densities are then converted to $H_2$ surface densities, $\Sigma_{H2}$, by subtracting the H I (also corrected for He and metals) in each spaxel. We carry out our analysis on all positive values of $\Sigma(H_2)$ that are larger than the combined uncertainty of total gas measurement and H I detection limit; we set this as our lower limit to $\Sigma(H_2)$, in order to avoid spurious values for the molecular gas surface density. Since $\Sigma(H_2)$ results from subtracting $\Sigma_{H1}$ from $\Sigma_{H1+H2}$, spurious or unphysical values can arise when the gas mass surface density is close in value, within the uncertainties, to the H I surface density. There are about 230 spaxels above this threshold (out of 393), or about 60% of the available ones. This fraction, ~60%, of usable spaxels remains the same also at larger binning scales. Figure 14 shows the $H_2/H1$ mass ratio as a function of the $H_2$ surface density; higher surface densities correspond to higher $H_2/H1$ ratios, in agreement with the small dynamic range in H I mass surface densities observed within galaxies (e.g., Kennicutt et al. 2007; Bigiel et al. 2008). The range of mass ratios is similar for both the smaller ($11'' \times 11''$) and larger ($44'' \times 44''$) regions and is in agreement with what is observed in most local galaxies (Saintonge et al. 2011).

# 8. Results and Discussion

## 8.1. The Relation between Star Formation and Gas

Following previous authors (e.g., Kennicutt 1998; Kennicutt et al. 2007; Bigiel et al. 2008; Rahman et al. 2012; Kennicutt & Evans 2012; Leroy et al. 2013), we produce scatter plots of the surface density of the SFR as a function of both the surface density of molecular gas (H2) and the total atomic+molecular gas (H1+H2). We fit linear relations to the scatter plots in logarithmic scale, using the ordinary least-squares (OLS) bisector linear fitting algorithm (Isobe et al. 1990). This algorithm is commonly employed in studies of the scaling relations of star formation, and we adopt it to facilitate comparisons with previous investigations.<fnref n="23" /> The best-fit lines through the data are expressed as

$$ \log(\Sigma_{SFR}) = \gamma_{H2} \times \log(\Sigma_{H2}) + A_{H2} \quad (7) $$

and

$$ \log(\Sigma_{SFR}) = \gamma_{H1+H2} \times \log(\Sigma_{H1+H2}) + A_{H1+H2}, \quad (8) $$

with $\gamma_{H2}$ and $\gamma_{H1+H2}$ the best-fit slopes and $A_{H2}$ and $A_{H1+H2}$ the intercepts. The widths of the distribution of the data about the best-fitting lines are also measured, in terms of standard deviations $\sigma_{H2}$ and $\sigma_{H1+H2}$.

The scatter plots, together with the best-fit lines, are shown in Figure 15 for the smallest and the largest spaxel sizes,

<footnote n="23">²³ The value of the linear fit slope $\gamma$ depends on the algorithm used. Calzetti et al. (2012) show from simulations that the bilinear regression fitting algorithm FITEXY (from the numerical recipes) yields systematically larger values of $\gamma$ than the OLS bisector fitting algorithm, by $\Delta\gamma \sim 0.2-0.5$ for an intrinsic linear slope $\gamma = 1.5$, when reasonable uncertainties are included in the data. The differences are larger or smaller for higher and lower slope, respectively.</footnote>
<figure><img src="image_16.png" /><figcaption>**Figure 12.** Maps of $f_{\text{PDR}}$, the fraction of dust luminosity emitted by star-forming regions (top panels), $\Sigma_{\text{dust}}$, the dust surface density (middle panels, in units of $M_{\odot} \text{pc}^{-2}$), and $\Sigma_{\text{SFR}}$, the SFR surface density (bottom panels, in units of $M_{\odot} \text{yr}^{-1} \text{kpc}^{-2}$), at both 11" × 11" (left column) and 44" × 44" (right column) resolution. $f_{\text{PDR}}$ and $\Sigma_{\text{dust}}$ are best-fitting parameters from the dust SED of each spaxel (Section 7). $\Sigma_{\text{SFR}}$ is from Equation (1).</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>**Figure 13.** Ratio of the dust surface densities in the S/N<sub>1100 µm</sub> ≥ 3.5 spaxels of the central region in NGC 4449, as determined from best fits of the IR SEDs including (Σ<sub>dust</sub>) or excluding (Σ<sub>dust,no 1100</sub>) the data points at 1100 µm. When excluding the 1100 µm data, the dust surface densities are systematically higher, with a median offset of about 60% and a slight trend for higher offsets in more strongly star-forming regions, i.e., those with higher values of Σ<sub>SFR</sub>.</figcaption></figure>

365 pc (11") and 1460 pc (44"), and listed in Table 4 for all four spaxel sizes considered here. The best-fit parameters for the Σ<sub>SFR</sub>−Σ<sub>H2</sub> relation are derived in two conditions: (1) the fits are limited to the H2 surface densities that are above the threshold discussed in the previous section, i.e., those to the right of the hatched yellow (cyan) region in Figure 14 (these fits are listed in the first line for each region size in Table 4); and (2) the fits include data down to 1/10 of the threshold (second line of Table 4 for each region's size, in parentheses). The second set of fits is performed to evaluate the stability of the fits with the more stringent limit.

The Kendall τ test is applied to each set of Σ<sub>SFR</sub> and Σ<sub>H2</sub>
(or Σ<sub>H1+H2</sub>) pairs, to evaluate the strength of each correlation,
and the resulting p<sub>τ</sub> values are listed in Table 4. For all spaxel
sizes, the correlation between Σ<sub>SFR</sub> and Σ<sub>H2</sub> is generally weak,
between 2σ and 3σ. Conversely, the correlations between Σ<sub>SFR</sub>
and Σ<sub>H1+H2</sub> are significant, with values between 4σ and 12σ,
for region sizes between 1.1 kpc and 360 pc, with only the
largest regions at 1.5 kpc showing a weak (2.5σ) correlation.
The largest spaxels suffer from small number statistics, which
account for the low correlation significance.

The smallest region sizes we analyze, ~360 pc on a side, are
sufficiently small to be affected by stochastic sampling of the
stellar IMF, in the sense that the IMF may not be fully
populated, especially at the lowest Σ<sub>SFR</sub> values (see Calzetti
2013, and references therein). We evaluate the scatter
introduced by stochastic IMF sampling by using the models
for the ionizing photon rate calculated by Cerviño et al. (2002).
We convert these authors’ values to our preferred Kroupa
(2001) IMF and recall that only 60% of the SFR is emerging
directly in Hα. The remaining ~40% is captured by dust and
measured through the 24 µm emission. This band receives
most of its heating contribution from the nonionizing UV
radiation, which is about 3–4 times less sensitive to stochastic

<figure><img src="image_18.png" /><figcaption>**Figure 14.** H<sub>2</sub>/H I surface density ratio as a function of the H<sub>2</sub> surface density, Σ(H<sub>2</sub>), for the 11" × 11" (black circles) and 44" × 44" (blue triangles) spaxels. Red circles are 11" × 11" spaxels with S/N<sub>1100 µm</sub> ≥ 3.5. The yellow region corresponds to the threshold defined in Section 7 at 11" × 11" resolution, while the cyan region is the same value at 44" × 44". Typical error bars for the faint Σ(H<sub>2</sub>) are also shown. There are no significant differences in the trend and dynamic range of the H<sub>2</sub>/H I mass ratio, for spaxels at the two resolutions. We request that all our H<sub>2</sub> values be larger than the 1σ uncertainty (i.e., be to the right of the color regions) to be considered valid measurements. All 44" × 44" spaxels are above the H<sub>2</sub> threshold (cyan region), and 95% of the 11" × 11" spaxels are (yellow region).</figcaption></figure>

sampling than the ionizing radiation (e.g., Andrews
et al. 2013). The two contributions of stochastic sampling to
the ionizing and nonionizing photons are shown in Figure 15
(top panels) as vertical bars: the effect of stochastic IMF
sampling decreases for increasing Σ<sub>SFR</sub>, and the scatter it
produces is generally smaller than the observed scatter in the
data. Thus, stochastic sampling of the IMF has a small effect on
our results, even at the smallest region sizes.

Our best-fit lines are compared with similar results from the
literature; where possible, we attempt to match the region sizes
used by the authors we compare our results against. Kennicutt
et al. (2007) measured the SFR–gas relation of ~500 pc size
star-forming regions in NGC 5194, for both the molecular gas,
obtaining γ<sub>H2</sub> = 1.37, and the total gas, obtaining
γ<sub>H1+H2</sub> = 1.56. We compare these results with those of our
365 pc spaxels (top panels of Figure 12). Bigiel et al. (2008)
and Rahman et al. (2012) measured the SFR–molecular gas
relation of samples of nearby galaxies with a spatial scale of
~0.75–1 kpc, obtaining γ<sub>H2</sub> = 0.95 and γ<sub>H2</sub> = 1.1, respec-
tively. Similarly, Leroy et al. (2013) obtain a slope of 1 ± 0.15,
although Shetty et al. (2014b), using a Bayesian approach to
the fits, finds a sublinear slope, γ<sub>H2</sub> ~ 0.7–0.8. We compare the
Bigiel et al. and Rahman et al. results with our 1460 pc sized
spaxels, although an equally good agreement would be
obtained if we compared those authors’ results with our
1 kpc spaxels. For the comparison with γ<sub>H1+H2</sub> at our largest
spaxel size, we use the mean value derived for whole galaxies
by Kennicutt (1998), γ<sub>H1+H2</sub> = 1.4. In contrast with this, we
should note that Bigiel et al. (2008) find γ<sub>H1+H2</sub> to vary from
<figure><img src="image_19.png" /></figure>

**Figure 15.** Scaling relations between the SFR and gas surface densities for two of the four spaxel sizes considered in this work: ~360 pc (top panels) and ~1.5 kpc (bottom panels). Both the relations for the molecular gas (left panels) and the total (atomic+molecular) gas (right panels) are shown. The data for each spaxel are shown as black circles; the red circles in the top panels are the spaxels with $S/N_{100} \mu m \ge 3.5$. Representative $1\sigma$ uncertainties for the data are reported in each panel. The yellow-shaded region is the lower limit for the $\Sigma_{H2}$ data points (see Figure 14 and text). The black vertical bars to the right of the top panels show the magnitude of the uncertainty due to stochastic sampling of the stellar IMF in the spaxels, which decreases for large SFR surface densities (Cerviño et al. 2002). Stochastic sampling is negligible for large galactic regions (bottom panels). The black lines show the OLS bisector linear best fits through the data, which are compared with best fits for other galaxies in the literature: Kennicutt et al. (2007, K+2007), Rahman et al. (2012, R+2012), Bigiel et al. (2008, B+2008), and Kennicutt (1998, K+1998). The slopes of the best linear fits, $\gamma$, and the Kendall $\tau$ correlation probabilities, $p_{\tau}$, are listed in each panel. See the text for details.

galaxy to galaxy in the range 1.1–2.7, when analyzing regions $\lesssim 1$ kpc in size. We briefly discuss these results in Section 8.3.

## 8.2. Simulations

We compare our results with the simulations of Calzetti et al. (2012), for the relation $\Sigma_{SFR}-\Sigma_{H2}$. Many details are provided in that paper, and here we summarize only the elements relevant to the current analysis. In the simulations, realizations of galaxy regions with sizes in the range 0.2–5 kpc are randomly populated with molecular clouds distributed according to an exponentially decreasing filling factor (to favor low-filling

factor regions, <15%, as observed). The clouds are extracted from a mass distribution with slope $\alpha$, in the range 500–$M_{\text{max}}$ $M_{\odot}$. We adopt here the same $\alpha = 2$ as those authors, but we lower their value of $M_{\text{max}}$ from $3 \times 10^7 M_{\odot}$ (appropriate for large late-type spirals like NGC 5194) down to $\sim 10^6 M_{\odot}$. This value of $M_{\text{max}}$ is what is observed in the LMC (Hughes et al. 2010), and we consider it more appropriate for low-mass dwarfs such as NGC 4449. The simulated molecular clouds, furthermore, obey the Larson laws (Larson 1981), which allow us to link their masses to sizes, and generate filling factors within spaxels; scatter in the parameters is added to the simulated data, as derived from observations.
<table><caption>**Table 4**
Star Formation Rate and Gas Surface Density Correlations</caption><thead><tr><th rowspan="2">Region Size<sup>a</sup><br>(pc)</th><th colspan="5">Σ<sub>SFR</sub> versus Σ<sub>H2</sub><sup>b</sup></th><th colspan="5">Σ<sub>SFR</sub> versus Σ<sub>H1+H2</sub><sup>b</sup></th></tr><tr><th>N<sup>c</sup></th><th>γ<sub>H2</sub><sup>c</sup></th><th>A<sub>H2</sub><sup>c</sup></th><th>σ<sub>H2</sub><sup>c</sup></th><th>p<sub>τ</sub><sup>c</sup></th><th>N<sup>d</sup></th><th>γ<sub>H1+H2</sub><sup>d</sup></th><th>A<sub>H1+H2</sub><sup>d</sup></th><th>σ<sub>H1+H2</sub><sup>d</sup></th><th>p<sub>τ</sub><sup>d</sup></th></tr></thead><tbody><tr><td>360</td><td>229<br>(247)</td><td>1.49 ± 0.13<br>(1.34 ± 0.09)</td><td>-3.81 ± 0.12<br>(-3.61 ± 0.08)</td><td>0.41 ± 0.05</td><td>2.3E-03</td><td>390</td><td>1.51 ± 0.05</td><td>-4.18 ± 0.06</td><td>0.37 ± 0.03</td><td>4.6E-34</td></tr><tr><td>730</td><td>55<br>(59)</td><td>1.35 ± 0.16<br>(1.17 ± 0.14)</td><td>-3.68 ± 0.15<br>(-3.41 ± 0.15)</td><td>0.33 ± 0.06</td><td>1.7E-02</td><td>90</td><td>1.57 ± 0.12</td><td>-4.25 ± 0.14</td><td>0.36 ± 0.04</td><td>1.5E-06</td></tr><tr><td>1100</td><td>26<br>(27)</td><td>1.06 ± 0.14<br>(1.02 ± 0.10)</td><td>-3.62 ± 0.23<br>(-3.54 ± 0.17)</td><td>0.28 ± 0.10</td><td>2.2E-02</td><td>43</td><td>1.42 ± 0.21</td><td>-4.17 ± 0.23</td><td>0.30 ± 0.07</td><td>4.3E-05</td></tr><tr><td>1460</td><td>13<br>(13)</td><td>1.18 ± 0.21<br>(1.18 ± 0.21)</td><td>-3.36 ± 0.15<br>(-3.36 ± 0.15)</td><td>0.20 ± 0.10</td><td>4.2E-02</td><td>22</td><td>1.53 ± 0.20</td><td>-4.06 ± 0.20</td><td>0.37 ± 0.06</td><td>1.1E-02</td></tr></tbody></table>

**Notes.**

<sup>a</sup> Size of each spaxel in pc. The size is derived from the square root of the inclination-corrected area of the spaxel.

<sup>b</sup> The SFR density Σ<sub>SFR</sub> is in units of M<sub>☉</sub> yr<sup>-1</sup> kpc<sup>-2</sup>. The gas surface density, both Σ<sub>H2</sub> and Σ<sub>H1+H2</sub>, is in units of M<sub>☉</sub> pc<sup>-2</sup>. Σ<sub>H1+H2</sub> is derived directly from Equation (4) and is directly proportional to the dust surface density.

<sup>c</sup> The best-fit slope γ<sub>H2</sub> and intercept A<sub>H2</sub> obtained by applying the OLS bisector linear fitting method to the N data in the Σ<sub>SFR</sub>-Σ<sub>H2</sub> scatter plots; N is the total number of spaxels used in the fit; σ<sub>H2</sub> is the scatter of the data about the best-fit line. The first line includes the best fits to the spaxels that have F(8)/F(24) ≤ 0.55 and molecular gas values above a threshold calculated from the combined uncertainty of the total gas measurement and H I detection limit (Section 7). The p<sub>τ</sub> values indicate the probability that the data are uncorrelated. On the second line for each region, in parentheses, we report also the number of spaxels used and the best-fit values derived by including spaxels that have Σ<sub>H2</sub> values as low as 1/10 of the combined uncertainty.

<sup>d</sup> The best-fit slope γ<sub>H1+H2</sub> and intercept A<sub>H1+H2</sub> obtained by applying the OLS bisector linear fitting method to the N data in the Σ<sub>SFR</sub>-Σ<sub>H1+H2</sub> scatter plots; N is the total number of spaxels used in the fit; σ<sub>H1+H2</sub> is the scatter of the data about the best-fit line. The p<sub>τ</sub> values indicate the probability that the data are uncorrelated. For the Σ<sub>SFR</sub>-Σ<sub>H1+H2</sub> fits we use all spaxels with F(8)/F(24) ≤ 0.55, since the total gas uncertainty is low enough that no spaxel is rejected.

Star formation is related to the molecular cloud mass as

$$SFR \propto M_{H2}^{\beta}, \qquad (9)$$

with β in the range [1, 2]. The relation between SFR and M<sub>H2</sub> has a built-in scatter expressed as a Gaussian in log–log space with a factor of 2 standard deviation. Changing this standard deviation to a factor of 4 or 0 has minimal impact on the results (Calzetti et al. 2012). The choice of a direct relation between SFR and molecular gas mass, instead of a relation between the volume densities of these quantities, stems from a model limitation: Calzetti et al. simulate clouds, but not their internal structure, which is necessary to obtain a relation between ρ<sub>SFR</sub> and ρ<sub>H2</sub> (Lada et al. 2010, 2012). However, as discussed in that paper, a power-law relation between volume densities translates into a power-law relation between SFR and gas mass. At the most basic level, β = 1 implies ρ<sub>SFR</sub> ∝ ρ<sub>H2</sub>. The data on molecular clouds within ~1 kpc of the Sun (Heiderman et al. 2010; Evans et al. 2014) imply SFR ∝ M<sub>H2</sub><sup>1.3</sup> and ρ<sub>SFR</sub> ∝ ρ<sub>H2</sub><sup>1.7</sup>.

The Schmidt–Kennicutt relation is then measured in the simulations following the same strategy as the observational approach, i.e., by deriving the parameters of Equations (7) and (8) from the simulated data. Calzetti et al. (2012) apply selection biases and detection “limits” that attempt to mimic as closely as possible typical observational conditions, including a dynamic range ≲1.5 in log(Σ<sub>H2</sub>). This is consistent with the dynamic range we obtain for NGC 4449, although, thanks to the IR/millimeter relations (Equations (4) and (5)), our range is about 0.5 dex wider than most previous studies.

We compare the data for NGC 4449 with simulations for β = 1 and β = 1.5, which bracket the range of observed slopes. We also consider a third model, with β = 1 and with a threshold for the star formation: the clouds in the lowest 18% bin, by mass, of the mass function do not form stars. This third model aims at simulating the case in which our method for removing a smooth background from the Hα and 24 µm images (Section 5), which averages to ~18% of the total

emission from the two images, artificially removes some of the true star formation.

### 8.3. Comparing Data with Simulations

The slopes and scatter predicted by the simulations are compared with the observed values for NGC 4449 as a function of increasing region size in Figure 16. The same algorithms are used to measure γ and σ in the simulated and actual data, to ensure that they can be compared. We use the OLS bisector fitting routine to measure the slope γ, and we fit the distributions of perpendicular distances of the data from the best-fit lines with Gaussians; we adopt the standard deviation σ of the Gaussian as our measure of the spread of the data. Among the three models described above, only the case of β = 1.5 consistently, albeit not perfectly, agrees with the observations. We do not attempt to improve the agreement between the β = 1.5 simulation and the data because of parameter degeneracy. The trends of both γ<sub>H2</sub> and σ<sub>H2</sub> can be made steeper toward smaller spaxel sizes by increasing (slightly) β, changing the slope of the cloud mass function, or increasing the maximum value of the molecular cloud mass by a factor of 2–3. The latter is within the regime of observed values for the LMC, which, as mentioned above, we take as being similar to NGC 4449.

The model comparisons in Figure 16 indicate that the observed trends for γ<sub>H2</sub> and σ<sub>H2</sub> are driven by the increasingly better sampling of the molecular cloud mass function at larger region sizes, as expected. They also suggest that the relation between SFR and cloud mass in NGC 4449 is superlinear, which implies a superlinear relation between volume densities of the same parameters. A linear relation between the SFR and cloud mass (or between ρ<sub>SFR</sub> and ρ<sub>H2</sub>) is excluded by the data at more than the 5σ level, from the combined trends of both γ<sub>H2</sub> and σ<sub>H2</sub>.

From Figure 16, the value of γ<sub>H2</sub> at 1100 pc is lower than what would be expected for a smooth trend between γ<sub>H2</sub> and region size, although still consistent within the uncertainties. A
<figure><img src="image_20.png" /></figure>

<figure><img src="image_21.png" /><figcaption>**Figure 16.** Slope, $\gamma_{H2}$, of the log($\Sigma_{SFR}$)–log($\Sigma_{H2}$) best-fit linear relation (left panel) and the scatter, $\sigma_{H2}$, of the data about the best-fit line (right panel), as a function of the linear size of the galactic regions in NGC 4449 used to perform the measurements. The data, with their $1\sigma$ uncertainties, are shown as black circles. Three simulations from Calzetti et al. (2012) are used for comparison, discriminated by the value of the exponent $\beta$ in the relation between SFR and molecular gas mass, SFR $\propto M_{H2}^{\beta}$; the case $\beta = 1.5$ (red crosses and line), the case $\beta = 1$ (blue crosses and line), and, finally, the case $\beta = 1$ with the lowest 18% of the mass in clouds being devoid of star formation (cyan line). Data from the literature are also shown in the panels: M51 and NGC 3521 are shown as dark- and light-gray circles, respectively, from Liu et al. (2011); magenta circles show the location of the data from Kennicutt et al. (2007, K+07, ~500 pc), Bigiel et al. (2008, B+08, ~750 pc), and Rahman et al. (2012, R+12, ~1000 pc).</figcaption></figure>

close inspection of the data and fits does not reveal any unusual circumstance, beyond the known sensitivity of the OLS bisector fitting algorithm to outliers at the edges of the dynamic range of the data. Data on other galaxies will be required to fully understand this deviation.

Our observed trends for the relation between $\Sigma_{SFR}$ and $\Sigma_{H2}$ in NGC 4449 are in agreement with those found by Liu et al. (2011). These authors analyze two star-forming galaxies, NGC 3521 and NGC 5194 (M51), and find that both the power-law exponent and the scatter about the best-fit lines of the molecular gas–SFR relation are a decreasing function of increasing region size, between 250 pc and 1.3 kpc. We report the data of Liu et al. (2011) in Figure 16. The data for M51 are systematically higher than those for NGC 4449, while the data for the slope $\gamma_{H2}$ of NGC 3521 are roughly consistent with those of our galaxy. While M51 is an almost face-on grand design spiral, NGC 3521 is highly inclined. The methodology of Liu et al. (2011) consists of identifying the H II regions before removing the diffuse stellar light from older stellar populations, in order to measure the SFRs; this method is subject to line-of-sight confusion for highly inclined galaxies. We also report in Figure 16 the results of Kennicutt et al. (2007), Bigiel et al. (2008), and Rahman et al. (2012). These authors only measure a single region size, common to all galaxies, as specified in the caption of Figure 16. Kennicutt et al. (2007) and Rahman et al. (2012) obtain results for $\gamma_{H2}$ (and, in the case of Kennicutt et al. [2007], also for $\sigma_{H2}$) that are consistent with our results for NGC 4449, once we take into account the appropriate region sizes for the measurements. These authors use a methodology for removing the diffuse stellar emission from the SFR indicators that are the closest to our method. Conversely, the data of Bigiel et al. (2008) are the most discrepant from our results; we speculate that the reason for the discrepancy is that the fluxes used to derive the SFR indicators in that work

include the contribution of diffuse emission from old stellar populations, which may yield artificially low values of the slopes and scatter. This effect is quantified in Calzetti et al. (2012). The data from the spirals show both similarities and differences with our observed trends for $\gamma_{H2}$ and $\sigma_{H2}$. Some of the discrepancies could be intrinsic, but others could be due to differences in the measurement techniques. All results from the literature discussed here, for instance, use CO emission to measure the molecular gas mass and do not include CO-dark molecular gas (Pineda et al. 2013). Conversely, our dust-based measurements include the CO-dark H$_2$. Thus, meaningful comparisons between galaxies will first need to address the issue of getting homogeneous measurements.

Our main conclusion from the analysis of Figure 16 is that trends of $\gamma_{H2}$ and $\sigma_{H2}$ with region size are present in both the current and previous data but had not been recognized before Liu et al. (2011).

In contrast, the slope $\gamma_{H1+H2}$ between the surface densities of the SFR and total gas remains roughly constant as a function of region size (Figure 17) and is in agreement with the slope observed for whole galaxies (Kennicutt 1998; Kennicutt & Evans 2012) and for the star-forming regions in M51 (Kennicutt et al. 2007). The spread of the data about the best-fit lines, $\sigma_{H1+H2}$, also remains roughly constant with region size (Table 4) and comparable in value to the largest values observed for $\sigma_{H2}$; this suggests that the spread in H I +H2 is driven predominantly by the scatter in the H I component. The roughly constant values of both $\gamma_{H1+H2}$ and $\sigma_{H1+H2}$, when combined with the high significance of the $\Sigma_{SFR}-\Sigma_{H1+H2}$ correlations, suggest that the onset of the superlinear trend between the SFR and the total gas surface densities (the SK law) occurs already at the subgalactic scale and that the atomic gas is an integral part in establishing the scaling between gas and SFR. In other words, the process
<figure><img src="image_22.png" /><figcaption>**Figure 17.** Slope, $\gamma_{H1+H2}$, of the log($\Sigma_{SFR}$)-log($\Sigma_{H1+H2}$) best-fit linear relation as a function of the size of the galactic regions in NGC 4449. The relation observed for the total (atomic+molecular) gas remains steep, around a slope of $\approx 1.5$, for all region sizes, and it agrees, within the uncertainties, with the slope of $\sim 1.4$ measured for whole galaxies (red horizontal line; Kennicutt 1998; Kennicutt & Evans 2012). Additional data from the literature are shown as magenta circles: Kennicutt et al. (2007, K+07, $\sim 500$ pc) and Bigiel et al. (2008, B+08, $\sim 750$ pc).</figcaption></figure>

that determines the balance between the atomic and molecular gas phases and between the gas and the SFR is established already at scales of a few hundred parsecs, much larger than the typical size of a single molecular cloud, but much smaller than the size of a galaxy.

Bigiel et al. (2008) find the value of $\gamma_{H1+H2}$ to vary considerably from galaxy to galaxy, with a range of 1.1–2.7, when performing their analysis using region sizes of 750 pc in each galaxy (Figure 17). The authors’ interpretation is that since H I is not directly related to star formation, its influence on $\gamma_{H1+H2}$ is to make it a galaxy-dependent measure. This interpretation appears in contradiction with our result that $\gamma_{H1+H2}$ is independent of the region’s size in NGC 4449. In addition to the already-mentioned difference in the method to derive SFR indicators, the authors measure the molecular mass surface density directly from CO (2–1). We speculate that the large range of values they find may arise from the combination of a region-dependent contribution of the old stellar population to the SFRs, spatially variable CO (2–1)/CO (1–0) ratios (Koda et al. 2012), and uncertain CO-to-H2 conversion values. As already noted above, meaningful comparisons between different samples and galaxies will first require that homogeneous measurements are performed across all diagnostics. A larger sample of galaxies observed in the cold dust regime, thus capturing the entire gas component with one measurement technique, is needed, in order to confirm or refute our speculations.

Hopkins et al. (2014) and Orr et al. (2017) perform high-resolution cosmological simulations of star-forming galaxies, finding that the Schmidt–Kennicutt relation arises naturally as an effect of feedback on local scales within the galaxies, irrespective of the specific star formation prescription used. The local stellar feedback, including supernovae, stellar winds,

radiation pressure, and photo heating, generates turbulence and acts as a self-regulating mechanism, which not only keeps star formation inefficient but also imparts a scaling between gas (both total and molecular) and SFR. Our results tend to support this scenario: in NGC 4449 we observe a balance between total gas and SFR, which occurs at all scales probed, from $\sim 350$ pc to 1.5 kpc. A scale-free mechanism like turbulence would be able to account for this result, although there may be additional favorable conditions unique to the starburst environment of NGC 4449, such as a high filling factor for star formation. This has been shown to be a common feature of starbursts, as opposed to low filling factors for normal star-forming disks (Elmegreen et al. 2014). The intense, localized star formation is likely to favor a scenario in which feedback-induced turbulence reaches equilibrium already at small scales of a few hundred parsecs. Observations of larger, and more varied, samples of galaxies will enable testing this and other models.

## 9. Summary and Conclusions

We have combined 1.1 mm (=1100 $\mu$m) maps of the nearby starburst galaxy NGC 4449 from the LMT/AzTEC with images at infrared wavelengths from the *Spitzer Space Telescope* and the *Herschel Space Observatory* to investigate the distribution and properties of the dust in the central $\sim 4.7 \times 12.5$ kpc$^2$ active region. Our main findings are as follows:

1. The dust emission in this dwarf has similar characteristics to those of other nearby starburst dwarfs, including a relatively high value of $U_{min} \sim 4–5$, indicating a relatively high mean dust temperature (Draine et al. 2007). We find that the dust masses derived from the Draine & Li (2007) models require a downward correction by a factor of $\lesssim 1.5$. This is in line with what has been found by the *Planck Observatory* for the dust in our MW, for which correction factors around 1.8–2.0 have been determined (Fanciullo et al. 2015; Planck Collaboration 2016b). The dust in NGC 4449 appears to require a smaller correction factor than the MW, likely a reflection of its higher mean starlight intensity.

2. The mid-IR ($F(8)/F(24)$) color is correlated with the far-IR/millimeter ($F(70)/F(1100)$ and $F(160)/F(1100)$) colors, for regions of sizes in the range from $\sim 360$ pc ($11'' \times 11''$) to $\sim 1.5$ kpc ($44'' \times 44''$) in NGC 4449. This IR-millimeter correlation, which is recovered also for whole galaxies, enables us to increase the dynamic range of our dust mass determinations. We use high-S/N data (mid-IR and far-IR data close to the peak of infrared emission) to derive the luminosity at 1100 $\mu$m for spaxels that are undetected or marginally detected by AzTEC. This extends our dust mass derivations down to a one order of magnitude fainter limit in 1100 $\mu$m luminosity than would be possible by using only the millimeter data points.

3. The total (H I+H$_2$) and molecular (H$_2$) mass surface densities, derived from the combination of dust and H I mass surface densities, are compared with the SFR surface density as a function of region size. The $\Sigma_{SFR}-\Sigma_{H2}$ scatter plots can be fit by power laws with exponents that decrease from $\sim 1.5$ to $\sim 1.2$ for increasing region size from 360 pc to 1.5 kpc; the same decreasing trend is observed for the scatter of the data about the
best-fit lines. The $\Sigma_{SFR}-\Sigma_{H2}$ trend with region size is consistent with random sampling of both molecular clouds and star-forming regions within randomly selected galactic regions (Calzetti et al. 2012). The presence of such trends had already been recognized by Liu et al. (2011).

4. Conversely, the power laws describing the $\Sigma_{SFR}-\Sigma_{H1+H2}$ scatter plots show a constant exponent with a value of ~1.5 and a constant scatter about the best-fit lines, independent of region size. The dominant factor in the spread of the data is likely the scatter in the H1 component. We also find that the $\Sigma_{SFR}-\Sigma_{H1+H2}$ correlation is significantly stronger than the $\Sigma_{SFR}-\Sigma_{H2}$ correlation. We suggest that the constant slope reflects the fact that both gas components play a role in determining the scaling with the SFR in this H1-dominated galaxy, and their balance is established at galactic scales that are larger than those of individual clouds, but considerably smaller than those of the whole galaxy. This is in agreement with a scenario in which the equilibrium between gas and SFR is set by feedback-induced turbulence (Hopkins et al. 2014; Orr et al. 2017).

The results for NGC 4449, albeit tantalizing, are still limited by having a single target. Analyses of larger samples of nearby galaxies, conducted with either the LMT or other millimeter facilities, will be required to place the results from this work on a more solid footing. The IR/millimeter color-color relation is of particular interest. If confirmed and calibrated for a significant number of star-forming galaxies and regions within galaxies, it will provide a venue for expanding the dynamic range of dust mass determinations. By using high-luminosity portions of the infrared SEDs of galaxies to infer their millimeter wavelength emission (which is intrinsically faint), this relation would increase the dynamic range of dust mass and surface density determinations by more than an order of magnitude relative to what is possible with millimeter data alone.

The authors would like to thank the anonymous referee for many constructive comments that have helped significantly improve the presentation of this paper.

Based on observations made with the Large Millimeter Telescope Alfonso Serrano, a binational project between Mexico and the USA, led by the Instituto Nacional de Astrofísica, Optica, y Electrónica (INAOE) and by the University of Massachusetts at Amherst, respectively.

This research has made use of the NASA/IPAC Extragalactic Database (NED), which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration.

This work is based in part on data downloaded from the *Spitzer* Heritage Archive. The *Spitzer* Space Telescope is operated by the Jet Propulsion Laboratory, California Institute of Technology, under a contract with NASA.

This work used data obtained from the ESA *Herschel* Science Archive. *Herschel* was an ESA space observatory with science instruments provided by European-led Principal Investigator consortia and with important participation from NASA.

Facility: LMT(AzTEC).

# Appendix

## The Cottingham Method Applied to AzTEC Data

In this appendix, we summarize the Cottingham (Cottingham 1987) method as applied to the AzTEC data.

The Cottingham method is a fitting approach that results in a maximum likelihood estimate of parameters in a model of the temporally varying atmospheric contamination. The approach begins with the despiked and calibrated time stream detector signals from each of AzTEC's detectors. The median value of each individual time stream is estimated and subtracted, an additional linear gain correction is applied to remove a long-term time drift—observed as variable tilt in the raw time streams—and then the atmosphere coefficients, $\alpha$, for a B-spline basis matrix, $B$ (the atmosphere model), are found by solving the matrix equation

$$AB\alpha = Ad^{\text{row}}, \quad (10)$$

where the matrix $A$ is defined as

$$A = B^T(1 - P\Pi). \quad (11)$$

$\Pi(P^T P)^{-1} P^T$ is known as the “pointing matrix,” and $d^{\text{row}} = (d_1, d_2, \dots, d_{\text{nbolo}})$ is a row vector holding the concatenation of the different detectors. These equations are written assuming uncorrelated noise among the detectors, i.e., $N = I$, where $N$ is the noise covariance matrix.

Because of AzTEC's relatively small number of detectors, and given that key detectors near the center of the array are inactive, the pointing matrix is not directly invertible for most of the data. To solve this, we perform a first-order polynomial fit to the detector time stream residuals in azimuthal coordinates of the form

$$P_j^i(a, e) = A + b a_j^i + C e_j^i, \quad (12)$$

where $a_j^i$ and $e_j^i$ are the azimuth and elevation coordinates of detector $i$ at sample $j$, respectively. This is equivalent to fitting an instantaneous plane to the array sky brightness in azimuth and elevation to each set of detector samples. Finally, the atmosphere-cleaned time streams are calculated by subtracting both the B-spline atmosphere template and the polynomial fit from the observed raw time streams. The reduction process then continues without modification from the standard AzTEC pipeline approach (e.g., Scott et al. 2012, and references therein).

### A.1. Iterative Flux Recovery

Although the Cottingham method is an unbiased surface brightness estimator, some of the reduction steps, such as the median estimate, the linear gain, and the polynomial fit, are susceptible to introducing some degree of bias. In order to recover a map that better represents the astronomical surface brightness distribution, the Flux Recovery Using Iterations Technique was used in a similar fashion to that described in Liu et al. (2010); the main difference relative to that paper is that we use the Cottingham method instead of the standard PCA-based atmospheric removal approach, together with some other minor modifications to the pixel selection criteria.

## ORCID iDs

D. Calzetti https://orcid.org/0000-0002-5189-8004

G. W. Wilson https://orcid.org/0000-0003-2705-9152
B. T. Draine [https://orcid.org/0000-0002-0846-936X](https://orcid.org/0000-0002-0846-936X)

H. Roussel [https://orcid.org/0000-0001-5617-7129](https://orcid.org/0000-0001-5617-7129)

K. E. Johnson [https://orcid.org/0000-0001-8348-2671](https://orcid.org/0000-0001-8348-2671)

W. F. Wall [https://orcid.org/0000-0001-6142-397X](https://orcid.org/0000-0001-6142-397X)

K. Grasha [https://orcid.org/0000-0002-3247-5321](https://orcid.org/0000-0002-3247-5321)

A. Battisti [https://orcid.org/0000-0003-4569-2285](https://orcid.org/0000-0003-4569-2285)

A. Kirkpatrick [https://orcid.org/0000-0002-1306-1545](https://orcid.org/0000-0002-1306-1545)

D. Rosa González [https://orcid.org/0000-0003-1327-0838](https://orcid.org/0000-0003-1327-0838)

O. Vega [https://orcid.org/0000-0002-2852-9737](https://orcid.org/0000-0002-2852-9737)

J. Puschnig [https://orcid.org/0000-0003-1111-3951](https://orcid.org/0000-0003-1111-3951)

M. Yun [https://orcid.org/0000-0001-7095-7543](https://orcid.org/0000-0001-7095-7543)

G. Östlin [https://orcid.org/0000-0002-3005-1349](https://orcid.org/0000-0002-3005-1349)

A. S. Evans [https://orcid.org/0000-0003-2638-1334](https://orcid.org/0000-0003-2638-1334)

Y. Tang [https://orcid.org/0000-0002-6581-3307](https://orcid.org/0000-0002-6581-3307)

J. Lowenthal [https://orcid.org/0000-0001-9969-3115](https://orcid.org/0000-0001-9969-3115)

D. Sánchez-Arguelles [https://orcid.org/0000-0002-7344-9920](https://orcid.org/0000-0002-7344-9920)

## References

Accurso, G., Saintonge, A., Catinella, B., et al. 2017, MNRAS, 470, 4750

Andrews, J. E., Calzetti, D., Chandar, R., et al. 2013, ApJ, 767, 51

Aniano, G., Draine, B. T., Calzetti, D., et al. 2012, ApJ, 756, 138

Aniano, G., Draine, B. T., Gordon, K. D., & Sandstrom, K. 2011, PASP, 123, 1218

Asplund, M., Grevesse, N., Sauval, A. J., & Scott, P. 2009, ARA&A, 47, 481

Austermann, J. E., Dunlop, J. S., Perera, T. A., et al. 2010, MNRAS, 401, 160

Berg, D. A., Skillman, E. D., Marble, A. R., et al. 2012, ApJ, 754, 98

Bigiel, F., Leroy, A. D., Walter, F., et al. 2008, AJ, 136, 2846

Bolatto, A. D., Wolfire, M., & Leroy, A. K. 2013, ARA&A, 51, 207

Böttner, C., Klein, U., & Heithausen, A. 2003, A&A, 408, 493

Calapa, M. D., Calzetti, D., Draine, B. T., et al. 2014, ApJ, 784, 130

Calzetti, D. 2013, in *Secular Evolution of Galaxies*, ed. J. Falcón-Barroso & J. H. Knapen (Cambridge: Cambridge Univ. Press), 419

Calzetti, D., Kennicutt, R. C., Bianchi, L., et al. 2005, ApJ, 633, 871

Calzetti, D., Kennicutt, R. C., Engelbracht, C. W., et al. 2007, ApJ, 666, 870

Calzetti, D., Lee, J. C., Sabbi, E., et al. 2015, AJ, 149, 51

Calzetti, D., Liu, G., & Koda, J. 2012, ApJ, 752, 98

Calzetti, D., Wu, S.-Y., Hong, S., et al. 2010, ApJ, 714, 1256

Cerviño, M., Valls-Gabaud, D., Luridiana, V., & Mas-Hesse, J. M. 2002, A&A, 381, 51

Cook, D. O., Dale, D. A., Johnson, B. D., et al. 2014, MNRAS, 445, 899

Cottingham, D. A. 1987, PhD thesis, Princeton Univ.

Daddi, E., Elbaz, D., Wlater, F., et al. 2010, ApJL, 714, L118

Dale, D. A., Cohen, S. A., Johnson, L. C., et al. 2009, ApJ, 703, 517

Dale, D. A., Cook, D. O., Roussel, H., et al. 2017, ApJ, 837, 90

Draine, B. T., Aniano, G., Krause, O., et al. 2014, ApJ, 780, 172

Draine, B. T., Dale, D. A., Bendo, G., et al. 2007, ApJ, 663, 866

Draine, B. T., & Li, A. 2007, ApJ, 657, 810

Elmegreen, D. M., Elmegreen, B. G., Adamo, A., et al. 2014, ApJL, 787, L15

Evans, N. J., Heiderman, A., & Vutisalchavakul, N. 2014, ApJ, 782, 114

Faesi, C. M., Lada, C. J., Forbrich, J., Menten, K. M., & Bouy, H. 2014, ApJ, 789, 81

Fanciullo, L., Guillet, V., Aniano, G., et al. 2015, A&A, 580A, 136

Feigelson, E. D., & Jogesh Babu, G. 2012, Modern Statistical Methods for Astronomy (Cambridge: Cambridge Univ. Press)

Galametz, M., Kennicutt, R. C., Albrecht, M., et al. 2012, MNRAS, 425, 763

Gao, Y., & Solomon, P. M. 2004, ApJ, 606, 271

Genzel, R., Tacconi, L. J., Gracia-Carpio, J., et al. 2010, MNRAS, 407, 2091

Genzel, R., Tacconi, L. J., Lutz, D., et al. 2015, ApJ, 800, 20

Grasha, K., Calzetti, D., Andrews, J. E., Lee, J. C., & Dale, D. A. 2013, ApJ, 773, 174

Gutermuth, R. A., Pipher, J. L., Megeath, S. T., et al. 2011, ApJ, 739, 84

Hao, C.-N., Kennicutt, R. C., Johnson, B. D., et al. 2011, ApJ, 741, 124

Heiderman, A., Evans, N., Allen, L. E., Huard, T., & Heyer, M. 2010, ApJ, 723, 1019

Helou, G., Roussel, H., Appleton, P., et al. 2004, ApJS, 154, 253

Hincks, A. D., Acquaviva, V., Ade, P. A. R., et al. 2010, ApJS, 191, 423

Hony, S., Gouliermis, D. A., Galliano, F., et al. 2015, MNRAS, 448, 1847

Hopkins, P. F., Keres, D., Óñorbe, J., et al. 2014, MNRAS, 445, 581

Huchtmeier, W. K., & Richter, O.-G. 1989, A General Catalog of HI Observations of Galaxies (New York: Springer)

Hughes, A., Wong, T., Ott, J., et al. 2010, MNRAS, 406, 2065

Hunter, D. A., Gallagher, J. S., Rice, W. L., & Gillett, F. C. 1989, ApJ, 336, 152

Hunter, D. A., Rubin, V. C., Swaters, R. A., Sparke, L. S., & Levine, S. E. 2002, ApJ, 580, 194

Hunter, D. A., Rubin, V. C., Swaters, R. A., Sparke, L. S., & Levine, S. E. 2005, ApJ, 634, 281

Hunter, D. A., van Woerden, H., & Gallagher, J. S. 1999, AJ, 118, 2184

Hunter, D. A., Wilcots, E. M., van Woerden, H., Gallagher, J. S., & Kohle, S. 1998, ApJL, 495, L47

Isobe, T., Feigelson, E. D., Akritas, M. G., & Babu, G. J. 1990, ApJ, 364, 10

Ivezic, Z., Connolly, A. J., VanderPlas, J. T., & Gray, A. 2014, Statistics, Data Mining, and Machine Learning in Astronomy (Princeton, NJ: Princeton Univ. Press)

Jameson, K. E., Bolatto, A. D., Leroy, A. K., et al. 2016, ApJ, 825, 12

Karachentsev, I. D., Sharina, M. E., Dolphin, A. E., et al. 2003, A&A, 398, 467

Karczewski, O. L., Barlow, M. J., Page, M. J., et al. 2013, MNRAS, 431, 2493

Kawamura, A., Mizuno, Y., Minamidani, T., et al. 2009, ApJS, 184, 1

Kendall, S., Kennicutt, R. C., Clarke, C., & Thornley, M. D. 2008, MNRAS, 387, 1007

Kennicutt, R. C., Armus, L., Bendo, G., et al. 2003, PASP, 115, 928

Kennicutt, R. C., Calzetti, D., Aniano, G., et al. 2011, PASP, 123, 1347

Kennicutt, R. C., Calzetti, D., Walter, F., et al. 2007, ApJ, 671, 333

Kennicutt, R. C., & Evans, N. J. 2012, ARA&A, 50, 531

Kennicutt, R. C., Lee, J. C., Funes, S. J., Sakai, S., & Akiyama, S. 2008, ApJS, 178, 247

Kennicutt, R. C., Jr. 1998, ApJ, 498, 541

Koda, J., Scoville, N., Hasegawa, T., et al. 2012, ApJ, 761, 41

Kohle, S. 1999, PhD thesis, Univ. Bonn

Kohle, S., Klein, U., Henkel, C., & Hunter, D. A. 1998, in *The Magellanic Clouds and Other Dwarf Galaxies*, ed. T. Richtler & J. M. Braun (Aachen: Shaker Verlag), 265

Kroupa, P. 2001, MNRAS, 322, 231

Kruijssen, J. M. D., & Longmore, S. N. 2014, MNRAS, 439, 3239

Lada, C. J., Forbrich, J., Lombardi, M., & Alves, J. F. 2012, ApJ, 745, 190

Lada, C. J., Lombardi, M., & Alves, J. F. 2010, ApJ, 724, 687

Lada, C. J., Lombardi, M., Roman-Zuniga, C., Forbrich, J., & Alves, J. F. 2013, ApJ, 778, 133

Larson, R. B. 1981, MNRAS, 194, 809

Lee, J. C., Gil de Paz, A., Kennicutt, R. C., et al. 2011, ApJS, 192, 6

Lee, J. C., Gil de Paz, A., Tremonti, C., et al. 2009, ApJ, 706, 599

Lelli, F., Verheijen, M., & Fraternali, F. 2014, MNRAS, 445, 1694

Leroy, A., Bolatto, A. D., Simon, J. D., & Blitz, L. 2005, ApJ, 625, 763

Leroy, A. K., Walter, F., Sandstrom, K., et al. 2013, AJ, 146, 19

Li, Y., Crocker, A. F., Calzetti, D., et al. 2013, ApJ, 768, 180

Liu, G., Calzetti, D., Yun, M. S., et al. 2010, AJ, 139, 1190

Liu, G., Koda, J., Calzetti, D., Fukuhara, M., & Momose, R. 2011, ApJ, 735, 63

Meidt, S. E., Schinnerer, E., Knapen, J. H., et al. 2012, ApJ, 744, 17

Mogotsi, K. M., de Blok, W. J. G., Caldú-Primo, A., et al. 2016, AJ, 151, 15

Momose, R., Okumura, S. K., Koda, J., & Sawada, T. 2010, ApJ, 721, 383

Oey, M. S., Meurer, G. R., Yelda, S., et al. 2007, ApJ, 661, 801

Onodera, S., Kuno, N., Tosaki, T., et al. 2010, AJL, 722, L127

Orr, M. E., Hayward, C. C., Hopkins, P. F., et al. 2017, MNRAS, submitted (arXiv:1701.01788)

Peng, C. Y., Ho, L. C., Impey, C. D., & Rix, H.-W. 2002, AJ, 124, 266

Peng, C. Y., Ho, L. C., Impey, C. D., & Rix, H.-W. 2010, AJ, 139, 2097

Perera, T. A., Chapin, E. L., Austermann, J. E., et al. 2008, MNRAS, 391, 1227

Pilyugin, L. S., Grebel, E. K., & Zinchenko, I. A. 2015, MNRAS, 450 3254

Pineda, J. L., Langer, W. D., Velusamy, T., & Goldsmith, P. F. 2013, A&A, 554A, 103

Planck Collaboration, Ade, P. A. R., Aghanim, N., Argüeso, F., et al. 2016a, A&A, 594A, 26

Planck Collaboration, Ade, P. A. R., Aghanim, N., Alves, M. I. R., et al. 2016b, A&A, 586A, 132

Rahman, N., Bolatto, A. D., Xue, R., et al. 2012, ApJ, 745, 83

Reines, A. E., Johnson, K. E., & Goss, W. M. 2008, AJ, 135, 2222

Rémy-Ruyer, A., Madden, S. C., Galliano, F., et al. 2014, A&A, 563A, 31

Rémy-Ruyer, A., Madden, S. C., Galliano, F., et al. 2015, A&A, 582A, 121

Rojas Bolivar, R., Calzetti, D., Dale, D. A., & Cook, D. O. 2016, AAS Meeting #227, 242.03

Roman-Duval, J., Bot, C., Chastenet, J., & Gordon, K. D. 2017, ApJ, 841, 72

Roman-Duval, J., Gordon, K. D., Meixner, M., et al. 2014, ApJ, 797, 86

Roussel, H. 2013, PASP, 125, 1126

Saintonge, A., Kauffmann, G., Kramer, C., et al. 2011, MNRAS, 415, 32
Sandstrom, K. M., Bolatto, A. D., Draine, B. T., Bot, C., & Stanimirović, S. 2010, ApJ, 715, 701

Schruba, A., Leroy, A. K., Walter, F., Sandstrom, K., & Rosolowsky, E. 2010, ApJ, 722, 1699

Scott, K. S., Austermann, J. E., Perera, T. A., et al. 2008, MNRAS, 385 2225

Scott, K. S., Wilson, G. W., Aretxaga, I., et al. 2012, MNRAS, 423, 575

Shetty, R., Kelly, B. C., & Klessen, R. S. 2014a, MNRAS, 442, 2208

Shetty, R., Kelly, B. C., Rahman, N., et al. 2014b, MNRAS, 437, L61

Smith, J. D. T., Draine, B. T., Dale, D. A., et al. 2007, ApJ, 656, 770

Srivastava, S., Kantharia, N. G., Basu, A., Srivastava, D. C., & Ananthakrishnan, S. 2014, MNRAS, 443, 860

Tully, R. B., Courtois, H. M., Dolphin, A. E., et al. 2013, AJ, 146, 86

Wall, W. F., Puerari, I., Tilanus, R., et al. 2016, MNRAS, 459, 1440

Walter, F., Brinks, E., de Blok, W. J. G., et al. 2008, AJ, 136, 2563

Whitaker, K. E., van Dokkum, P. G., Brammer, G., & Franx, M. 2012, ApJL, 754, L29

Whitmore, B. C., Brogan, C., Chandar, R., et al. 2014, ApJ, 795, 156

Wilson, G. W., Austermann, J. E., Perera, T. A., et al. 2008, MNRAS, 386 807