Gas kinematics in FIRE simulated galaxies compared to spatially unresolved H i observations

Abstract

The shape of a galaxy's spatially unresolved, globally integrated 21-cm emission line depends on its internal gas kinematics: galaxies with rotationally supported gas discs produce double-horned profiles with steep wings, while galaxies with dispersion-supported gas produce Gaussian-like profiles with sloped wings. Using mock observations of simulated galaxies from the FIRE project, we show that one can therefore constrain a galaxy's gas kinematics from its unresolved 21-cm line profile. In particular, we find that the kurtosis of the 21-cm line increases with decreasing V/σ and that this trend is robust across a wide range of masses, signal-to-noise ratios, and inclinations. We then quantify the shapes of 21-cm line profiles from a morphologically unbiased sample of ∼2000 low-redshift, H i-detected galaxies with M_star = 10⁷⁻¹¹ M_⊙ and compare to the simulated galaxies. At M_star ≳ 10¹⁰ M_⊙, both the observed and simulated galaxies produce double-horned profiles with low kurtosis and steep wings, consistent with rotationally supported discs. Both the observed and simulated line profiles become more Gaussian like (higher kurtosis and less-steep wings) at lower masses, indicating increased dispersion support. However, the simulated galaxies transition from rotational to dispersion support more strongly: at M_star = 10⁸⁻¹⁰ M_⊙, most of the simulations produce more Gaussian-like profiles than typical observed galaxies with similar mass, indicating that gas in the low-mass simulated galaxies is, on average, overly dispersion supported. Most of the lower-mass-simulated galaxies also have somewhat lower gas fractions than the median of the observed population. The simulations nevertheless reproduce the observed line-width baryonic Tully–Fisher relation, which is insensitive to rotational versus dispersion support.

galaxies: dwarf, galaxies: irregular, galaxies: kinematics and dynamics

1 INTRODUCTION

The shape of a galaxy's globally integrated 21-cm line profile encodes information about the kinematics of neutral hydrogen (H i) in its interstellar medium (ISM). Because the 21-cm transition is usually optically thin, an observed line profile can be directly translated into a line-of-sight velocity distribution (LOSVD) of H i gas in the galaxy. The shape of the 21-cm line thus depends on a galaxy's internal kinematics, and it will in general be different for discy, rotationally supported galaxies than for irregular galaxies in which gas is primarily supported by dispersion.

Most large spiral galaxies produce characteristic ‘double-horned’ unresolved 21-cm line profiles, which arise naturally from differential Doppler effects in rotating discs (Epstein 1964; Shostak 1977). In fact, double-horned line profiles are not unique to 21-cm emission from H i in galaxies but are a generic consequence of rotating disc geometries and are produced in a wide range of astrophysical settings, including lines formed in accretion discs in some broad-line AGN (Stella 1990; Storchi-Bergmann et al. 2017), protostellar discs (Yamada et al. 2009), and accretion discs around massive stars (Doazan 1965) and white dwarfs (Tylenda 1981). In many cases, the precise shape of such double-peaked emission lines can be used to constrain the rotation curve and/or spatial distribution of the emitting material.

Unlike large spiral galaxies, many low-mass galaxies do not produce double-horned profiles, instead exhibiting single-peaked, Gaussian-like profiles with sloping wings (e.g. Geha et al. 2006; Begum et al. 2008; Haynes et al. 2011; Bradford, Geha & Blanton 2015). Single-peaked, Gaussian-like profiles are typically interpreted as evidence of lower-rotational velocities and more disordered gas kinematics than is typical for more massive disc galaxies (Singhal 2008; Stewart, Blyth & de Blok 2014; Papastergis, Adams & van der Hulst 2016). This suggests that in the absence of spatially resolved observations, the shape of a galaxy's unresolved H i profile can be used as a rough proxy for V/σ (i.e. the ratio of ordered rotation to dispersion) of the galaxy's gas.

In contrast to resolved H i kinematic maps, which are expensive to obtain and have only been produced of the order of 100 galaxies in the local Universe (e.g. Walter et al. 2008; Hunter et al. 2012; Ott et al. 2012), unresolved global 21-cm line profiles can be efficiently obtained with short, single-dish observations. Unresolved H i profiles are already cataloged for a few ×  10⁴ low-redshift galaxies (e.g. Paturel et al. 2003; Meyer et al. 2004; Begum et al. 2008; Haynes et al. 2011; Giovanelli & Haynes 2015), and next-generation radio surveys are expected to increase the size of unresolved H i samples by a factor of 10 within a decade (Dewdney et al. 2009; Duffy et al. 2012; Li et al. 2018).

Besides their availability, another advantage of using unresolved 21-cm line profiles as probes of H i kinematics is that most unresolved surveys are morphologically blind, providing a homogeneous and unbiased sample of the local gas-rich galaxy population against which to test the predictions of theoretical models. This is in contrast to spatially resolved studies, which often preferentially target ‘well-behaved’ disc galaxies. Mock unresolved H i observations are also straightforward to produce from simulations with minimal post-processing and can be compared directly to observations, without the need to fit kinematic models.

Both observational (McGaugh et al. 2000; Haynes et al. 2011; Bradford, Geha & Blanton 2015; Bradford, Geha & van den Bosch 2016) and theoretical (Obreschkow et al. 2009a,b; Brook & Shankar 2016; Macciò et al. 2016; Brooks et al. 2017) works have previously interpreted H i profiles from unresolved data. However, these works have typically not used the shape of line profiles to extract kinematic information, but have instead quantified observed 21-cm lines with only two numbers: the integrated flux, which can be translated to total H i mass, and the line width, which to first order traces the depth of the gravitational potential. These parameters are not sensitive to whether the emitting gas is supported primarily by rotation or dispersion.

In this paper, we compare unresolved 21-cm line profiles from simulations to profiles of real galaxies, with the goal of comparing the simulated and observed galaxies’ H i kinematics. We first show how the information encoded in the shapes of unresolved 21-cm line profiles can be leveraged to statistically constrain the gas kinematics of the low-redshift galaxy population. We use mock H i profiles produced from simulated galaxies, for which the internal gas kinematics are known a priori, as a guide for interpreting the line profiles of real galaxies. We then compare the line profiles of real galaxies to the mock profiles in order to compare the degree of rotational versus dispersion support in the simulated galaxies to that in the observed population.

The remainder of this paper is organized as follows. In Section 2, we describe our procedure for mock-observing simulated galaxies and show how unresolved H i line shape is related to internal rotational versus dispersion support. In Section 3, we quantify the shapes of observed H i line profiles and compare to our simulated galaxies. We summarize our results and discuss their implications in Section 4. In Appendix A, we describe a simple analytic model for predicting the shape of unresolved H i line profiles. We discuss the relation between rotational/dispersion support and other galaxy properties in Appendix B.

2 MOCK OBSERVATIONS

We analyse the same suite of cosmological hydrodynamic zoom-in simulations from the FIRE project (Hopkins et al. 2014) that was studied by El-Badry et al. (2018, hereafter E18). The simulations were run with the GIZMO hydrodynamics code (Hopkins 2015) in the Lagrangian ‘meshless finite mass’ (MFM) mode, using the FIRE-2 model for galaxy formation and feedback¹ (Hopkins et al. 2017). For details regarding the physical processes modelled in these simulations and their numerical implementation, we refer the reader to Hopkins et al. (2017, 2018).

This simulation suite contains 24 isolated galaxies, with z = 0 stellar masses spanning M_star = 10^6.3−11.1 M_⊙, residing in haloes with z = 0 masses M_200m = 10^9.9−12.2 M_⊙. Here, M_200m is the total mass within R_200m, which is the radius within which the matter density is 200× the mean matter density at z = 0. E18 studied the gas properties of these galaxies in detail, showing that they exhibit a diverse range of gas kinematics and morphologies: thin, rotationally supported disc galaxies dominate at M_star ≳ 10¹⁰ M_⊙, while puffy, dispersion-supported irregular galaxies dominate at lower masses. The properties of the individual simulated galaxies studied in this work are summarized in E18 (their table 1).

We mock-observe these simulated galaxies in order to produce synthetic unresolved H i profiles that can be compared to observations. The observational sample to which we compare is the homogeneously reduced sample of unresolved H i observations produced by Bradford et al. (2015, hereafter B15); our mock observations are designed to emulate these data. Most of the lower-mass galaxies (M_star ≲ 10^8.5 M_⊙) in the observational sample were observed with a typical velocity resolution of Δ|$v$| = 0.65 km s⁻¹ as part of the low-mass H i programme introduced in Geha et al. (2006); see fig. 3 of B15. Most of the higher-mass galaxies were originally observed with a velocity resolution of Δ|$v$| = 5.5 km s⁻¹ through the ALFALFA survey (Haynes et al. 2011) and subsequently re-reduced by B15. Spectra collected through the low-mass H i programme were Hanning smoothed to a resolution of 5 km s⁻¹ by B15; those from the ALFALFA survey, to a resolution of 10 km s⁻¹. In this work, we only analyse the subsample of the B15 data with a smoothed signal-to-noise ratio (SNR) > 15 (see equation 2); this subsample contains 2002 galaxies with a median SNR of 25.

For our main analysis, we mock-observe all simulated galaxies at z = 0 from 100 random viewing angles distributed uniformly on the unit sphere. For a given viewing angle and desired SNR, we apply the following procedure to create a mock-H i profile:

• We consider all H i gas within a circular aperture of diameter 80 kpc centred on the galaxy. This is comparable to the 3.5 arcmin Arecibo beam diameter for an object at a distance of 70 Mpc, which is typical for our sample of observations. Our results are not sensitive to this choice, because the aperture is larger than the region containing significant H i for all the simulated galaxies in our sample. We calculate the H i mass of each gas element in the simulation as in E18, excluding ionized gas and cold gas likely to be molecular.² The neutral hydrogen fraction of each gas element is computed by GIZMO, which incorporates ionization from a spatially uniform, redshift-dependent UV background (Faucher-Giguère et al. 2009) and approximates the radiation field due to local sources (see Hopkins et al. 2017).

• We construct the H i mass-weighted LOSVD of all H i gas in this aperture, binned at a velocity resolution Δ|$v$|⁠. To account for thermal broadening, we treat each gas element as a Gaussian in velocity space with |$\sigma = c_s=\sqrt{kT/(\mu m_p)}$|⁠, where c_s is the isothermal sound speed. We choose Δ|$v$| to match the observations against which we compare: for galaxies with M_star < 10^8.5 M_⊙, we use |$\Delta v = 0.65\, \rm km\,s^{-1}$|⁠; for more massive galaxies, we take |$\Delta v = 5.5\, \rm km\,s^{-1}$|⁠.

• For mock-observations of a target at a distance D, we compute the 21-cm flux density in the ith velocity channel using the formula  

  \begin{eqnarray} \left(\frac{S_{{21,i}}}{{\rm Jy}}\right)\!=\!\frac{1}{2.36\times 10^{5}}\left(\frac{m_{{\rm {\rm H}\,\small {I}},i}}{M_{{\odot }}}\right)\!\left(\frac{\Delta v}{{\rm km\,s^{-1}}}\right)^{-1}\left(\frac{D}{{\rm Mpc}}\right)^{-2}, \end{eqnarray}

  (1)

  implicitly assuming the 21-cm transition is optically thin (e.g. Haynes & Giovanelli 1984). Here, |$m_{{\rm {\rm H}\,\small {I}},i}$| is the H i mass in the ith velocity channel. In practice, we solve for D for each observation in order to yield the desired SNR, as quantified in equation (2).

• We add Gaussian noise to the spectrum with a fixed dispersion of 0.004 Jy, which is typical for the observations against which we compare. We then Hanning smooth the spectrum to a resolution consistent with the observations: for galaxies with M_star < 10^8.5 M_⊙, we smooth to a resolution of |$5\, \rm km\,s^{-1}$|⁠; for higher-mass galaxies, to a resolution of |$10\, \rm km\,s^{-1}$|⁠.

• We identify the H i emission region of the spectrum in velocity space as follows. We construct a copy of the noisy spectrum that is Hanning smoothed to twice the resolution adopted in (iv). We define the edges of the emission region as the first point on either side of the emission peak where this doubly-smoothed spectrum drops below zero flux; we set the emission to zero outside this region.

   Denoting the peak flux density f_peak, we fit a quadratic polynomial to the wings on either side of the line profile between 0.5f_peak and the edges of the emission region, replacing the flux in this region with the polynomial interpolation. This procedure was found by B15 to minimize spurious signal due to noise in the wings of the profile; we adopt it here for consistency.

• We calculate the total integrated flux of the H i line as |$S_{21,{\rm tot}}=\int _{-\infty }^{\infty }S_{{\rm 21}}\left(v\right)\,{\rm d}v$|⁠, and W₅₀, the velocity width at 50 per cent of the maximum flux, as the difference between the velocities above and below the peak flux at which the interpolated emission profile reaches 0.5f_peak. We then calculate the SNR of the integrated emission using the same definition used by B15:  

  \begin{equation} {\rm SNR}=\left(\frac{S_{21,{\rm tot}}}{W_{50}}\right)\frac{w_{{\rm smo}}^{1/2}}{\sigma _{{\rm rms}}}. \end{equation}

  (2)

  Here, |$w_{\rm smo}=W_{50}/(10\rm \,km\,s^{-1})$| and σ_rms is the standard deviation of the spectrum outside the emission region, after the smoothing in (iv) is applied. Note that this definition represents the total SNR of the H i line integrated over all channels; the SNR in any one channel will be lower than the value yielded by equation (2).

2.1 Quantifying unresolved line shapes

We use two distinct but related statistics to quantify the shape of galaxies’ integrated H i lines.

2.1.1 Kurtosis

The kurtosis, which we denote κ, is a measure of the ‘peakedness’ of the line; i.e. how centrally concentrated the emission is in velocity space. Single-peaked, Gaussian-like H i lines have relatively high kurtosis, while line profiles with more weight in the edges, such as the characteristic double-horned profile, have lower kurtosis. The kurtosis of galaxies’ integrated H i lines was previously used as a proxy for V/σ by Papastergis et al. (2016), who used a cut of κ ≤ −1.2 to select rotationally supported galaxies in a study of the local baryonic Tully–Fisher relation (BTFR). For an observed line profile S₂₁(⁠|$v$|⁠), we define the kurtosis as  

\begin{equation} \kappa =\frac{\mu _{4}}{\mu _{2}^{2}}-3. \end{equation}

(3)

Here, μ_n is the nth central moment:  

\begin{equation} \mu _{n}=\frac{\int _{-\infty }^{\infty }S_{{\rm 21}}\left(v\right)\left(v-\overline{v}\right)^{n}\,{\rm d}v}{\int _{-\infty }^{\infty }S_{{\rm 21}}\left(v\right)\,{\rm d}v}, \end{equation}

(4)

where |$\overline{v}$| is the flux-weighted mean velocity:  

\begin{equation} \overline{v}=\frac{\int _{-\infty }^{\infty }v S_{{\rm 21}}\left(v\right)\,{\rm d}v}{\int _{-\infty }^{\infty }S_{{\rm 21}}\left(v\right)\,{\rm d}v}. \end{equation}

(5)

The constant −3 in equation (3) sets the kurtosis scale so that a Gaussian has κ = 0. We find that for both our mock observations and real galaxies, the majority of H i lines have −1.5 ≲ κ ≲ 0.

2.1.2 Steepness of profile wings

To quantify the steepness of a profile's wings, we first calculate W₂₀ and W₅₀, the width of the line profile at 20 per cent and 50 per cent of its peak value. We then measure their scaled difference, ΔW ≡ (W₂₀ − W₅₀)/W₅₀. This quantity is a measure of the slope of the line profile wings: for profiles with steep wings, W₂₀ and W₅₀ are very similar, and ΔW is small. For profiles with sloping wings, W₂₀ is significantly larger than W₅₀, so ΔW is larger; for a Gaussian, ΔW = 0.52. Most of the line profiles in both our observed and simulated samples have 0 ≲ ΔW ≲ 1.

2.2 Relation between line shape and kinematics

In Fig. 1(a), we show how a galaxy's gas kinematics change the shape of its unresolved 21-cm profile. We show column density and velocity maps of three simulated galaxies spanning the range of gas kinematics and morphology in our sample, ranging from irregular and dispersion supported (m11a, with M_star = 10⁸ M_⊙) to rotationally supported but with significant disordered motion (m11h, with M_star = 4 × 10⁹ M_⊙), to rotationally supported and morphologically discy (m11b, with M_star = 10⁸ M_⊙). Details regarding the construction of these maps, as well as comparable maps for all the galaxies in our sample, can be found in E18. Fig. 1(a) shows that the simulated galaxies’ unresolved H i profiles become single peaked, with increasingly sloping wings and higher kurtosis, as they become more dispersion supported.

Figure 1.

Left: The effect of internal gas kinematics on the shape of galaxies’ unresolved 21-cm lines. The first (second) row shows the resolved zeroth (first) H i moment map for three low-mass galaxies in our sample, with increasing rotational support from left to right. Each map is produced with inclination i = 60 deg. The third row shows the global H i mass-weighted LOSVD. The bottom row shows the mock 21-cm spectrum (grey) and identified emission region (black; see Section 2). Unresolved spectra encode information about the galaxies’ gas kinematics: rotationally supported galaxies produce double-horned line profiles with steep wings, while dispersion-supported galaxies produce single-peaked profiles with shallower wings. Right: The bottom panels show unresolved H i data for two low-mass galaxies (M_star ∼ 10⁸M_⊙), one with a Gaussian-like profile and one with a double-horned profile. The top panels show corresponding SDSS images.

In Fig. 1(b), we show the unresolved H i profiles of two galaxies in our observed sample, both with M_star ∼ 10⁸ M_⊙. One of these galaxies exhibits a double-horned profile with low kurtosis (qualitatively similar to m11b), while the other produces a single-peaked profile with higher kurtosis (qualitatively similar to m11a). Spatially resolved H i kinematics are not available for our observed sample, so we instead show the galaxies’ optical images from the NASA-Sloan Atlas, a reprocessing of SDSS DR8 (Blanton et al. 2011). The optical morphology of the target with the double-peaked 21-cm profile is clearly discy, while that of the target with the Gaussian-like profile appears clumpy and irregular. Stellar morphology from optical images is an imperfect proxy for H i morphology, and we do not carry out any quantitative analysis on optical images. However, we find that on average, galaxies with double-peaked H i profiles exhibit discier optical morphologies upon visual inspection than those with single-peaked profiles.

To determine quantitatively how the shape of the unresolved H i line is related to a galaxy's resolved gas kinematics, we mock-observe all our simulations from 100 random viewing angles and calculate the kurtosis, κ, wing steepness, ΔW, and the degree of rotational versus dispersion support, V/σ, for each viewing angle. To calculate V/σ for a given viewing angle, we first construct resolved mock H i velocity maps with spatial resolution 0.1 kpc (see E18 for details) and then define  

\begin{equation} V/\sigma =\frac{\sum _{j}m_{j}\left|V_{{\rm los},j}\right|}{\sum _{j}m_{j}\sigma _{{\rm los},j}}. \end{equation}

(6)

Here, m_j is the H i mass in a given spatial pixel, V_{los, j} and σ_{los, j} are the H i mass-weighted mean line-of-sight velocity and velocity dispersion in a pixel, and the sum is over all spatial pixels in the data cube. We take this definition from Obreja et al. (2016). Although V/σ as calculated in equation (6) is dependent on the size of spatial pixels; i.e. the beam size in resolved observations, it has the advantage of being defined in terms of quantities measured directly from the simulation, without requiring any kinematic model fitting. Like κ and ΔW, V/σ is a line-of-sight dependent quantity. Projection effects can cause rotational velocities to contribute significantly to σ, so V/σ values calculated from equation (6) are typically lower than values calculated from cylindrically averaged 3D quantities (see E18).

In Fig. 2, we plot κ and ΔW for our simulated galaxies as a function of V/σ. For all quantities, we show the median and middle 68 per cent values across 100 random viewing angles. κ and ΔW are measured from mock unresolved H i profiles with an SNR of 25, the median value for our observed sample. In agreement with the qualitative trend apparent in Fig. 1(a), we find that κ and ΔW decrease nearly monotonically with increasing V/σ, indicating that these quantities can be used as a proxy for V/σ. We find that κ and ΔW are also both anticorrelated with the intrinsic, three-dimensional V_rot/σ, as defined in E18. However, the correlation is somewhat weaker than the correlation with V/σ as defined in equation (6) because the unambiguously rotationally supported galaxies all have κ ∼ −1.4 and ΔW ∼ 0.1 but have a wider range of V_rot/σ, with 6 ≲ V_rot/σ ≲ 12.

Figure 2.

Kurtosis (top) and steepness of the wings (bottom) of unresolved 21-cm line profiles of simulated galaxies versus V/σ (from mock-IFU measurements; equation 6). Error bars show 68 per cent scatter across 100 different viewing angles. Discy, rotationally supported galaxies with higher V/σ have systematically lower kurtosis (κ ≲ −1) and steeper profile wings ((W₂₀ − W₅₀)/W₅₀ ≲ 0.2) than dispersion supported galaxies. Rotationally supported galaxy and dispersion-supported galaxy can thus be separated statistically based only on their unresolved line profile shape.

Fig. 2 also shows that, for the simulated galaxies, V/σ, κ, and ΔW are all strong functions of mass, with higher-mass galaxies having higher V/σ and lower κ and ΔW. However, the trend of decreasing kurtosis and increasingly steep profile wings at higher V/σ is not simply a result of the joint dependence on mass: the lowest-mass galaxy in our sample that is strongly rotationally supported (m11b, with M_star = 10⁸M_⊙) still has low κ and ΔW. For more detailed discussion of trends in gas kinematics with mass for this sample of simulated galaxies, we refer to E18.

Although the non-negligible 68 per cent scatter in Fig. 2 shows that κ and ΔW do vary across viewing angles, rotationally supported system and dispersion-supported system are still separated on average. In Fig. 3, we show explicitly how κ varies with SNR (left) and inclination (right). For each inclination and SNR value, we compute the kurtosis of 10⁴ mock H i profiles for m11a and m11b, two of the galaxies shown in Fig. 1(a). Both of these galaxies have M_star ≈ 10⁸ M_⊙ and M_200m ≈ 10^10.7 M_⊙, but m11b is rotationally supported while m11a is dispersion supported. The left-hand panel of Fig. 3 shows that, at SNR > 15 (the minimum SNR of our observed sample), kurtosis measurements for a given system are robust: the typical κ values for the two galaxies remain well separated, and the spread of κ values across different SNR is less than that across different viewing angles at fixed SNR.

Figure 3.

Left: Effect of SNR on 21-cm line profile kurtosis measurements for a simulated discy, rotationally supported galaxy (blue) and an irregular, dispersion supported galaxy (red). Dashed line shows |$\rm SNR=25$|⁠, the median value for our observed sample. At each SNR, we mock-observe each galaxy from 100 random viewing angles and generate 100 different noise realizations for each viewing angle. We plot the median and middle 68 per cent scatter for the resulting 10⁴ profiles generated at each SNR. Right: Here, we vary inclination while fixing |$\rm SNR=25$|⁠, the median for the observed sample in Figs 4 and 5. κ is robust to changes in SNR down to SNR ∼ 15. κ does vary with inclination for discy systems, but the effects are large only for rare, nearly face-on viewing angles.

The right-hand panel of Fig. 3 shows how kurtosis varies with inclination angle, which is defined in terms of the net angular momentum vector of H i. For the discy, rotationally supported galaxy, kurtosis varies significantly with inclination: it has κ ∼ −1.2 when viewed edge on, but κ ∼ 0 when viewed face on. This is unsurprising, since the face-on disc has no rotationally projected along the line of sight; generically, rotationally supported disc galaxies will produce single-peaked, Gaussian-like profiles if viewed completely face on.

Fortunately, most rotationally supported galaxies will nevertheless produce 21-cm profiles with low kurtosis, because face-on viewing angles are intrinsically rare: for random disc orientations, the probability distribution of inclinations is p(i) di = sin i di, meaning that far more galaxies are observed nearly edge on than nearly face on. Because inclination effects on kurtosis measurements only become significant at i ≲ 40 deg (representing 23 per cent of viewing angles), rotationally supported galaxies will produce unresolved H i profiles with low kurtosis for most viewing angles.

For the dispersion-supported galaxy, there is no significant trend in κ with inclination, though there is still significant scatter across different viewing angles at fixed inclination. This is unsurprising, since dispersion-supported systems have no unambiguously preferred kinematic axis. Most lines of sight for such systems produce single-peaked line profiles roughly resembling Gaussians, but large-scale inflows and outflows of H i, as well as irregular H i distributions, can produce a variety of irregular and asymmetric profiles.

Although it is not shown in Fig. 3, we find that ΔW has a similar scaling with SNR and inclination to κ. Due to viewing angle effects, the mapping from line shape as quantified by κ or ΔW to V/σ is not unique. However, since real galaxies should have random orientations, the observed distribution of κ and ΔW values places a statistical constraint on the observed population's distribution of V/σ values, even though it is not always possible to unambiguously measure the kinematics of a particular galaxy from its unresolved line profile.

2.2.1 Effect of rotation curve shape

The detailed shape of a galaxy's unresolved 21-cm line profile depends not only on rotational versus dispersion support but also on the shape of its rotation curve. We investigate this dependence in detail in Appendix A. Our primary conclusion is that while a double-peaked 21-cm profile is almost always indicative of rotational support, a single-peaked, Gaussian-like profile can arise either due to dispersion support or due to a linearly rising rotation curve: if a galaxy's circular velocity curve rises sufficiently slowly and does not begin to flatten in the region occupied by the galaxy, a single-peaked profile can be produced even in the limit where H i is completely rotationally supported and follows the circular velocity curve, |$v_c(r)=\sqrt{r\times {\rm d\Phi }/{\rm d}r}$|⁠. For the simulated galaxies, we can measure the H i rotation curve and |$v$|_c(r) directly, and we find that Gaussian-like profiles are due to dispersion support, not linearly rising rotation curves. However, for observed galaxies with Gaussian-like H i profiles, it is generally not possible to distinguish between a dispersion-supported system and a rotationally supported system with a linearly rising rotation curve.

3 COMPARISON TO OBSERVATIONS

3.1 Kurtosis

We now compare the unresolved H i profiles generated from our simulations to those of real galaxies. Our observed sample consists of 2002 21-cm H i profiles with SNR > 15. Almost all of these profiles were presented in B15; we also include in our sample observations of few dozen low-mass galaxies that were obtained since the publication of that work and were reduced following the same procedure described in B15. The majority of this observational sample is completely morphologically blind. 210 of the profiles first presented by B15 were selected with a weak selection criterion against galaxies with inclinations near 0. However, we find the full observational sample has a distribution of inclinations (computed from optical axis ratios by B15) consistent with random orientations.

We compute the kurtosis of the smoothed emission regions of these profiles, as defined in Section 2, using equation (3). Stellar masses were taken from the NASA-Sloan atlas (Blanton et al. 2011). The median SNR of these data is 25, so we adopt SNR = 25 in generating mock H i profiles for the simulated galaxies. In this observational sample, there is no significant correlation between SNR and κ or M_star.

In Fig. 4, we show the distribution of kurtosis values for these galaxies as a function of mass, overplotting measurements for our simulated galaxies. At M_star ≳ 10¹⁰ M_⊙, most observed galaxies have κ ≲ −1, indicating that they have rotationally supported kinematics. This is also true for the simulated galaxies, which fall on top of the median observed relation. At this mass scale, the scatter in the kurtosis measurements of different profiles for individual galaxies, which is due primarily to different viewing angles, is almost as large as all the scatter in the observed sample. This suggests that gas-rich galaxies at M_star ≳ 10¹⁰ M_⊙ likely all have similar rotationally supported H i kinematics.

Figure 4.

Top: 21-cm line kurtosis measurements for simulated galaxies compared to ∼2000 observed galaxies from B15. Error bars show 68 per cent scatter across 100 random viewing angles; dashed lines show median and 68 per cent scatter for observed galaxies in each mass bin. At M_star > 10¹⁰ M_⊙, where all the simulated galaxies have rotationally supported H i discs, the simulated galaxies produce line profiles with kurtosis values very similar to the observed galaxy population. At M_star = 10⁸⁻¹⁰ M_⊙, most of the simulated profiles have higher κ than the observed population, corresponding to lower V/σ (see Fig. 2). Bottom: Gas fractions. Most low-mass-simulated galaxies with higher κ than the median of the observed sample also have somewhat lower gas fractions than typical observed galaxies.

At lower masses, particularly M_star = 10⁸⁻¹⁰ M_⊙, most of the simulated galaxies have significantly higher kurtosis than at Milky Way (MW) masses, with κ ∼ −0.5. A similar trend towards higher κ at lower masses is evident in the observations, but the trend is less pronounced than for the simulated galaxies: more than half of the observed galaxies have κ < −1 over this mass range, indicating rotationally supported kinematics. Of the eight simulated galaxies with M_star = 10⁸⁻¹⁰ M_⊙, seven have median κ values above the median observed κ–M_star relation. The only exception is m11b, which has a rotationally supported disc (Fig. 1a) and falls below the median observed relation.

3.2 Gas fractions

In the bottom panel of Fig. 4, we compare the gas fractions of the simulated galaxies to observations. For both the observed and simulated galaxies, we define |$M_{{\rm gas}} \equiv 1.4M_{{\rm {\rm H}\,\small {I}}}$|⁠, where the factor of 1.4 accounts for helium, and then plot f_gas = M_gas/(M_gas + M_star). For both the observed and simulated samples, f_gas decreases with increasing M_star. The simulations with M_star > 10¹⁰ M_⊙ have gas fractions similar to the median of the observed population, but at lower masses, most of the simulated galaxies have somewhat lower f_gas than the median observed relation.

In particular, most of the low-mass simulated galaxies with higher-than-average κ (i.e. enhanced dispersion support) also have lower gas fractions than the median of the observed population. The lowest-mass simulated galaxy with a clear rotationally supported disc (m11b, with M_star = 10⁸ M_⊙) also has a higher-than-average gas fraction ( f_gas = 0.94). This suggests that for the simulated galaxies, increased dispersion supported is likely connected to reduced gas content. Intriguingly, we find no correlation between f_gas and κ at fixed mass in the observed sample: at all masses, observed galaxies with lower-than-average gas fractions appear to have similar gas kinematics to those with higher-than-average gas fractions. We return to this issue in Section 4.1.

3.3 Line widths and BTFR

The width of a galaxy's unresolved 21-cm line, as quantified by W₂₀ or W₅₀, is typically used as a diagnostic of its maximum circular velocity, and thus, the depth of its gravitational potential. Velocity widths are used to construct the BTFR for gas-rich galaxies (e.g. McGaugh et al. 2000). In this context, individual line widths are to first order agnostic to rotational versus dispersion support, since the range of line-of-sight velocities in a virialized system depends primarily on the depth of the gravitational potential. Different BTFR slopes can be constructed depending on whether W₂₀ or W₅₀ is used; see Bradford et al. (2016) for a review.

To compare the BTFR for our simulated galaxies with observations, we follow B15 in defining baryonic masses for the simulated galaxies as |$M_{\rm bar} = M_{\rm star}+1.4M_{\rm {\rm H}\,\small {I}}$|⁠. We correct the observed W₂₀ and W₅₀ values for redshift broadening by dividing them by (1 + z); this correction is always small because z ≪ 1. Most previous works on the BTFR have also attempted to correct for inclination effects by dividing observed line widths by sin i, where inclination is calculated from galaxies’ projected optical axis ratios with the assumption of an intrinsic thin-disc morphology. Such corrections can be problematic for low-mass galaxies, both because optical axis ratios are an imperfect tracer of H i morphology and because there is no guarantee that galaxies’ H i morphologies are intrinsically discy. We therefore do not apply any inclination correction to the observed or simulated line widths, instead accounting for the random orientations of the observed galaxies by mock-observing the simulated galaxies along many random lines of sight. As a result, the intrinsic scatter in the observed line width–mass relation is larger than that found when an inclination correction is applied (see B15).

In Fig. 5, we compare the BTFR predicted by our simulations as quantified by W₂₀ (top) and W₅₀ (middle) to that of the observational sample from B15. At M_bar > 10^10.5 M_⊙, the simulated H i profiles have both W₂₀ and W₅₀ in good agreement with the observed sample. At lower masses, the simulated W₂₀ values remain in reasonable agreement with the observed values, while the simulated W₅₀ values are on average lower than the median of the observed values. Because our simulated galaxies – and the observed galaxies – to a lesser extent preferentially produce 21-cm profiles with sloping wings at lower masses, the slope of the BTFR measured by W₅₀ is steeper than that measured by W₂₀.

Figure 5.

Spatially unresolved BTFR as probed by W₂₀ (top) and W₅₀ (middle), without inclination corrections. Error bars show 68 per cent scatter across 100 viewing angles. Bottom panel shows ΔW = (W₂₀ − W₅₀)/W₅₀, an estimate of the slope of the wings of the H i profile. The simulations reproduce the observed BTFR as probed by W₂₀ well. At ∼MW masses (M_bar > 10^10.5 M_⊙), the simulated galaxies also have similar W₅₀ and ΔW to the observed galaxies. At lower masses, most simulated galaxies have W₂₀ similar to the observed population but lower W₅₀ and thus higher ΔW, indicating less-steep profile wings and higher dispersion than average for the observed galaxies.

The bottom panel of Fig. 5 shows ΔW. At M_bar ≳ 10^10.5 M_⊙, values for the simulated and observed galaxies agree well. At M_bar = 10^8.5−10 M_⊙, most of the simulated galaxies have higher ΔW, corresponding to shallower H i profile wings, than typical observed galaxies. Consistent with Fig. 4, the higher ΔW for low-mass simulated galaxies suggests that they have enhanced dispersion support compared to the observed galaxies.

Measured from unresolved line widths, the BTFR is relatively insensitive to rotational versus dispersion support. This is in contrast to measurements from spatially resolved gas kinematics (e.g. Simons et al. 2015; Bloom et al. 2017), which can distinguish between rotation and dispersion. For unresolved observations, simulated galaxies can be expected to match the observed BTFR as long as they fall on a realistic M_bar–M_halo relation and are approximately virialized, irrespective of the mode of kinematic support. The simulated galaxies reproduce the observed BTFR as measured by W₂₀ well, despite the apparent excess dispersion support in low-mass galaxies. On the other hand, the dispersion-supported simulated galaxies fall slightly below the observed BTFR as measured by W₅₀. This occurs because W₅₀ is significantly less than W₂₀ for a dispersion-supported system with a Gaussian LOSVD, but the two quantities are similar for a rotationally supported systems with steep wings. This suggests that the scatter in the BTFR can be reduced by measuring W₂₀ as opposed to W₅₀, or by inflating W₅₀ values for targets with large ΔW and κ. Previous studies with observational data (Bradford et al. 2016; Papastergis et al. 2016) have indeed found this to be the case. Incorporating measurements of line shape as quantified by κ and ΔW in the observed BTFR will likely prove useful in minimizing the relation's scatter for the purpose of using it as a redshift-independent distance estimator (e.g. Obreschkow & Meyer 2013).

4 SUMMARY AND DISCUSSION

We have investigated the connection between the shape of galaxies’ spatially unresolved 21-cm line profiles and the kinematics of their gas. Rotationally supported galaxies typically produce double-peaked line profiles with steep wings, while dispersion-supported galaxies produce single-peaked Gaussian-like profiles with more sloping wings (see Fig. 1a). Using a suite of simulations from the FIRE project, we showed that the kurtosis and wing steepness of galaxies’ unresolved H i lines can thus be used as a proxy for V/σ (Fig. 2). Although the mapping between line shape and V/σ also depends on viewing angle (Fig. 3) and on the shape of galaxies’ rotation curves (Appendix A), marginalization over these parameters with a large observational sample provides statistical constraints on the degree of rotational support in the local galaxy population, without requiring spatially resolved data.

Unresolved 21-cm line profiles provide a homogeneous, unbiased observational sample with which to test the predictions of simulations. Because mock observations of unresolved H i are straightforward to carry out with minimal assumptions and no model fitting, they allow for robust comparison of simulations and observational data on equal terms. We compared mock 21-cm line profiles of our simulated galaxies to those of a sample of ∼2000 low-redshift observed galaxies, quantifying their shape through the measured kurtosis (Fig. 4) and the difference between their measured 20 per cent and 50 per cent velocity widths (Fig. 5). The primary results of these comparisons are as follows:

• Kinematics in MW-mass galaxies: At M_star = 10¹⁰⁻¹¹ M_⊙, the simulated galaxies all have rotationally supported gas discs and produce H i profiles with shapes in excellent agreement with the observed population. At this mass scale, the scatter in observed kurtosis values is small enough to be explained almost entirely by inclination effects, indicating that most H i-rich galaxies with M_star = 10¹⁰⁻¹¹ M_⊙ have similar, rotationally supported gas kinematics.

• Kinematics in lower-mass galaxies: At M_star = 10⁸⁻¹⁰ M_⊙, most of the simulated galaxies have higher kurtosis values and less-steep wings than more massive galaxies. This is true for the observed sample as well. However, the simulated lower-mass galaxies appear to have somewhat more dispersion support than the observed systems (Figs 4 and 5). This is not universally true: one low-mass galaxy (m11b, with M_star = 10⁸ M_⊙) does have a kinematically cold, rotationally supported gas disc and produces a characteristic double-horned H i profile. However, seven of the eight simulated galaxies in our sample with M_star = 10⁸⁻¹⁰ M_⊙ produce H i profiles with higher kurtosis and shallower wings than the median of the observed sample.³

• Gas fractions: The simulated MW-mass galaxies have gas fractions similar to the observed galaxies. At lower masses (M_star = 10⁸⁻¹⁰ M_⊙), the gas fractions of most simulated galaxies are somewhat lower than the median of the observed population (Fig. 4). For the simulated low-mass galaxies, f_gas and V/σ are correlated at fixed mass: galaxies with lower-than-average gas fractions have higher-than-average κ and ΔW, indicating enhanced dispersion support. We do not find such a correlation for the observed galaxies.

• Baryonic Tully–Fisher relation: At all masses, the simulated galaxies have 20 per cent line widths in good agreement with the observed BTFR, irrespective of whether they are supported by rotation or dispersion (Fig. 5). Global line widths are primarily probes of galaxies’ virial masses and are agnostic to the mode of kinematic support.

At the lowest mass scales probed by our simulations (M_star ≃ 10^6−7.5 M_⊙), the observed sample does not contain enough observed galaxies for a robust comparison between the simulations and observations. Like the simulated galaxies, the lowest-mass observed galaxies preferentially have single-peaked profiles with high κ and ΔW, but this could be explained by either increased dispersion support or linearly rising circular velocity curves (see Appendix A).

4.1 Possible selection effects

We note that although no explicit morphological selection criteria were used in constructing our observed sample, some selection effects could be present. In particular, our observed galaxies are all detected in both H i and optical light, which could select against (a) low surface brightness or (b) H i-poor galaxies. Many of our low-mass galaxies do have low surface brightness (Chan et al. 2017). We do not, however, expect a selection effect against low surface brightness galaxies in our observed sample to significantly bias our observationally inferred kinematics, as observed low surface brightness galaxies typically have ordinary levels of rotational support in their gas kinematics (de Blok, McGaugh & Rubin 2001; Oh et al. 2015).

Because most of the dispersion-supported galaxies in our simulation sample also have somewhat lower gas fractions than typical observed galaxies with similar mass, one might expect selection effects against H i-poor galaxies to introduce biases against dispersion-supported systems in the observed sample. However, the lack of any correlation between gas fraction and kinematics in the observed sample suggests that the excess dispersion found in the simulated galaxies is not primarily due to poor observational completeness at low f_gas: if it were, observed galaxies with lower-than-average f_gas would show enhanced dispersion support.

We also consider the possibility of biases in the sample of simulated galaxies. The primary selection effect in our simulation sample is that we target haloes that are isolated at z = 0, meaning that they have no more massive neighbours within at least 3 × R_200m. We have verified that the sample spans a normal range of halo spin, concentration and formation time E18, and previous work has also found that targeting isolated haloes does not introduce significant biases in most halo properties (Oñorbe et al. 2014).

To test explicitly whether our selection of isolated haloes introduces biases in the galaxy properties studied here, we repeated our analysis with only the 923 galaxies in our observational sample (out of 2002 total) that were found to be isolated by B15. We found this isolated subsample to have similar f_gas, κ, and ΔW to the full observational sample.⁴ We therefore conclude that simulating haloes that are isolated at late times does not introduce significant biases for the galaxy properties studied here.

4.2 Excess dispersion in low-mass galaxies

The resolved gas properties of the simulated galaxies studied here were previously investigated and compared to observations by E18. In qualitative agreement with our results, that work found most galaxies with M_star = 10⁸⁻¹⁰ M_⊙ to be somewhat less rotationally supported than typical observed galaxies at the same mass scale (e.g. their fig. 15). Due to imperfectly understood selection effects in the observed samples to which they compared, E18 concluded that there were hints of enhanced dispersion support in the simulations at low masses, but they could not verify the disagreement to be robust. Our results also find evidence of excess dispersion support at low masses, but with comparison to a larger, morphologically unbiased observational sample. We emphasize that although our simulated low-mass galaxies are on average more dispersion-supported than the observed population, they still fall in regions of parameter space populated with real galaxies; i.e. the simulations do not produce galaxies in observationally ‘forbidden’ regions of parameter space.

The mass scale at which we find the most enhanced gas dispersion in the simulated galaxies compared to observations is precisely the scale at which feedback-driven outflows have been shown to most strongly alter galaxies’ gravitational potentials (Di Cintio et al. 2014; Chan et al. 2015; El-Badry et al. 2016; Bullock & Boylan-Kolchin 2017). In the FIRE simulations, dynamical heating through a feedback-driven, time-variable potential at this mass scale has been shown to produce dark matter (DM) cores (Chan et al. 2015), low central stellar surface densities (Chan et al. 2017), and dispersion-supported stellar dynamics with primarily radial orbits (El-Badry et al. 2017).

The same large-scale outflows that lead to core formation (a) continually inject energy into low-mass galaxies’ ISM, making it difficult for gas to settle into a rotationally supported disc E18 and (b) produce very bursty star formation histories (SFHs) at low masses. This suggests that the excess gas dispersion we find at this mass scale is likely connected to tensions between the simulations and observations identified by previous work. Sparre et al. (2017) showed that FIRE galaxies with M_star ≲ 10^9.5 M_⊙ likely have SFHs that at z = 0 are on average burstier than those of typical low-mass field galaxies [though observations also imply non-trivial SFH burstiness at this mass scale; see Weisz et al. 2012], and Chan et al. (2017) showed that galaxies from the FIRE simulations with M_star = 10⁷⁻⁹ M_⊙ have, on average, lower surface brightnesses than typical observed galaxies at the same stellar mass.

One might therefore expect that low-mass galaxies with rotationally supported gas kinematics would have calmer SFHs and smaller cores than similar-mass galaxies with dispersion-supported gas; we investigate this further in Appendix B. We find that the lowest-mass galaxy in our simulation sample with a clear rotationally supported disc has an SFH that is substantially less bursty at late times than the SFHs of the dispersion-supported simulated galaxies at the same mass scale. It nevertheless has a DM core at z = 0, likely due to burstier star formation at early times, but the core is smaller than in dispersion-supported systems at the same mass scale.

The apparent excess of dispersion support in the low-mass simulated galaxies can serve as a useful constraint on feedback models and may motivate modifications to the FIRE model or the incorporation of additional physical processes in the future. We first consider the possibility that the FIRE model produces feedback that is in some sense ‘too strong’ at low masses. Because feedback energetics are taken directly from stellar evolution calculations, there is little room to arbitrarily reduce or increase the overall feedback energy budget.⁵ It is possible, however, that changing details in the implementation of physical processes that are currently modelled approximately (e.g. adding on-the-fly radiative transport or more explicit treatment of H ii regions) could change the burstiness of star formation. For example, if feedback can create hot ‘chimneys’ surrounding star-forming regions such that localized outflows can vent energy out of the disc (e.g. Ceverino & Klypin 2009; Hopkins, Quataert & Murray 2012), the effects of feedback on the bulk of a galaxy's gas may be less catastrophic than when energy and momentum are primarily deposited in cold gas. Alternatively, if radiative feedback can suppress star formation locally on time-scales shorter than the free-fall time in bound gas clouds (e.g. Faucher-Giguère 2018), star formation can become less spatially and temporally clustered, leading to weaker outflows and less dispersion-supported gas.

Simply turning off some feedback sources (e.g. radiative feedback or stellar winds) does not seem to reduce SFH burstiness, as it allows more stars to form before collapse is regulated by supernovae (Hopkins et al. 2017). As long as star formation is restricted to self-gravitating gas clouds, varying the local star formation efficiency or the critical density above which star formation is allowed to occur also does not significantly change galaxy SFHs in the default FIRE model, since star formation proceeds until feedback produces enough energy to regulate it (Hopkins et al. 2017). It may thus if fact be easier to facilitate disc formation by incorporating additional feedback processes, albeit processes that act less violently than supernovae, than by reducing the types of feedback sources or their energy budget. For example, Chen, Bryan & Salem (2016) found that including cosmic ray diffusion in galaxy formation simulations caused their low-mass galaxies to transition from being dispersion supported, with bursty SFHs, to rotationally supported, with smoother SFHs. Further investigation of the effects of changes in the details of feedback implementations on SFH burstiness and galaxy morphology represents a promising avenue for future work. We note, for example, that feedback resulting in weaker large-scale outflows in low-mass galaxies could also alleviate the difficulty identified by Anglés-Alcázar et al. (2017) of growing massive black holes in the FIRE simulations at early times.

The fact that our simulations do produce some low-mass galaxies with rotationally supported discs also suggests that at fixed mass, morphology is driven by a galaxy's formation history and/or environment. This warrants further investigation of the assembly histories of our low-mass galaxies and the primary factors that drive them to become rotationally or dispersion supported, as well as simulation of low-mass galaxies in a wider range of environments.

Acknowledgements

We thank the referee, Danail Obreschkow, for a constructive report and Alyson Brooks and Pascal Elahi for helpful discussions. K.E. acknowledges support from a Berkeley graduate fellowship, a Hellman award for graduate study, and an NSF Graduate Research Fellowship. E.Q. and K.E. are supported by a Simons Investigator Award from the Simons Foundation and by an NSF grant AST-1715070. M.B.K. acknowledges support from NSF grant AST-1517226 and from NASA grants NNX17AG29G and HST-AR-13888, HST-AR-13896, HST-AR-14282, HST-AR-14554, HST-GO-12914, and HST-GO-14191 from STScI. D.R.W. is supported by a fellowship from the Alfred P. Sloan Foundation. A.W. was supported by NASA through grants HST-GO-14734 and HST-AR-15057 from STScI. Support for P.F.H. was provided by an Alfred P. Sloan Research Fellowship, NASA ATP grant NNX14AH35G, NSF Collaborative Research grant #1411920, and CAREER grant #1455342. D.K. and T.K.C. were supported by NSF grants AST-1412153 and AST-1715101 and by the Cottrell Scholar Award from the Research Corporation for Science Advancement. CAFG was supported by NSF through grants AST-1412836, AST-1517491, AST-1715216; by CAREER award AST-1652522; and by NASA through grant NNX15AB22G. We ran numerical calculations on the Caltech compute cluster ‘Wheeler’, allocations TG-AST130039 & TG-AST150080 granted by the Extreme Science and Engineering Discovery Environment (XSEDE) supported by the NSF, and the NASA HEC Program through the NAS Division at Ames Research Center and the NCCS at Goddard Space Flight Center. The analysis in this paper relied on the python packages numpy (Van Der Walt, Colbert & Varoquaux 2011), matplotlib (Hunter 2007), and astropy (Astropy Collaboration et al. 2013).

Footnotes

1

See the FIRE project website at http://fire.northwestern.edu. A public version of the GIZMO code is available at http://www.tapir.caltech.edu/∼phopkins/Site/GIZMO.html.

2

Gas with |$T<300\, \rm K$| and n_H > 10 cm⁻³ is treated as molecular. Our results are not sensitive to excluding molecular gas, as only a small fraction (∼10 per cent by mass in MW–mass galaxies, and less at lower masses) of the ISM passes these molecular criteria at z = 0. We note that these molecular fractions are likely somewhat lower than typical values for observed galaxies.

3

We also note that the m11b halo has a relatively high spin parameter, with λ = 0.077, which is 2σ above the mean for all dark matter haloes (Bullock et al. 2001); this may be related to its discier-than-average morphology and kinematics.

4

This is true only because the observed galaxies all retain enough cold gas to be detected in H i. Stripped galaxies in dense cluster environments likely do have different gas properties than isolated galaxies, but such targets do no enter our observational sample.

5

In Hopkins et al. (2017) and Hopkins et al. (in preparation), we have nevertheless carried out tests of the effects of modified feedback recipes on the burstiness of star formation in dwarfs. While reducing the feedback energy budget can make galaxies somewhat less bursty and morphologically discier at fixed halo mass, doing so increases galaxy stellar masses, such that at fixed stellar mass, burstiness is not significantly changed. Significantly changing the feedback energy budget also produces a M_star–M_halo relation in worse agreement with constraints from abundance matching than the default FIRE model.

6

We generally treat the rotation curve as arising from a spherically symmetric mass distribution, but we emphasize that the model is agnostic to the precise form of the rotation curve; for example, accounting for the potential due to the disc simply requires substituting the appropriate form of |$v$|_c(r) for |$v$|_{ϕ, 0}(r).

7

In practice, the integrand is undefined for |$v_{\phi }^2 < v_{\rm los}^2$|⁠, so we only integrate over r_min < r < r_max, where r_min and r_max bound the region where |$v$|_ϕ(r) > |$v$|_los. This makes sense physically, because annuli where |$v$|_ϕ(r) < |$v$|_los clearly will not have any material moving at velocity |$v$|_los.

8

These curves are presented for all the simulations analysed in this work in the appendix of E18.

9

Here, the SFH is defined ‘archaeologically’ for all stars within 20 kpc at z = 0, and the SFR at a given time represents the average over the previous 100 Myr. The SFH defined in this way differs slightly from that of the main progenitor, since a small fraction (<10 per cent) of the stars in the galaxies at z = 0 formed ex situ.

10

We multiply the density profile from the DM-only simulation by (1 − Ω_b/Ω_m) to correct for the lack of baryons.

APPENDIX A: TOY MODEL FOR UNRESOLVED H i

In this section, we use a simple analytic model to show how the shape of a galaxy's integrated H i profile depends on its rotation curve, gas velocity dispersion, and gas surface density profile. We begin by considering a thin rotating gas disc in the x–y plane, with H i surface density profile Σ(r) and an intrinsic rotation curve |$v$|_{ϕ, 0}(r), viewed at inclination angle i; for a rotation-supported, spherically symmetric system in equilibrium, |$v_{\phi ,0}(r) = v_c(r) = \sqrt{GM({<}r)/r}$|⁠.⁶ We first consider an infinitesimal parcel of gas at (r, ϕ), defining the coordinate system such that ϕ = 0 points towards the observer, who is at an angle |$\theta = {\pi} /2-i$| above the disc mid-plane. The parcel's projected radial velocity along a line of sight to the observer is 

\begin{eqnarray} v_{\rm los}(r, \phi )=v_{\phi ,0}(r)\sin \phi \sin i=v_{\phi }(r)\sin \phi , \end{eqnarray}

(A1)

where we define |$v$|_ϕ(r) = |$v$|_{ϕ, 0}(r)sin i. The mass in this gas parcel is dm = Σ(r) × dr × rdϕ, so the gas mass at radius r in an infinitesimal velocity bin d|$v$|_los follows  

\begin{eqnarray} \frac{{\rm d}m}{{\rm d}v_{\rm los}}(r, v_{\rm los})=2\frac{{\rm d}m}{{\rm d}\phi }\frac{{\rm d}\phi }{{\rm d}v_{\rm los}}=2r\Sigma \left(r\right){\rm d}r\frac{{\rm d}\phi }{{\rm d}v_{\rm los}}, \end{eqnarray}

(A2)

where the extra factor of 2 arises because material at ϕ has the same |$v$|_los as material at |${\pi} - \phi$|⁠. Differentiating equation (A1) with respect to ϕ and inverting the derivative, we find  

\begin{eqnarray} \frac{{\rm d}\phi }{{\rm d}v_{\rm los}}(r, v_{\rm los})=\frac{1}{v_{\phi }\left(r\right)\cos \phi }=\frac{1}{\sqrt{v_{\phi }^{2}\left(r\right)-v_{\rm los}^{2}}}. \end{eqnarray}

(A3)

The total H i mass in the disc per velocity interval d|$v$|_los can then be found by integrating equation (A2) over all⁷ radii:  

\begin{eqnarray} M(v_{\rm los})=\frac{{\rm d}M}{{\rm d}v_{\rm los}}=2\int _{0}^{\infty }\frac{\Sigma (r)r\,{\rm d}r}{\sqrt{v_{\phi }^{2}(r)-v_{\rm los}^{2}}}. \end{eqnarray}

(A4)

Thus far, we have assumed that the projected rotation velocity |$v$|_ϕ(r) is the same for all particles in a given annulus, as will be the case if the gas is supported purely by rotation. However, if there is intrinsic dispersion in the rotation velocity at a given radius, this assumption will not hold. If we take the azimuthal velocity distribution in a given annulus to be a Gaussian centred on |$v$|_ϕ(r) with dispersion σ(r), then equation (A3) for a single annulus must be convolved with a Gaussian. Again integrating over all radii, we obtain  

\begin{eqnarray} M(v_{{\rm los}})=\sqrt{\frac{2}{{\pi} }}\int _{0}^{\infty }\int _{-\infty }^{\infty }\frac{\exp \left[-\frac{\left(v_{{\rm los}}-\tilde{v}\right)^{2}}{2\sigma ^{2}(r)}\right]\Sigma (r)r}{\sigma (r)\sqrt{v_{\phi }^{2}(r)-\tilde{v}^{2}}}\,{\rm d}\tilde{v}\,{\rm d}r, \end{eqnarray}

(A5)

which obeys the expected normalization condition that |$M_{{\rm {\rm H}\,\small {I}}}=\int _{-\infty }^{\infty }M\left(v_{{\rm los}}\right)\,{\rm d}v_{{\rm los}}=2{\pi} \int _{0}^{r}\Sigma \left(r\right)r\,{\rm d}r$|⁠. Assuming negligible H i self-absorption, the observable 21-cm flux density is then  

\begin{eqnarray} S_{{\rm 21}}(v_{{\rm los}})=4.24\times 10^{-6}\left(\frac{M(v_{{\rm los}})}{M_{{\odot }}}\right)\left(\frac{D}{{\rm Mpc}}\right)^{-2}\,{\rm Jy}. \end{eqnarray}

(A6)

We note that a very similar model was used by Obreschkow et al. (2009a) to predict the unresolved H i and CO lines of galaxies for their semi-analytic model.

In Fig. A1, we show how the H i profile shapes predicted by equation (A5) depend on |$v$|_ϕ(r) (top), Σ(r) (middle), and σ(r) (bottom). In each row, the black line shows the same fiducial model: a Hernquist (1990) rotation curve with total mass M = 10¹¹ M_⊙ and scale length a = 20 kpc, an exponential H i surface density with scale length r_d = 2.5 kpc, and a constant gas dispersion of |$\sigma = 5\, \rm km\,s^{-1}$|⁠. We then show the effect of varying each of these parameters, considering constant and linearly rising rotation curves (top row), an exponential surface density with r_d = 5 kpc, a constant density disc (middle row), and gas dispersions of σ = 15 and 30 |$\rm km\,s^{-1}$| (bottom row).

Figure A1.

Schematic effect of the rotation curve, gas surface density, and gas velocity dispersion on a galaxy's observable H i profile as predicted by our toy model (Section A). The black line shows the same fiducial model in all three rows. We vary the rotation curve (top), surface density profile (middle), and dispersion (bottom), using equation (A6) to predict the H i profile for each configuration. Numbers in the upper right-hand panels indicate the profile's kurtosis. We assume |$M_{\rm {\rm H}\,\small {I}} = 10^9\,{\rm M}_{{\odot }}$| and a distance of 30 Mpc in all panels. Increasing the velocity dispersion or decreasing the rotation velocity in the region occupied by most of the gas causes the line profile to transition from its characteristic double-horned profile to a single-peaked profile with broader wings.

Fig. A1 shows that double-horned H i profiles arise naturally for rotating geometries when the gas velocity dispersion is low compared to the maximum rotation velocity and a significant fraction of the gas mass (M ∼ r²Σ(r)) is rotating at near the maximum rotation velocity. A single-peaked H i profile with high kurtosis and shallow wings can be produced either by a high velocity dispersion (dashed red line, bottom row) or by a low rotation velocity in the region where most of the gas is found (dashed red line, top row). Low-mass galaxies are thus not guaranteed to produce double-peaked H i profiles even if their gas is rotationally supported; haloes with solid-body rotation curves extending to the outskirts of the H i disc will generically produce single-peaked unresolved profiles.

For the simulated galaxies, we can measure |$v$|_c(r) and |$v$|_ϕ(r) directly.⁸ In most of the simulated galaxies with M_star ≲ 10¹⁰ M_⊙, the circular velocity |$v$|_c(r) is qualitatively similar to the default rotation curve adopted in Fig. A1, but the true gas rotation velocity |$v$|_ϕ(r) is significantly lower; for these systems, the lack of double-peaked profiles is verifiably due to enhanced dispersion support and non-circular motions, not a linearly-rising circular velocity curve.

For the simple model considered here, double-peaked line profiles can only be produced by rotation. The model will break down, however, when dispersion support begins to dominate or outflows drive large non-circular motions. In such cases, it is in principle possible to produce a double-horned profile purely due to radial kinematics, i.e. when a galaxy with a large bipolar outflow containing most of the H i mass is viewed head on. Such outflows are common in the simulated low-mass galaxies we study (Muratov et al. 2015; El-Badry et al. 2016, 2018), but the simulated galaxies that are not rotationally supported generally do not produce double-peaked lines. This, coupled with the ubiquity of double-peaked profiles in the observed sample, suggests that rotation is the primary driver of the double-peaked profiles and lower kurtosis values found in the observed sample.

Several resolved observational studies (Swaters et al. 2009; Oh et al. 2015; Iorio et al. 2017) have found some low-mass galaxies to have solid-body rotation curves out to ∼2 disc scale radii; as illustrated in the top row of Fig. A1, such systems can produce single-peaked unresolved profiles even if they have rotationally supported discs. We therefore caution that the increase in the kurtosis values of the observed sample at low masses could be attributed either to increased dispersion support or to linearly rising rotation curves.

In Fig. A2, we show the relation predicted by this analytic model between V_rot/σ, unresolved line kurtosis, and the steepness of the line's wings. Here, we assume a Hernquist rotation curve, exponential surface density profile, and constant gas dispersion, but our predictions are qualitatively similar for a wide range of models. We define V_rot as the maximum value of |$v$|_ϕ(r); we note that V_rot/σ defined in this way will generally be higher than the observationally motivated V/σ defined in equation (6) due to projection effects. In agreement with the trend found for our simulated galaxies (Fig. 2), the model predicts that κ and (W₂₀ − W₅₀)/W₅₀ increase monotonically as galaxies transition from rotational to dispersion support.

Figure A2.

Relation between H i line kurtosis, (W₂₀ − W₅₀)/W₅₀, and V_rot/σ as predicted by our toy model, assuming a Hernquist rotation curve and exponential disc surface density profile (see Section A). Kurtosis and (W₂₀ − W₅₀)/W₅₀ are both measures of the steepness of the unresolved H i profile's wings and are inversely correlated with V_rot/σ.

Although we have not done so in this work, we note that this model can also be applied to fit the observed H i profiles of real galaxies; i.e. given parametrized models for |$v$|_ϕ(r), σ(r), and Σ(r), it is a straightforward exercise in optimization or sampling to determine what set of parameters produce an unresolved profile in best agreement with an observed profile. In this case, the tight correlation between observed galaxies’ H i disc sizes and H i masses (e.g. Wang et al. 2015) places a prior on Σ(r), so that fitting an observed 21-cm profile strongly constrains |$v$|_ϕ(r) and σ(r). Previous works (e.g. Stewart et al. 2014; Westmeier et al. 2014) have presented other flexible fitting functions for unresolved 21-cm profiles that capture most of the observed diversity in line shapes. The advantage of using a model such as the one introduced here is that the parameters returned from fitting are directly interpretable in terms of physical properties of the galaxy. The realization of next-generation radio facilities such as WALLABY and the SKA will produce H i observations of more than half a million galaxies, the vast majority of which will be spatially unresolved or marginally resolved (Dewdney et al. 2009; Duffy et al. 2012). Interpreting these galaxies’ H i profiles with physically motivated models will place strong statistical constraints on the distribution and kinematics of gas throughout the galaxy population.

APPENDIX B: RELATION BETWEEN MORPHOLOGY, SFH, AND CORES

As discussed in Section 4, the dispersion-supported kinematics and irregular morphologies of many of our simulated low-mass galaxies are likely related to their bursty SFHs. To investigate this further, we compare in Fig. B1 the SFHs⁹ of m11a and m11b, the ‘irregular’ and ‘discy’ galaxies shown in Fig. 1(a). Both galaxies have bursty SFHs at early times, but the SFH of m11b becomes calmer at late times, while that of m11a remains quite bursty until z = 0. In m11a, strong bursts of star formation and subsequent temporary quenched episodes occur when a large fraction of the galaxy's gas falls into the galactic centre and accumulates at high densities. Because most of the gas in m11b is angular momentum supported at late times, it cannot coherently fall into the centre, and star formation is regulated more locally in the disc (see E18 for details).

Figure B1.

Top: SFHs of two of our simulated galaxies with similar stellar and halo masses but different morphologies. The discy galaxy has a substantially less bursty SFH at late times than the irregular galaxy. Bottom: DM density profiles at z = 0 of the two haloes from our standard simulations with baryons and feedback compared to profiles of the same haloes in DM only simulations. Both galaxies have a cored density profile, but the irregular galaxy with a burstier SFH forms a somewhat larger core.

In the bottom panels of Fig. B1, we show the z = 0 DM density profiles of the two haloes. To highlight the effects of feedback-driven potential fluctuations in reducing the central density of DM, we compare the DM density profiles from the main simulations analysed here to those from collisionless DM-only simulations of the same zoom-in region.¹⁰ Consistent with previous work on coring at this mass scale (e.g. Chan et al. 2015), we find that both galaxies have substantially reduced DM densities in the central ∼few kpc in the simulations with baryons, indicating that feedback-driven potential fluctuations have redistributed DM from the central regions to larger radii.

The core is somewhat larger in m11a, the dispersion-supported galaxy. This is not surprising, since the galaxy's burstier SFH at late times indicates that it has undergone larger outflows and rapid potential fluctuations at late times (see e.g. Pontzen & Governato 2014 for discussion of the conditions required for feedback-driven core formation). The presence of a core in m11b can likely be attributed to its bursty SFH at early times (e.g. Madau, Shen & Governato 2014; Oñorbe et al. 2015). Indeed, we find that the two galaxies have qualitatively similar morphologies and evolutionary histories until z ≈ 0.5, when m11b goes through a merger. This brings in a large quantity of high-angular momentum gas, and the galaxy quickly forms a rotation-supported disc thereafter.