KINEMATICS OF EXTREMELY METAL-POOR GALAXIES: EVIDENCE FOR STELLAR FEEDBACK

For the study of the kinematics of XMPs, we selected a sample of nine galaxies, eight of which were previously studied in Paper I (see Table 1), and thus have spatially resolved measurements of oxygen abundance. This time the galaxies were observed with the ISIS⁸ instrument on the 4.2 m William Herschel Telescope (WHT), using a high-resolution grating.⁹ The slit of the spectrograph was positioned along the major axis of the galaxies (see Figure 1), overlaying the slit used in Paper I. The slit width was selected to be 1 arcsec, according to the seeing at the telescope and also matching the value of our previous observation of these galaxies in Paper I. The spectrograph has two arms, the red one and the blue one, each with its own detector. All galaxies were observed in both arms, except for HS0822, where we have only the red spectrum. A 2 × 2 binning of the original images was performed to reduce oversampling. After binning, the dispersion turns out to be 0.52 Å pix⁻¹ for the red arm and 0.46 Å pix⁻¹ for the blue arm. Similarly, the spatial scale is 0.44 arcsec pix⁻¹ for the red arm and 0.40 arcsec pix⁻¹ for the blue arm. The wavelength coverage for the blue arm is 940 Å centered at Hβ (from 4390 Å to 5330 Å), whereas the read arm covers 1055 Å around Hα (from 6035 Å to 7090 Å). The total exposure time per target was approximately two hours (Table 1), divided into snapshots of 20 minutes as a compromise to simultaneously minimize cosmic-ray contamination and readout noise.

Zoom In Zoom Out Reset image size

Table 1.  Parameters of the Observed XMP Galaxies

Our ID	Galaxy Name	Obs. Date	Seeinga	WUib	Dc	Scalec	Exp. Time	$\mathrm{log}{M}_{\star }$ d
			(arcsec)	(km s−1)	(Mpc)	(pc arcsec−1)	(s)	(M${}_{\odot }$)
ML4	${\rm{J}}030331.26-010947.0$	2012 Aug 16	0.5	38.1 ± 2.7	125 ± 9	608	3600	8.33 ± 0.33
HS0822	${\rm{J}}082555.52+353231.9$	2013 Feb 02	1.4	37.8 ± 1.6	9.6 ± 0.7	47	7200	6.04 ± 0.03
SBS0	${\rm{J}}094416.59+541134.3$	2015 Mar 15	0.9	41.1 ± 3.2	23.1 ± 1.6	112	4800	7.05 ± 0.05
SBS1	${\rm{J}}113202.41+572245.2$	2015 Mar 15	1.2	40.4 ± 3.6	22.5 ± 1.6	109	7200	7.53 ± 0.48
ML16	${\rm{J}}114506.26+501802.4$	2015 Mar 15	1.0	40.1 ± 4.0	23.6 ± 1.7	114	7200	6.71 ± 0.07
UM461	${\rm{J}}115133.34-022221.9$	2015 Mar 15	1.2	40.6 ± 4.0	12.6 ± 0.9	61	6400	6.7 ± 1
N241e	${\rm{J}}144412.89+423743.6$	2015 Mar 15	1.0	40.3 ± 4.0	10.2 ± 0.7	49	6000	6.52 ± 0.07
SBS2104	${\rm{J}}210455.31-003522.3$	2012 Aug 16	0.5	39.0 ± 2.2	20.1 ± 1.7	97	7200	6.19 ± 0.05
ML32	${\rm{J}}230209.98+004939.0$	2012 Aug 16	0.5	42 ± 5	138 ± 10	666	7200	8.39 ± 0.33

Notes.

^aMean value of the seeing parameter during the observation of each object, measured by RoboDIMM@WHT (http://catserver.ing.iac.es/robodimm/). ^bMean FWHM of the telluric lines measured in each spectrum. ^cDistance and scale at the position of the source, from NED (NASA/IPAC Extragalactic Database). ^dTotal stellar mass of the galaxy from MPA-JHU (see Paper I). ^eTarget not included in Paper I.

Download table as:  ASCIITypeset image

The spectral resolution was measured from the FWHM of several telluric lines, using the calibrated data described in the next section. They yield a value of approximately 40 km s⁻¹ (Table 1), which renders a spectral resolution close to 7500 at Hα. This resolution is in good agreement with the value expected for the 1 arcsec slit, once the anamorphic magnification of ISIS is taken into account (Schweizer 1979).

The angular resolution during observations ranges from fair (1.3 arcsec) to good (0.5 arcsec), as measured by the seeing monitor at the telescope (Table 1). Additional information on the galaxies and the observations is included in Table 1.

The spectra were reduced using standard routines included in the IRAF¹⁰ package. The process comprises bias subtraction, flat-field correction using both dome and sky flat-fields, removal of cosmic rays (van Dokkum 2001), and wavelength calibration using spectral lamps. In order to evaluate the precision of the wavelength calibration, we measured the variation of the centroids of several sky emission lines along the spatial direction (if the calibration were error-free, the wavelength of the telluric line should not vary). For different galaxies and emission lines, the rms variation of the centroids ranges from 0.036 Å to 0.072 Å, which correspond to 1.6 km s⁻¹ and 3.3 km s⁻¹, respectively. It drops even further if only strong telluric lines are considered.

The sky subtraction, the geometrical distortion correction, and the combination of the spectra of each single source were carried out using custom-made Python routines. For the sky subtraction, we performed a linear fit, sampling the sky on both sides of the galaxy spectrum. This model sky emission was then subtracted at each wavelength. The effect of the geometrical distortion is that the spatial and wavelength direction are not exactly perpendicular on the detector. We determine the angle between the two directions as the angle between the continuum emission of the brightest point in the galaxy and the telluric lines. A linear fit allowed us to model the geometric transformation needed to make the spatial and wavelength axes perpendicular. The spectra from different exposures of the same galaxy were realigned before combining them. The left panel in Figure 2 illustrates the final result of the reduction procedure. It contains the spectrum around Hα of the galaxy N241. The signal-to-noise ratio (S/N) in continuum pixels has an average value of approximately 5 in the inner 10 arcsec around the brightest H ii region of the galaxy. The S/N in Hα ranges from 200 to 900 at the brightest point in the galaxy, decreasing toward the outskirts. Because our main goal is the kinematic analysis, we only performed a rough flux calibration by comparison with the calibrated Sloan Digital Sky Survey spectra of the galaxies (SDSS DR-12; Alam et al. 2015), rescaling the flux by accounting for the difference between our 1 arcsec slit and the 3 arcsec diameter of the SDSS fiber. Specifically, the flux integrated over the central 3 arcsec of the galaxies is scaled to be equal to 0.42 times the SDSS flux, where 0.42 is the ratio between the areas covered by the slit and the fiber. SDSS spectra of SBS1 and UM461 are not available. Therefore, we use for calibration the spectra of these galaxies analyzed in Paper I.

Zoom In Zoom Out Reset image size

Since the photometric calibration is only approximate, we have preferred not to correct the spectra for internal reddening in the H ii region. This correction is based on comparing the observed Balmer decrement with its theoretical value (e.g., Emerson 1996, p. 472), which in practice implies comparing fluxes of spectra taken with the red arm (Hα) and the blue arm (Hβ) of the spectrograph. Rather than carrying out an uncertain correction, we have preferred to neglect the small expected reddening (e.g., Sánchez Almeida et al. 2016). If one uses the SDSS spectra of our targets to estimate the Balmer decrement, the ensuing reddening correction implies increasing Hα by 50% or less. Neglecting this correction does not modify any of the conclusions drawn in the paper.

Mean velocities, velocity dispersions, and rotation curves of the emitting gas are estimated following Sánchez Almeida et al. (2013, Section 3), and we refer to this work for details. For the sake of completeness, however, the main assumptions and the resulting equations are provided in the following.

The velocity, U, was calculated from the global wavelength displacement of the Hα emission. It was derived from the centroid of the line, and also from the central wavelength of a Gaussian fitted to the emission line. The two values generally agree, with differences $\lt 1$ km s⁻¹ in the regions of interest. The error in velocity was estimated by assuming the rms of the continuum variations to be given by noise, and then propagating this noise into the estimated parameter (e.g., Martin 1971). In the case of the centroid, the expression is analytic. For the Gaussian fit, it is the standard error given by the covariance matrix of the fit (e.g., Press et al. 1986, chapter 14). The spatially unresolved velocity dispersion, W_U, was calculated from the width of the Gaussian model fitted to the line, but also directly from the line profile. We correct the observed FWHM, W_U0, for the instrumental spectral resolution, W_Ui (≃40 km s⁻¹; see Section 2.1 and Table 1), for thermal motions, W_Ut (≃25 km s⁻¹, corresponding to H atoms at $1.5\times {10}^{4}$ K, a temperature typical of the H ii regions in our galaxies; see Table 2), and for the natural width of Hα, W_Un ($\simeq 7$ km s⁻¹), so that

$$\begin{eqnarray}&&{W}_{U}^{2}={W}_{U0}^{2}-{W}_{{Ui}}^{2}-{W}_{{Ut}}^{2}-{W}_{{Un}}^{2}.\end{eqnarray} \tag{ 1 }$$

Table 2.  Physical Parameters for the Main Star-forming Region in the Galaxies

ID	H ii size	${{\mathscr{L}}}_{H\alpha }$	SFR	$\langle 12+\mathrm{log}({\rm{O}}/{\rm{H}})\rangle$ a	$\langle {T}_{e}\rangle$ a	$\langle {n}_{e}\rangle$ b	${R}_{e}$ c	${{\rm{W}}}_{U}$ d	$\mathrm{log}{M}_{\mathrm{dyn},\mathrm{turb}}$ e
	(arcsec)	(1038 erg s−1)	(M${}_{\odot }$ yr−1)		(104 K)	(cm−3)	(kpc)	(km s−1)	(${M}_{\odot }$)
ML4	1.4	30 ± 4	0.016 ± 0.002	7.76 ± 0.06	1.57	13	0.435 ± 0.012	62 ± 2	8.3 ± 0.03
HS0822	1.1	1.9 ± 0.3	0.00101 ± 0.00015	7.69 ± 0.07	1.56	38	0.027 ± 0.006	36 ± 2	6.61 ± 0.1
SBS0	1.0	7.1 ± 1	0.0038 ± 0.0005	7.63 ± 0.06	1.67	129	0.052 ± 0.011	38 ± 3	6.95 ± 0.12
SBS1	1.1	2.3 ± 0.3	0.00121 ± 0.00017	7.61 ± 0.09	1.44	84	0.059 ± 0.018	51 ± 3	7.27 ± 0.14
ML16	2.4	7.5 ± 1.1	0.004 ± 0.0006	7.92 ± 0.04	1.46	48	0.139 ± 0.005	35 ± 5	7.3 ± 0.12
UM461	1.2	1.32 ± 0.19	0.00071 ± 0.0001	7.38 ± 0.09	2.47	212	0.036 ± 0.009	42 ± 4	6.88 ± 0.14
N241	1.4	87 ± 13	0.047 ± 0.007	7.11 ± 0.28f	1.5g	73	0.035 ± 0.003	48 ± 3	6.99 ± 0.07
SBS2104	1.4	4.1 ± 0.7	0.0022 ± 0.0004	7.48 ± 0.08	1.58	113	0.069 ± 0.002	18 ± 5	6.4 ± 0.2
ML32	1.3	190 ± 30	0.103 ± 0.015	7.7 ± 0.05	1.68	13	0.434 ± 0.013	55 ± 4	8.19 ± 0.06

Notes.

^aMean value taken from Paper I. ^bInferred from the ratio between [S ii] λ6716 and [S ii] λ6731 (e.g., Osterbrock 1974). ^cHalf-light radius of the clump, measured from the FWHM of the Hα flux distribution. ^dVelocity dispersion of the clump, calculated from the Hα line of the clump-integrated spectra. ^eDynamical mass of the clump from Equation (5). Errors do not include the uncertainty in the subtraction of thermal motion. ^fFrom Sánchez Almeida et al. (2016). ^gAssumed to be similar to the value in the other objects.

Download table as:  ASCIITypeset image

The velocity at each spatial pixel along the major axis of the galaxy was fitted with the universal rotation curve by Salucci et al. (2007), i.e.,

$$\begin{eqnarray}&&U(d)={U}_{0}+{U}_{1}\displaystyle \frac{d-{d}_{0}}{\sqrt{{{\rm{\Delta }}}^{2}+{(d-{d}_{0})}^{2}}},\end{eqnarray} \tag{ 2 }$$

where U(d) stands for the observed velocity at a distance d, U₀ represents the systemic velocity, i.e., the velocity at the dynamic center $d={d}_{0}$, Δ gives a spatial scale for the central gradient of the rotation curve, and, finally, U₁ provides the amplitude of the rotational velocity.

The dynamical mass enclosed within a distance d, M(d), follows from the balance between the centrifugal force and the gravitational pull, and it is given by

$$\begin{eqnarray}&&M(d){\sin }^{2}i=(2.33\times {10}^{5}\,{M}_{\odot })(d-{d}_{0}){U}^{2}(d),\end{eqnarray} \tag{ 3 }$$

where i stands for the inclination of the disk along the line of sight, and where distances are in kiloparsecs and velocities in km s⁻¹. The velocity U in Equation (3) should be the component along the line of sight of the circular velocity, ${v}_{c}\,\sin i$. In the case of a purely stellar system, $U/\sin i$ and v_c differ because the stellar orbits are not circular, with the difference given by the so-called asymmetric drift (e.g., Hinz et al. 2001; Binney & Tremaine 2008). In our case, where we measure the velocity of the gas, they also differ because gas pressure gradients partly balance the gravitational force, and this hydrodynamical force has to be taken into account in the mechanical balance (e.g., Dalcanton & Stilp 2010). In a first approximation (e.g., Dalcanton & Stilp 2010; Read et al. 2016),

$$\begin{eqnarray}&&{v}_{c}^{2}\simeq {U}^{2}/{\sin }^{2}i+0.06\,{W}_{U}^{2}(d-{d}_{0})/{R}_{\star },\end{eqnarray} \tag{ 4 }$$

where the coefficient that gives the correction has been inferred assuming the gas density to drop exponentially with the distance to the center ($| d-{d}_{0}|$), with a length scale around three times the length scale of the stellar disks ${R}_{\star }$. This difference between the gaseous and the stellar disks is typical of XMPs (Filho et al. 2013). Since $| d-{d}_{0}| \sim {R}_{\star }$, Equation (4) predicts that, even for a large turbulent velocity (${W}_{U}\simeq U/\sin i$), the difference between v_c and $U/\sin i$ is small. Since we do not know the inclination angle of the galaxies, we cannot apply the correction for pressure gradients in Equation (4). Thus, we employ U in Equation (3), knowing that the resulting dynamical masses slightly underestimate the true masses by an amount given by ${W}_{U}/\sin i$ and Equation (4).

Our targets have distinct star-forming regions (see Figure 1), also denoted here as star-forming clumps or starbursts. We refer to the brightest clump in each galaxy as the main star-forming region, to distinguish it from other, fainter clumps that appear in some galaxies. They have been highlighted in Figure 1. The dynamical masses of the individual star-forming clumps are calculated from the velocity dispersion, assuming virial equilibrium,

$$\begin{eqnarray}&&{M}_{\mathrm{dyn},\mathrm{turb}}=(1.20\times {10}^{5}\,{M}_{\odot }){R}_{e}\,{{W}_{U}}^{2}.\end{eqnarray} \tag{ 5 }$$

The half-light radius, R_e, is calculated as in Paper I, from the FWHM of a 1D Gaussian fitted to the Hα flux of the star-forming clump, which is subsequently corrected for seeing (Table 1). W_U in Equation (5) is calculated from the Hα line that results from integrating all the spectra along the spatial axis, from $-{R}_{e}$ to $+{R}_{e}$, with the interval centered at the clump. R_e in Equation (5) is given in kiloparsecs and W_U in km s⁻¹. The use of W_U to estimate dynamical masses implies excluding thermal motions in the virial equilibrium equation (Equation (1)). This approximation may be more or less appropriate depending on the structure of the velocity field in the emitting gas. In practice, however, thermal motions are always smaller than W_U, so that including them does not significantly modify the masses estimated through Equation (5).

Our long-slit spectra only provide cuts across the star-forming regions. In order to obtain the total luminosity, ${{\mathscr{L}}}_{{\rm{H}}\alpha }$, we assume the region to be circular, with a radius R_e. Then the flux obtained by integrating the observed spectra from $-{R}_{e}$ to $+{R}_{e}$, ${F}_{{\rm{H}}\alpha }$, is corrected according to the ratio between the observed area, $2\,{R}_{e}\times 1\,\mathrm{arcsec}$, and the circular area, $\pi \,{R}_{e}^{2}$, which yields

$$\begin{eqnarray}&&{{\mathscr{L}}}_{{\rm{H}}\alpha }=2{\pi }^{2}\,{D}^{2}\,\displaystyle \frac{{R}_{e}}{1\,\mathrm{arcsec}}\,{F}_{{\rm{H}}\alpha },\end{eqnarray} \tag{ 6 }$$

with D the distance to the source. We take D from NED, as listed in Table 1. SFRs are inferred from ${{\mathscr{L}}}_{{\rm{H}}\alpha }$ using the prescription by Kennicutt & Evans (2012) given in Equation (20).

When the S/N is sufficiently good, the observed Hα profiles show faint broad emission, often in the form of multiple pairs of separate emission-line components (see Figure 3). The existence of pairs of components is particularly interesting from a physical standpoint, since they may reveal the feedback of star formation on the surrounding medium or the presence of a black hole (BH, to be discussed in Section 5.2).

Zoom In Zoom Out Reset image size

In order to determine the properties of the secondary lobes (number, wavelength shifts, relative fluxes, etc.), we performed a multi-Gaussian fit to Hα using the Python package LMFIT (Newville et al. 2014). The fitting is not unique unless we constrain the free parameters, which we do by attempting to reproduce Hα in the brightest pixel of the galaxy as a collection of discrete components, with widths similar to that of the central component. The procedure is detailed below. Several examples of the fit in various S/N conditions are given in Figure 3. We initialize the fit in the spatial pixel with highest S/N in continuum, where the secondary components clearly stand out. This pixel corresponds to the brightest star-forming region of the galaxy (the bright blue knots in Figure 1). Then the rest of the pixels are fitted, using as starting values the free parameters obtained in the fit to the adjacent pixel. The number of secondary components is set at the central brightest pixel, where the S/N needed to detect them is highest, and we fix this number for each galaxy. In order to follow each component across the spatial direction, we put two constraints on its spatial variation, namely, (1) the center does not deviate from the previous value by more than twice the wavelength sampling (∼1 Å), and (2) the velocity dispersion is bound between a minimum set by the spectral resolution and twice the value in the previous pixel. These constraints take account of the fact that neighboring pixels are not independent, since our sampling grants at least two pixels per spatial resolution element. As a sanity check, we verified that the constraint on the center of the components does not affect the inferred velocities. The maximum wavelength displacement between two adjacent pixels for all secondary components in the clumps is approximately half the wavelength sampling. Although the S/N decreases with the distance to the brightest pixel, we checked that each detected secondary component has a signal at least twice the noise in continuum in all the spectra used to characterize the star-forming regions.

The uncertainties in each parameter can be derived from the covariance matrix of the fit (e.g., Press et al. 1986). However, when the fitted parameters are near the values set by the constraints, the covariance matrix cannot be computed and the uncertainties provided by the Python routine turn out to be unreliable. Therefore, errors are estimated by running a Monte Carlo simulation. The best-fit profile is contaminated with random noise that has the observed S/N. This simulated observation is processed as a true observation, and the exercise is repeated 80 times with independent realizations of the noise. The rms fluctuations of the parameters thus derived are quoted as 1σ errors.

In order to estimate the mass, the mass loss rate, and the kinetic energy corresponding to the secondary components observed in the wings of Hα (Figure 3), we follow a procedure similar to that described in Carniani et al. (2015). It is detailed here both for completeness and to tune the approach to our specific needs.

We assume that the small emission-line components in the wings of the main Hα profile are produced by the recombination of H. Then the luminosity in Hα, ${L}_{{\rm{H}}\alpha }$, is given by

$$\begin{eqnarray}&&{L}_{{\rm{H}}\alpha }=\int \,f\,{n}_{e}\,{n}_{p}\,{j}_{{\rm{H}}\alpha }\,{dV},\end{eqnarray} \tag{ 7 }$$

where the integral extends to all the emitting volume, and n_e, n_p, and ${j}_{{\rm{H}}\alpha }$ represent the number densities of electrons and protons, and the emission coefficient, respectively. The symbol f stands for a local filling factor that accounts for a clumpy medium, so that only a fraction f is contributing to the emission. In a tenuous plasma, the emission coefficient is given by (Osterbrock 1974; Lang 1999)

$$\begin{eqnarray}&&{j}_{{\rm{H}}\alpha }\simeq 3.56\times {10}^{-25}\,{t}_{4}^{-1}\,\,\mathrm{erg}\,{\mathrm{cm}}^{3}\,{{\rm{s}}}^{-1},\end{eqnarray} \tag{ 8 }$$

where t₄ is the temperature in units of 10⁴ K. The mass of emitting gas turns out to be

$$\begin{eqnarray}&&{M}_{g}=\int \,f\,{m}_{{\rm{H}}}\,{n}_{{\rm{H}}}\,{X}^{-1}\,{dV},\end{eqnarray} \tag{ 9 }$$

where ${n}_{{\rm{H}}}$ is number density of H (neutral and ionized), ${m}_{{\rm{H}}}$ stands for the atomic mass unit, and X is the fraction of gas mass in H. If the gas is fully ionized¹¹ (${n}_{{\rm{H}}}={n}_{p}$), and the electron temperature and the H mass fraction are constant, then Equations (7)–(9) give

$$\begin{eqnarray}&&{M}_{g}=\displaystyle \frac{{m}_{{\rm{H}}}\,{L}_{{\rm{H}}\alpha }}{X\,{j}_{{\rm{H}}\alpha }\,\langle {n}_{e}{\rangle }_{1}},\end{eqnarray} \tag{ 10 }$$

with the mean electron density defined as

$$\begin{eqnarray}&&\langle {n}_{e}{\rangle }_{1}=\displaystyle \frac{\int \,f\,{n}_{p}\,{n}_{e}\,{dV}}{\int \,f\,{n}_{p}\,{dV}}.\end{eqnarray} \tag{ 11 }$$

Using astronomical units, and assuming X = 0.75 (which corresponds to the solar composition; e.g., Asplund et al. 2009), Equation (10) can be written as

$$\begin{eqnarray}&&{M}_{g}=3.15\times {10}^{3}\,{M}_{\odot }\,{t}_{4}\,\displaystyle \frac{{L}_{{\rm{H}}\alpha }/{10}^{38}\,\mathrm{erg}\,{{\rm{s}}}^{-1}}{\langle {n}_{e}{\rangle }_{1}/{10}^{2}\,{\mathrm{cm}}^{-3}}.\end{eqnarray} \tag{ 12 }$$

Assume that the motions associated with the weak emission-line components in the wings of the main profile trace some sort of expansion. The mass loss rate, ${\dot{M}}_{g}$, is defined as the total mass carried away per unit time by these motions. It can be computed as the flux of mass across a closed surface around the center of expansion,

$$\begin{eqnarray}&&{\dot{M}}_{g}=\int f\,{\rho }_{g}\,{\boldsymbol{v}}\cdot {\boldsymbol{d}}{\boldsymbol{\Sigma }},\end{eqnarray} \tag{ 13 }$$

where ${\rho }_{g}$ is the gas density at the surface and ${\boldsymbol{v}}$ represents the velocity of the flow. Assuming that motions are radial, and that the density and filling factor are constant at the surface, the previous integral can be simplified considering a spherical surface of radius R,

$$\begin{eqnarray}&&{\dot{M}}_{g}=4\pi \,{R}^{2}\,f\,{\rho }_{g}{v}_{\mathrm{out}},\end{eqnarray} \tag{ 14 }$$

with the radial speed ${v}_{\mathrm{out}}=| {\boldsymbol{v}}|$. If the moving mass forms a shell of width ${\rm{\Delta }}R$ (Figure 4), then ${M}_{g}=4\pi \,{R}^{2}\,f{\rho }_{g}\,{\rm{\Delta }}R$, so that the mass loss rate turns out to be

$$\begin{eqnarray}&&{\dot{M}}_{g}\simeq {M}_{g}\,{v}_{\mathrm{out}}\,{\rm{\Delta }}{R}^{-1}.\end{eqnarray} \tag{ 15 }$$

If the moving mass is distributed on a sphere, then ${M}_{g}=4\pi \,{R}^{3}\,f{\rho }_{g}/3$, and the mass loss rate is still given by Equation (15) provided that ${\rm{\Delta }}R/R\simeq 1/3$. In other words, by using ${\rm{\Delta }}R/R\simeq 0.3$ in Equation (15) one approximately considers the full range of geometries between shells and spheres (for further discussion on the assumptions made to estimate ${\dot{M}}_{g}$ see, e.g., Maiolino et al. 2012; Keller et al. 2014; or Carniani et al. 2015).

Zoom In Zoom Out Reset image size

Once the expanding shell geometry has been assumed, the radius of the shell is set by the gas mass and the electron density, since Equation (9) can be written as

$$\begin{eqnarray}&&{M}_{g}\simeq \displaystyle \frac{8\pi \,{m}_{H}}{1+X}\,\langle {n}_{e}{\rangle }_{2}\,\langle f{\rangle }_{3}\,{R}^{2}\,{\rm{\Delta }}R,\end{eqnarray} \tag{ 16 }$$

with the mean values given by

$$\begin{eqnarray*}\langle {n}_{e}{\rangle }_{2} & = & \displaystyle \frac{\displaystyle \int \,f\,{n}_{e}\,{dV}}{\displaystyle \int \,f\,{dV}},\\ \langle f{\rangle }_{3} & = & \displaystyle \frac{\displaystyle \int \,f\,{dV}}{\displaystyle \int \,{dV}}.\end{eqnarray*}$$

Equation (16) assumes the gas to be composed of fully ionized H and He, so that

$$\begin{eqnarray}&&{n}_{e}\simeq {n}_{{\rm{H}}}+2{n}_{\mathrm{He}}\simeq {n}_{{\rm{H}}}\displaystyle \frac{1+X}{2X},\end{eqnarray} \tag{ 17 }$$

which is a reasonable estimate considering that we are dealing with metal-poor gas, such that the metals do not contribute to the electron density. Equation (16) provides a way to estimate the size of the emitting region in terms of the electron density, the filling factor, and the relative width of the shell; explicitly,

$$\begin{eqnarray}&&{\left[\displaystyle \frac{R}{14.7\mathrm{pc}}\right]}^{3}=\left[\displaystyle \frac{{M}_{g}}{{10}^{4}\,{M}_{\odot }}\right]{\left[\displaystyle \frac{\langle {n}_{e}{\rangle }_{2}}{{10}^{2}{\mathrm{cm}}^{-3}}\displaystyle \frac{\langle f{\rangle }_{3}}{0.3}\displaystyle \frac{{\rm{\Delta }}R/R}{0.3}\right]}^{-1}.\end{eqnarray} \tag{ 18 }$$

The expansion velocity and the radius allow us to estimate the age of the expanding shell, or the time lag from blowout,

$$\begin{eqnarray}&&{\rm{Age}}=R/{v}_{\mathrm{out}},\end{eqnarray} \tag{ 19 }$$

where ${v}_{\mathrm{out}}$ is assumed to be constant in time. Since the expansion decelerates, Age in Equation (19) represents an upper limit to the true age.

A convenient (and intuitive) way of expressing the mass loss rate in Equation (15) is to rewrite the expression in terms of the star formation rate (SFR) estimated using the prescription by Kennicutt & Evans (2012),¹²

$$\begin{eqnarray}&&\displaystyle \frac{\mathrm{SFR}}{1\ {M}_{\odot }\ {\mathrm{yr}}^{-1}}\simeq \displaystyle \frac{{{\mathscr{L}}}_{{\rm{H}}\alpha }}{1.86\times {10}^{41}\,\mathrm{erg}\,{{\rm{s}}}^{-1}},\end{eqnarray} \tag{ 20 }$$

which depends only on the total Hα flux of the region, ${{\mathscr{L}}}_{{\rm{H}}\alpha }$. We parameterize the Hα flux in a weak emission-line component, ${L}_{{\rm{H}}\alpha }$, in terms of total Hα flux,

$$\begin{eqnarray}&&{L}_{{\rm{H}}\alpha }=\varepsilon {{\mathscr{L}}}_{{\rm{H}}\alpha },\end{eqnarray} \tag{ 21 }$$

with ε the scaling factor inferred from the fits described in Section 3.2. We do not correct ${{\mathscr{L}}}_{{\rm{H}}\alpha }$ for extinction in Equation (20) because the photometric calibration needed to infer it from the Balmer decrement is uncertain (Section 2.2); this is because the main Hα lobe and the weak emission components do not necessarily have the same extinction, and because the expected extinction is very low (reddening coefficient around 0.1; Sánchez Almeida et al. 2016). Using Equations (20) and (21), Equation (15) can be expressed as

$$\begin{eqnarray}&&\displaystyle \frac{{\dot{M}}_{g}}{\mathrm{SFR}}\simeq 0.20\displaystyle \frac{\varepsilon }{{10}^{-3}}\displaystyle \frac{{v}_{\mathrm{out}}}{{10}^{2}\,\mathrm{km}\,{{\rm{s}}}^{-1}}\,{t}_{4}{\left[\displaystyle \frac{{\rm{\Delta }}R}{3\mathrm{pc}}\displaystyle \frac{\langle {n}_{e}{\rangle }_{1}}{{10}^{2}{\mathrm{cm}}^{-3}}\right]}^{-1}.\end{eqnarray} \tag{ 22 }$$

Another important parameter derived from the gas mass and velocity is the kinetic energy involved in the outflow,

$$\begin{eqnarray}&&{E}_{k}=\displaystyle \frac{1}{2}{M}_{g}\,{v}_{\mathrm{out}}^{2}\simeq {10}^{51}\,\mathrm{erg}\,\displaystyle \frac{{M}_{g}}{{10}^{4}\,{M}_{\odot }}{\left[\displaystyle \frac{{v}_{\mathrm{out}}}{{10}^{2}\mathrm{km}{{\rm{s}}}^{-1}}\right]}^{2}.\end{eqnarray} \tag{ 23 }$$

The equation has been written in units of ${10}^{51}\,\mathrm{erg}$, which is the typical kinetic energy released in a core-collapse supernova (SN; e.g., Leitherer et al. 1999; Kasen & Woosley 2009), corresponding to the explosion of a massive star.

This energy range is also within reach of an accreting BH of intermediate mass, like the ${10}^{5}\,{M}_{\odot }$ BH seeds needed to explain the existence of supermassive BHs in the early universe (e.g., Volonteri 2010; Ferrara et al. 2014; Pacucci et al. 2016). Assuming that the BH releases kinetic energy at a fraction of the Eddington limit, ψ, and that the accretion event has a duration, ϒ, then

$$\begin{eqnarray}&&{E}_{\mathrm{BH}}\simeq {10}^{51}\,\mathrm{erg}\,\displaystyle \frac{\psi }{0.1}\,\displaystyle \frac{{\rm{\Upsilon }}}{25\,\mathrm{yr}}\,\displaystyle \frac{{M}_{\mathrm{BH}}}{{10}^{5}\,{M}_{\odot }},\end{eqnarray} \tag{ 24 }$$

with M_BH the BH mass (e.g., Frank et al. 2002, p. 398).

In principle, the observed faint emission-line pairs in the wings of Hα might also be generated in a rotating disk around a massive object (e.g., Epstein 1964). Then the approximation used to estimate the radius of the shell in Equation (18) still holds, except that

$$\begin{eqnarray}&&{\left[\displaystyle \frac{R}{27.6\mathrm{pc}}\right]}^{3}=\left[\displaystyle \frac{{M}_{g}}{{10}^{4}\,{M}_{\odot }}\right]{\left[\displaystyle \frac{\langle {n}_{e}{\rangle }_{2}}{{10}^{2}{\mathrm{cm}}^{-3}}\displaystyle \frac{\langle f{\rangle }_{3}}{0.3}\displaystyle \frac{{\rm{\Delta }}R/R}{0.3}\displaystyle \frac{{\rm{\Delta }}z/R}{0.3}\right]}^{-1},\end{eqnarray} \tag{ 25 }$$

where ${\rm{\Delta }}R$ represents the difference between the inner and outer radii, and ${\rm{\Delta }}z$ stands for the thickness of the disk. This size, together with the circular velocity, v_c, directly inferred from the spectra, allow us to estimate the mass of the central object, ${M}_{\bullet }$, from the balance between gravity and centrifugal force,

$$\begin{eqnarray}&&{M}_{\bullet }\simeq 7.0\times {10}^{7}\,{M}_{\odot }\,\displaystyle \frac{R}{30\,\mathrm{pc}}{\left(\displaystyle \frac{{v}_{c}}{{10}^{2}\mathrm{km}{{\rm{s}}}^{-1}}\right)}^{2}.\end{eqnarray} \tag{ 26 }$$

A separate, but important, result from the analysis of the spectra deals with the metallicity. Most of the XMPs that we study are characterized by having a drop in gas-phase metallicity associated with the main starburst (Paper I). We searched for possible relationships between such inhomogeneities in metallicity and the kinematic properties described above but found no correlation. What we found, however, is a notable lack of correlation between the metallicity (as traced by O/H) and the ratio N/O. This is apparent in Figure 8, where the variations of the metallicity, SFR, and N/O ratio are displayed. The N/O ratio was estimated from [N ii]λ5583 and [S ii] λλ6716, 6731 following the calibration by Pérez-Montero & Contini (2009), which is almost independent of dust reddening (it was not estimated in Paper I because of the difficulty in separating [N ii] from Hα, given the limited wavelength resolution of those spectra). Figure 8 shows that the N/O ratio remains fairly constant across the galaxy, even at the position of the main starbursts, where the Hα flux peaks and the metallicity drops. The lack of correlation between the ratios N/O and O/H represents the general behavior. This observational fact is consistent with the gas accretion scenario. If the accretion of metal-poor gas is triggering the observed starburst, then the fresh gas would locally reduce the O/H ratio in the star-forming region relative to the rest of the galaxy. However, mixing with external gas cannot modify the pre-existing ratio between metals, thus leaving the N/O ratio unchanged. This argument has been put forward before to support the gas accretion scenario (Amorín et al. 2010, 2012; Sánchez Almeida et al. 2014a). The value for $\mathrm{log}({\rm{N}}/{\rm{O}})$ that we infer, around −1.5, is also consistent with the N/O ratio found in metal-poor α-enhanced stars in the solar neighborhood (e.g., Israelian et al. 2004; Spite et al. 2005). Since the star formation is enhanced where the O/H ratio drops, the lack of correlation between the ratios N/O and O/H reveals a lack of correlation between N/O and the SFR. This property seems to be common among star-forming galaxies at all redshifts (e.g., Pérez-Montero et al. 2013), and supports the metal-poor nature of the gas that forms stars.

Zoom In Zoom Out Reset image size

The first significant result is the mere existence of multiple components in the wings of Hα, which are often paired, so that for each red component there is a blue component with a similar wavelength shift and amplitude (see the example in Figure 3(a)). This fact eliminates spatially unresolved, optically thin, and uniformly expanding shells as the source of emission in the line wings, since they produce a top-hat line profile without individual emission peaks (e.g., Zuckerman 1987, pp. 185–209; Cid Fernandes & Terlevich 1994; Tenorio-Tagle et al. 1996; and also the Appendix). However, the two-hump emission can be caused by a non-spherically symmetric shell or a shell having internal absorption (the other alternatives are analyzed and discarded in Section 5.2). Even though the observed amplitudes are similar, the blue component of each pair tends to be larger than the red component. Figure 9 shows the ratio between the areas of the red and blue components, and in all but one case, this ratio is one or smaller, considering error bars. Moreover, if one of the components is missing, it tends to be the red one. These two facts suggest a dusty expanding shell model, where the receding cap of the shell is obscured by the approaching cap (see the Appendix and Figure 17). These issues are discussed in detail in Section 5.2. The color code employed in Figures 9–15 assigns a color to each galaxy, and then the different symbols of the same color correspond to the various pairs of secondary components in that galaxy.

Zoom In Zoom Out Reset image size

Using the equations in Section 3.3, we have estimated the mass, mass loss rate, and kinetic energy associated with all the weak components in the Hα wings (Section 3.2, Figure 3). The mean values, as well as the errors, are given in Table 4. The velocities and the ratio between the secondary components and the central component, ε, are provided directly by the multi-Gaussian fit described in Section 3.2. In the case of paired components, the velocity is the mean of the absolute values of the two velocities, and ε is the total ε of the pair. The estimates also depend on the temperature and electron density of the emitting plasma, as well as on the relative thickness of the shell, ${\rm{\Delta }}R/R$, and the filling factor, f.

Table 4.  Parameters for the Weak Components in the Hα Wings^a

IDb	Velocityc	$\varepsilon$ d	Red/bluee	Mg	R	$\dot{{M}_{g}}$	Ek	Age
	(km s−1)			(103 ${M}_{\odot })$	(pc)	$({M}_{\odot }$ yr−1)	(1051 erg)	(104 yr)
ML4-1.f	98 ± 8	0.031 ± 0.012	⋯	37 ± 15	30 ± 4	0.41 ± 0.17	3.5 ± 1.5	30 ± 5
HS0822-1.0	114 ± 6	0.019 ± 0.005	0.6 ± 0.15	0.47 ± 0.13	4.9 ± 0.5	0.037 ± 0.011	0.061 ± 0.018	4.2 ± 0.5
HS0822-2f	−248 ± 9	0.0024 ± 0.0008	⋯	0.06 ± 0.02	2.5 ± 0.3	0.02 ± 0.008	0.036 ± 0.013	0.97 ± 0.12
SBS0-1.0	16 ± 8	0.28 ± 0.13	0.12 ± 0.03	8 ± 4	8.5 ± 1.3	0.05 ± 0.04	0.02123 ± 0.02	50 ± 30
SBS0-2.0	211 ± 0.4	0.0061 ± 0.0005	1.3 ± 0.6	0.18 ± 0.03	2.35 ± 0.13	0.054 ± 0.009	0.079 ± 0.013	1.09 ± 0.06
SBS0-3f	−370 ± 20	0.00015 ± 0.00013	⋯	0.004 ± 0.004	0.68 ± 0.2	0.008 ± 0.007	0.006 ± 0.005	0.18 ± 0.05
SBS1-1.0	98 ± 3	0.11 ± 0.03	0.51 ± 0.03	1.3 ± 0.4	5.3 ± 0.6	0.08 ± 0.03	0.12 ± 0.04	5.3 ± 0.6
ML16-1.0	110 ± 40	0.036 ± 0.006	3 ± 4	2.6 ± 0.6	8 ± 0.6	0.12 ± 0.05	0.3 ± 0.2	7 ± 2
UM461-1.0	116 ± 3	0.0165 ± 0.002	1.2 ± 0.5	0.08 ± 0.015	1.53 ± 0.09	0.021 ± 0.004	0.011 ± 0.002	1.29 ± 0.09
UM461-2.0	248 ± 3	0.00297 ± 7e-05	1.1 ± 0.5	0.014 ± 0.002	0.86 ± 0.04	0.014 ± 0.002	0.0088 ± 0.0013	0.34 ± 0.017
UM461-3.0	372 ± 4	0.00082 ± 7e-05	1.2 ± 1.4	0.0039 ± 0.0007	0.56 ± 0.03	0.0089 ± 0.0016	0.0054 ± 0.0009	0.148 ± 0.008
N241-1.0	115 ± 4	0.017 ± 0.005	0.5 ± 0.3	10 ± 3	10.7 ± 1.2	0.35 ± 0.12	1.3 ± 0.4	9.2 ± 1
N241-2.0	236 ± 9	0.0043 ± 0.0011	2 ± 3	2.4 ± 0.7	6.8 ± 0.6	0.28 ± 0.09	1.3 ± 0.4	2.8 ± 0.3
N241-3.0	395 ± 7	0.0012 ± 0.0005	0.2 ± 1.4	0.6 ± 0.3	4.4 ± 0.7	0.2 ± 0.1	1 ± 0.5	1.08 ± 0.17
SBS2104-1.0	112 ± 2	0.025 ± 0.002	1.1 ± 0.3	0.46 ± 0.09	3.4 ± 0.2	0.052 ± 0.01	0.057 ± 0.011	2.96 ± 0.19
SBS2104-2.0	246 ± 4	0.0038 ± 0.0002	1 ± 0.5	0.07 ± 0.013	1.8 ± 0.11	0.032 ± 0.006	0.042 ± 0.008	0.72 ± 0.04
SBS2104-3f	−372 ± 10	0.0013 ± 0.0007	⋯	0.024 ± 0.014	1.3 ± 0.2	0.024 ± 0.014	0.033 ± 0.019	0.33 ± 0.06
ML32-1.0	93 ± 9	0.062 ± 0.014	1.7 ± 0.4	490 ± 130	71 ± 6	2.2 ± 0.7	42 ± 14	74 ± 10

Notes.

^aThe errors in this table represent the standard deviation across the starburst (for velocity, ε, and area ratio) or errors for the value obtained from the mean Hα line of the starburst (for the rest). ^bThe number added to the galaxy name identifies each of the pairs of components in the line profiles, with "1.0" representing the pair of lowest velocity, "2.0" the second lowest velocity, and so on. ^cVelocity of the center of the component. It is negative when there is no red counterpart to the blue component. In paired components, it is the mean of the absolute value of the pair. ^dArea ratio between the secondary component and the main component. In paired components, it is the total ε of the pair. ^eArea ratio between the red and blue emissions in the case of paired components. ^fUnpaired component.

Download table as:  ASCIITypeset image

The electron temperature is taken from the data used to determine oxygen abundance in Paper I. Specifically, we use the mean over the main star-forming region of the galaxy. The electron density has been inferred from the ratio of the emission lines [S ii]λ6716 and [S ii] λ6731, following the standard procedure (e.g., Osterbrock 1974), using the mean value over the whole star-forming region. These two lines are present within the spectral range observed with the WHT. Electron densities and temperatures are listed in Table 2. For the thickness of the shell, we use ${\rm{\Delta }}R/R=0.3$ (see Section 3.3). The filling factor, f, is assumed to be one.

The velocities of the weak secondary components are typically in excess of 100 km s⁻¹ (see Table 4). These velocities are very large compared with the rotational and turbulent velocities of the galaxies (Section 4), and, hence, are almost certainly larger than the escape velocity of the galaxy, which is expected to be of the order of the rotational plus turbulent velocity in the disk (e.g., Binney & Tremaine 2008). Assuming that the components in the wings correspond to outflows in expanding shells, it is still unclear whether or not the involved mass will finally escape from the galaxy. It all depends on the ISM of the host galaxy and the pressure it exerts on the expanding shell. The issue is analyzed in Section 5.2, with the conclusion that the gas will likely escape. The mass involved in these motions is given in Table 4 and summarized in Figure 10. The masses are in the range between $10\,{M}_{\odot }$ and ${10}^{5}\,{M}_{\odot }$, and they represent a few per cent of the total gas mass in the H ii region (the parameter ε in Table 4 gives the flux ratio between the secondary and the main emission lobe, and it has been taken as a proxy for the mass ratio; see Equations (21) and (22)). These masses are very small compared to the stellar masses and the dynamical masses of the entire galaxies (Tables 1 and 3).

Zoom In Zoom Out Reset image size

Figure 10 shows that the mass loss rate scales with the gas mass, since the velocities and radii of the different components and galaxies are similar. The weak components have associated mass loss rates between 10⁻² and 1 M${}_{\odot }$ yr⁻¹. They are very significant compared with the SFRs, as is attested to in Figure 11, which shows the mass loading factor, ${\dot{M}}_{g}/\mathrm{SFR}$. Most mass loading factors are in excess of one, meaning that for every mass unit transformed into stars, several mass units are swept away. Whether this mass escapes the H ii region, and subsequently the galaxy, we cannot know from our observations. However, cosmological numerical simulations predict large mass loading factors in low-mass galaxies (10 is not uncommon; see, e.g., Peeples & Shankar 2011; Davé et al. 2012; Oppenheimer et al. 2012; Bournaud et al. 2014; Thompson & Krumholz 2016). Mass loading factors observed in high-redshift objects are generally lower than the ones we infer (say, of the order of 2, e.g., Martin et al. 2012; Newman et al. 2012), likely because they refer to massive galaxies ($\mathrm{log}({M}_{\star }/{M}_{\odot })\gt 10$). However, it is not uncommon to find factors up to 10 in local dwarfs (e.g., Martin 1999; Veilleux et al. 2005). A caveat is in order, though. We are assuming the density of the shell to be the same as the mean density of the star-forming region. If the expanding shell builds up mass by sweeping gas around the explosion site, its density is expected to be higher than the mean density of the region. This increase in density would reduce our estimates of the mass loading factor (see Equation (22)).

Zoom In Zoom Out Reset image size

Figure 12 shows the kinetic energy associated with the weak components. Similar to the mass loss rate, it scales with the mass of gas that is involved in the motion. Figure 12 also shows the energy typically associated with a core-collapse SN explosion, of the order of 10⁵¹ erg (e.g., Leitherer et al. 1999). The energies involved in the observed motions are usually around this value, albeit smaller; ML32 is the only exception. Therefore, energy-wise, the components may be tracing shells produced by individual SN explosions. ML32 would require some 10 SNe exploding simultaneously.

Zoom In Zoom Out Reset image size

Depending on the initial mass function, population synthesis models predict a SN rate between 0.02 and 0.005 yr⁻¹ for a constant SFR = 1 M${}_{\odot }$ yr⁻¹ (Leitherer et al. 1999). By using a value in between these two extremes, say 0.01 yr⁻¹/$({M}_{\odot }{\mathrm{yr}}^{-1})$, the SFRs of our star-forming regions (Table 2) allow us to estimate the expected SN rates if we assume the observed SFRs to be constant. The predicted SN rate turns out to be between 10⁻³ and ${10}^{-5}\,{\mathrm{yr}}^{-1}$, which corresponds to a time lag between SN explosions (${\tau }_{\mathrm{SN}}$, defined as the inverse of the SN rate) from 1000 yr to 0.1 Myr. These rates suffice to maintain the observed signals. Figure 13 provides ${\rm{Age}}/{\tau }_{\mathrm{SN}}$ for individual shells. ${\rm{Age}}\gt {\tau }_{\mathrm{SN}}$ in half of the cases, which implies that we should expect more than one SN signal at a time, as we observe.

Zoom In Zoom Out Reset image size

Ages are estimated from the measured velocity and the estimated radius of the shell (Equation (19)). These radii are small compared to the size of the star-forming regions. Figure 13 shows the radius relative to the effective radius of the region, with the typical ratio being around 10% or smaller. Therefore, it is important to keep in mind, for interpretation of the weak emission components in the wings of Hα, that there are one or several (probably non-concentric) expanding shells embedded in a large star-forming region. This substructure is common in giant H ii regions, such as 30 Doradus in the Large Magellanic Cloud (e.g., Meaburn 1984; Sabbi et al. 2013; Camps-Fariña et al. 2015). Shells and holes appear in Hubble Space Telescope images of the large starburst in Kiso 5639, a tadpole galaxy very similar to the XMPs studied here (Elmegreen et al. 2016).

We do not see a dependence of the properties of the weak emission-line components on the mean metallicity of the starburst, except perhaps for the mass loading factor. Figure 14 shows a scatter plot of metallicity versus mass loading factor, and it hints at an increase in the factor with increasing metallicity. The trend may be coincidental, due to the small number of points in the plot. However, it may also indicate the contribution of stellar winds to the energy and momentum driving the observed expansion. The winds of massive stars can provide a total kynetic energy comparable to the energy released during their SN explosions (e.g., Fierlinger et al. 2016), and these winds driven by radiation pressure have mass loss rates that increase with increasing metallicity (e.g., Martins 2015). The metallicities represented in Figure 14 correspond to the metallicity averaged over the starburst and, thus, due to the variation of metallicity across the star-forming region; their values are slightly higher than the values quoted in Paper I. ML16, the point of highest metallicity in Figure 14, is also the galaxy with no detected inhomogeneity in metallicity in Paper I.

Zoom In Zoom Out Reset image size

One can think of several alternatives to explain the existence of multiple paired emission-line components in the wings of Hα. The previous section analyzes the possibility of the components being produced by spatially unresolved expanding shells driven by SN explosions. Here we discuss other options, such as (1) spatially unresolved rotating structures, (2) spatially resolved expanding shells, (3) bipolar outflows, (4) blown-out spherical expansion leading to expanding rings, and (5) spatially unresolved expanding shells driven by BH accretion feedback.

• (1)  
  Spatially unresolved rotating disk-like structures. Rotation naturally produces double-horned profiles (e.g., Epstein 1964). Therefore, one of the possibilities to explain the existence of blue and red paired emission peaks is the presence of a gaseous structure rotating around a massive object. The mass of the required central object can be estimated from the size of the rotation structure and the velocity (Section 3.4). Figure 15 shows that the mass of the central object needs to be around three orders of magnitude larger than the gas mass involved in the motions, with central masses, ${M}_{\bullet }$, ranging from 10⁶–10⁸ M${}_{\odot }$. This mass turns out to be of the order of the stellar mass of the whole galaxy (Table 1), which rules out the stellar nature of the required central concentration. If, on the other hand, it were a BH, then it would have to have ${M}_{\bullet }\sim {M}_{\star }$, which is completely out of the so-called Magorrian relationship between the central BH mass and the stellar mass of the galaxy, for which ${M}_{\bullet }/{M}_{\star }\sim {10}^{-3}$ (e.g., Kormendy & Ho 2013). Thus, the exotic nature of the required central mass disfavors the interpretation of the weak emission in the wings of Hα as being caused by a disk orbiting around a massive central object, such as a BH.
• (2)  
  Spatially resolved expanding shells. The inner part of a spatially resolved expanding shell also produces double-peak emission-line profiles, since the outer parts of the shell, contributing to the line core, are excluded (e.g., Zuckerman 1987, pp. 185–209). We often see several pairs of emission peaks (Figure 3(a)); therefore, this interpretation requires the existence of several resolved concentric shells, which is unrealistic since the first shell would have already swept the ISM around the center of the expansion. Moreover, the size of the shell is estimated to be smaller than the typical spatial resolution, so that the expanding regions cannot be spatially resolved.
• (3)  
  Bipolar outflows. Insofar as they comply with the energetics required to explain the observed motions, bipolar flows cannot be easily discarded. They create two-horned profiles, and several of them coexisting in the resolution element could explain the presence of multiple paired components. In fact, bipolar flows are extreme cases of non-spherically symmetric shells, where double-horned line profiles appear naturally (e.g., Gill & O'Brien 1999).
• (4)  
  Blown-out spherical expansion leading to expanding rings. If the radius of the shell exceeds the scale height of the disk of the galaxy, then the shell may break up, leading to a ring-like expanding structure that produces two-horned profiles when observed edge-on. This possibility can be discarded since the ring would have to be as large as the disk thickness, which is incompatible with the small size of the expanding structures (Table 4).
• (5)  
  Spatially unresolved expanding shells driven by BH accretion feedback. Energy-wise, intermediate-mass BHs of ${10}^{5}\,{M}_{\odot }$ can also power the expanding shells that we characterize in Section 5.1 (see Equation (24)). However, the observations also disfavor this option for a number of reasons. It is unclear how a primordial BH seed can coincide spatially with a star-forming region unless there are many such BHs lurking in the galaxy or the BH is born in situ. ${10}^{5}\,{M}_{\odot }$ BHs can be produced by direct collapse in the early universe, but this explanation can be discarded, since they are very rare (e.g., van Wassenhove et al. 2010) and many of them are needed. On the other hand, $100\,{M}_{\odot }$ BHs are expected to be more common, being left as remnants from Population III stars formed in the early universe or even today (e.g., Kehrig et al. 2015). However, this alternative is also very unlikely, since H ii regions are short-lived (a few tens of millions of years; e.g., Leitherer et al. 1999), and a BH requires several hundred million years of continuous feeding to grow from 100 to ${10}^{5}\,{M}_{\odot }$ (see, e.g., Volonteri 2010). Finally, accreting BHs are expected to be luminous X-ray sources. Five of our XMPs were observed with the X-ray satellite Chandra: HS0822, SBS1, N241, SBS2109, and UM 461 (Kaaret et al. 2011; Prestwich et al. 2013). Only one of them shows a signal above the noise, namely SBS1, with an X-ray luminosity around ${10}^{39}\,\mathrm{erg}\,{{\rm{s}}}^{-1}$. This level of emission is consistent with young SNe (e.g., Filho et al. 2004; González-Martín et al. 2006; Dwarkadas & Gruszko 2012).

Zoom In Zoom Out Reset image size