Spectral energy distribution and radio halo of NGC 253 at low radio frequencies

TL;DR

This study presents low-frequency radio observations of NGC 253, modeling its spectral energy distribution, revealing a complex interplay of starburst and extended emissions, and detecting its large-scale radio halo with characteristic structures.

Contribution

First detailed broadband radio spectral modeling of NGC 253 at low frequencies, identifying absorption features and halo morphology using MWA data.

Findings

Central starburst region modeled as free-free absorbed synchrotron plasma.

Extended emission shows spectral flattening at low frequencies.

Detected large-scale X-shaped radio halo extending ~8 kpc.

Abstract

We present new radio continuum observations of NGC253 from the Murchison Widefield Array at frequencies between 76 and 227 MHz. We model the broadband radio spectral energy distribution for the total flux density of NGC253 between 76 MHz and 11 GHz. The spectrum is best described as a sum of central starburst and extended emission. The central component, corresponding to the inner 500pc of the starburst region of the galaxy, is best modelled as an internally free-free absorbed synchrotron plasma, with a turnover frequency around 230 MHz. The extended emission component of the NGC253 spectrum is best described as a synchrotron emission flattening at low radio frequencies. We find that 34% of the extended emission (outside the central starburst region) at 1 GHz becomes partially absorbed at low radio frequencies. Most of this flattening occurs in the western region of the SE halo, and may be indicative of synchrotron self-absorption of shock re-accelerated electrons or an intrinsic low-energy cut off of the electron distribution. Furthermore, we detect the large-scale synchrotron radio halo of NGC253 in our radio images. At 154 - 231 MHz the halo displays the well known X-shaped/horn-like structure, and extends out to ~8kpc in z-direction (from major axis).

Tables

Table 1. Table 1 Range of angular resolution and noise values of the GLEAM data.

$\nu$	synthesised beam		rms noise	background
(MHz)	bmaj$\times$bmin	PA	(mJy beam-1)	noise
	(arcmin2)	(deg)		(mJy beam-1)
76	$5.03\times 4.72$	$-18.8$	$107$	$-44$
227	$1.73\times 1.67$	$-26.0$	$12.8$	$-3.3$

Table 2. Table 2 Flux density measurements of background sources embedded in NGC 253 radio emission. All measurements are in the same absolute flux density scale of Baars et al. ( 1977 ) .

Frequency	background	background	background	Angular	References
$[$MHz$]$	source 1	source 2	source 3	resolution
	(mJy)	(mJy)	(mJy)	(arcsec2)
74	$263\pm {76}^{†}$			$75\times 75$	a, ※
150	$212\pm 19$	${15.6}^{\mathrm{up}}$	$162\pm 17$	$36\times 24$	b, ※
200	$233\pm 23$			$138\times 126$	c
330	$190\pm 10$			$72\times 72$	d
610	$84\pm 15$			$114\times 24$	e
843	$97\pm 10$			$47\times 43$	f
1465	$63.0\pm 2.5$			$66\times 38$	g
1465	$53.0\pm 3.0$	$17.2\pm 1.0$	$26\pm 2$	$30\times 30$	d
1490	$5.6\pm 0.5$		$53.0\pm 0.5$	$19.0\pm 0.5$	h
4850	$19.0\pm 1.0$	$6.6\pm 0.3$	$7.9\pm 0.4$	$30\times 30$	d
8350	$9.5\pm 0.5$			$84\times 84$	d

Table 3. Table 3 Results of the spectral modeling of the background sources.

Background	polynomial	$S_{1GHz}$	$\alpha$	$c_1$	dof	${\chi }^2$
source	order	(mJy)
1	2	$80.7\pm 1.9$	$0.79\pm 0.02$	$-0.25\pm 0.04$	7	19.2
2	2	$21.9\pm 1.0$	$0.32\pm 0.07$	$-0.63\pm 0.11$	1	3.5
3	1	$33.1\pm 1.4$	$0.88\pm 0.03$	–	1	2.8

Table 4. Table 4 Total integrated flux density measurements of NGC 253. Uncertainties associated with the GLEAM, MWA/EoR and TGSS measurements include the fitting error and the uncertainty of the absolute flux density calibration, for the remaining data the errors quoted in references are adopted. All measurements are brought to the same absolute flux density scale of Baars et al. ( 1977 ) . If no information on the largest angular scale of the final images was reported, we state ‘no info’, and provide a theoretical value (upper limit, marked as ‘th’) for those measurements for which information on interferometer’s shortest baseline was provided, and assuming that adequate u ​ v 𝑢 𝑣 uv coverage has been achieved and no minimum- u ​ v 𝑢 𝑣 uv cut has been applied.

Frequency	total	Angular	Largest	References
$[$MHz$]$	(Jy)	resolution	angular scale
		(arcmin2)	(deg)
76	$23.9\pm 2.0$	$5.1\times 4.7$	29	a
80	$23.7\pm 4.0$	$3.7\times 3.7$	no info	b,c
84	$22.6\pm 1.9$	$5.1\times 4.7$	27	a
92	$20.9\pm 1.7$	$5.1\times 4.7$	24	a
99	$20.0\pm 1.6$	$5.1\times 4.7$	23	a
107	$20.3\pm 1.7$	$5.1\times 4.7$	21	a
115	$20.0\pm 1.6$	$5.1\times 4.7$	19	a
122	$18.8\pm 1.6$	$5.1\times 4.7$	18	a
130	$17.8\pm 1.5$	$5.1\times 4.7$	17	a
143	$17.9\pm 1.5$	$5.1\times 4.7$	16	a
151	$16.6\pm 1.4$	$5.1\times 4.7$	15	a
158	$16.1\pm 1.3$	$5.1\times 4.7$	14	a
166	$15.7\pm 1.3$	$5.1\times 4.7$	13	a
169	$15.0\pm 1.2$	$5.1\times 4.7$	13	a
174	$16.4\pm 1.3$	$5.1\times 4.7$	13	a
181	$15.7\pm 1.3$	$5.1\times 4.7$	12	a
189	$15.3\pm 1.2$	$5.1\times 4.7$	12	a
197	$15.4\pm 1.2$	$5.1\times 4.7$	11	a
204	$15.8\pm 1.3$	$5.1\times 4.7$	11	a
212	$15.1\pm 1.2$	$5.1\times 4.7$	11	a
220	$14.7\pm 1.2$	$5.1\times 4.7$	10	a
227	$15.0\pm 1.2$	$5.1\times 4.7$	10	a
330	$16.5\pm 1.9$	$1.2\times 1.2$	1.2	d, e
408	$15.7\pm 1.9$	$2.9\times 2.86$	no info	f
468	$15.1\pm 1.5$	$2.1\times 2.1$, $5.2\times {5.2}^{♢}$	extr. single dish♣	c, g
610	$9.4\pm 0.6$	$1.9\times 0.4$	no info	h
843	$9.0\pm 0.9$	$0.72\times 0.78$	no info (th: 1.1)	i
960	$8.0\pm 0.12$	$20.2\times 20.2$	single dish	c, g
1100	$6.7\pm 0.08$	$4.2\times 4.2$	no info, dense core	j
1200	$6.68\pm 0.10$	$3.4\times 3.4$	no info, dense core	j
1300	$6.22\pm 0.07$	$3.2\times 3.2$	no info, dense core	j
1400	$5.89\pm 0.16$	$2.9\times 2.9$	no info, dense core	j
1410	$6.12\pm 0.12$	$15.5\times 15.5$	single dish	c, g
1430	$5.7\pm 0.5$	$0.91\times 0.83$	no info (th: 0.8)	h
1465	$5.9\pm 0.1$	$1.1\times 0.63$	extr. single dish♣	k
1465	$6.3\pm 1.1$	$0.5\times 0.5$	0.25	d, e
1490	$5.6\pm 0.5$	$0.9\times 0.9$	0.27	m
2650	$3.85\pm 0.12$	$8.3\times 8.3$	single dish	c, g
2695	$4.26\pm 0.14$	$4.9\times 4.9$	single dish	n
2700	$3.49\pm 0.12$	$8.0\times 8.0$	single dish	c, g
4850	$2.93\pm 0.13$	$2.7\times 2.7$	single dish	n
4850	$2.71\pm 0.14$	$0.5\times 0.5$	single dish	e
4850	$2.69\pm 0.10$	$4.2\times 4.2$	single dish	p
5009	$2.50\pm 0.23$	$4.0\times 4.0$	single dish	c, g
5009	$2.12\pm 0.09$	$4.0\times 4.0$	single dish	c, r
8350	$1.66\pm 0.08$	$1.4\times 1.4$	single dish	e
8700	$2.06\pm 0.12$	$1.5\times 1.5$	single dish	n
10550	$1.98\pm 0.18$	$1.2\times 1.2$	single dish	s
10700	$1.95\pm 0.15$	$1.2\times 1.2$	single dish	t

Table 5. Table 5 Flux density measurements of NGC 253 central starburst region. All measurements are in the same absolute flux density scale of Baars et al. ( 1977 ) .

Frequency	central	Angular	References
$[$MHz$]$	component	resolution
	(Jy)	(arcsec2)
150	$2.16\pm 0.15$	$36\times 24$	a, ※
330	$2.67\pm 0.16$	$33\times 21$	b
610	$2.3\pm 0.2$	$114\times 24$	c
1413	$2.33\pm 0.14†$	$3\times 1.8$	d
1450	$2.07\pm 0.04$	$33\times 21$	b
1465	$2.04\pm 0.10$	$30\times 30$	e
1660	$1.96\pm 0.04$	$33\times 21$	b
4520	$1.36\pm 0.04$	$33\times 21$	b
4850	$1.27\pm 0.06$	$30\times 30$	e
4890	$1.30\pm 0.04$	$33\times 21$	b
6700	$1.13\pm 0.04$	$37\times 37$	f
7000	$1.04\pm 0.04$	$35\times 35$	f
8090	$0.93\pm 0.03$	$33\times 21$	b
8350	$0.98\pm 0.05$	$84\times 84$	e
8470	$0.89\pm 0.03$	$33\times 21$	b

Table 6. Table 6 Radio spectral indices ( α 𝛼 \alpha ) measured for selected regions of NGC 253 (Figure  6 ). Radio spectral index maps are used at 169 MHz (MWA/EoR0; S 169 ​ M ​ H ​ z subscript 𝑆 169 M H z S_{\rm 169MHz} ), 200 MHz (GLEAM; S 200 ​ M ​ H ​ z subscript 𝑆 200 M H z S_{\rm 200MHz} ), and 1.465 GHz convolved to the 169 MHz ( S 1.4 ​ GHz 169 ​ r ​ e ​ s superscript subscript 𝑆 1.4 GHz 169 r e s S_{\rm 1.4GHz}^{\rm 169res} ) and separately to 200 MHz ( S 1.4 ​ GHz 200 ​ r ​ e ​ s superscript subscript 𝑆 1.4 GHz 200 r e s S_{\rm 1.4GHz}^{\rm 200res} ) angular resolution.

Region	${\alpha }_{1.4\mathrm{GHz}}^{169MHz}$	$S_{169MHz}$	$S_{1.4\mathrm{GHz}}^{169res}$	${\alpha }_{1.4\mathrm{GHz}}^{200MHz}$	$S_{200MHz}$	$S_{1.4\mathrm{GHz}}^{200res}$
		(mJy)	(mJy)		(mJy)	(mJy)
1	$0.54\pm 0.06$	$242\pm 26$	$75\pm 10$	$0.57\pm 0.06$	$233\pm 25$	$75\pm 10$
2	$0.63\pm 0.04$	$586\pm 59$	$149\pm 13$	$0.68\pm 0.03$	$581\pm 51$	$149\pm 13$
3	$0.57\pm 0.02$	$883\pm 89$	$257\pm 11$	$0.62\pm 0.02$	$884\pm 73$	$257\pm 11$
4	$0.67\pm 0.04$	$415\pm 42$	$97\pm 9$	$0.72\pm 0.04$	$412\pm 36$	$97\pm 9$
5		$6558\pm 656$	$3341\pm 14$		$6589\pm 528$	$3341\pm 14$
6	$0.60\pm 0.07$	$224\pm 24$	$61\pm 10$	$0.70\pm 0.07$	$246\pm 25$	$61\pm 10$
7	$0.58\pm 0.03$	$488\pm 49$	$139\pm 10$	$0.63\pm 0.03$	$487\pm 42$	$139\pm 10$
8	$0.50\pm 0.02$	$616\pm 62$	$211\pm 9$	$0.51\pm 0.02$	$578\pm 49$	$211\pm 9$
9	$0.53\pm 0.06$	$216\pm 23$	$69\pm 9$	$0.51\pm 0.06$	$192\pm 21$	$69\pm 9$

Table 7. Table 7 Projected maximum vertical extent of the NGC 253 disk and halo at 169 MHz ( h ℎ h ) measured at a distance x 𝑥 x from the nuclear region along the major axis, where 0 is centred on the nucleus of the galaxy. The NE direction along the major axis is negative and SW is positive. The extent is measured perpendicular to the major axis (PA = − 38 ∘ absent superscript 38 =-38^{\circ} ) separately for the North and South side of the disk and halo. The uncertainties on the measurements are 21.6 arcsec (equivalent to 0.4 kpc). See § V.1 for discussion.

$x$	$x$	$h$ North	$h$ North	$h$ South	$h$ South
(arcsec)	(kpc)	(arcsec)	(kpc)	(arcsec)	(kpc)
-658	-11.8	113	2.0	100	1.8
-526	-9.5	213	3.8	211	3.8
-395	-7.1	262	4.7	327	5.9
-263	-4.7	396	7.1	440	7.9
-132	-2.4	416	7.5	367	6.6
0	0.0	311	5.6	291	5.2
132	2.4	270	4.9	291	5.2
263	4.7	297	5.5	323	5.8
395	7.1	279	5.0	260	4.7
526	9.5	201	3.6	171	3.1
658	11.8	56	1.0	100	1.8

Equations

$$\text{log}(S_{\nu })$$

$$=A_0+A_1\text{log}\left(\frac{\nu }{{\nu }_0}\right)+...+A_n{\text{log}}^n\left(\frac{\nu }{{\nu }_0}\right),$$

$$S_{\nu }=\prod _{i=0}^{n}10^{A_i{\text{log}}^i(\nu /{\nu }_0)},$$

$$S_{\nu }=S_{\tau =1}{\left(\frac{\nu }{{\nu }_{\tau =1}}\right)}^{-\alpha }\left(\frac{1-e^{-\tau (\nu )}}{\tau (\nu )}\right),$$

$$\tau (\nu )=(\nu /{\nu }_{\tau =1}{)}^{-(\alpha +2.5)},$$

$$S_{\nu }=S_{\tau =1}{\left(\frac{\nu }{{\nu }_{\tau =1}}\right)}^2(1-e^{-{\tau }_{\mathrm{ff}}(\nu )}),$$

$${\tau }_{\mathrm{ff}}(\nu )=(\nu /{\nu }_{\tau =1}{)}^{-2.1}.$$

$$S_{\nu }=S_0{\left(\frac{\nu }{{\nu }_0}\right)}^{-\alpha }\left(\frac{1-e^{-{\tau }_{\mathrm{ff}}(\nu )}}{{\tau }_{\mathrm{ff}}(\nu )}\right).$$

$${\chi }^2=\sum _i^n{\left(\frac{{\text{M}}_i-{\text{data}}_i}{{\text{error}}_i}\right)}^2$$

$$\text{ln}\mathcal{L}(\theta )=-\frac{1}{2}\sum _i^n\left[\frac{({\text{M}}_i-{\text{data}}_i{)}^2}{{\text{error}}_i^2}+\text{ln}(2\pi {\text{error}}_i^2)\right].$$

$$\Delta \text{ln}(Z)=\text{ln}(Z_2)-\text{ln}(Z_1)$$

$$Z_{1,2}=\int \int ...\int \mathcal{L}(\theta )\Pi (\theta )d(\theta ).$$

$$S_{\nu }^{\mathrm{mod}}=(S_{\nu }^{\mathrm{E}}+S_{\nu }^{\mathrm{C}}{)}^{\mathrm{tot}}+(S_{\nu }^{\mathrm{C}}{)}^{\mathrm{cor}},$$

$$S_{\nu }^{\mathrm{E}}=S_0^{\mathrm{ext}}{\left(\frac{\nu }{{\nu }_0}\right)}^{-{\alpha }^{\mathrm{ext}}},$$

$$S_{\nu }^{\mathrm{C}}=S_0^{\mathrm{cor}}{\left(\frac{\nu }{{\nu }_0}\right)}^{-{\alpha }^{\mathrm{cor}}}\left(\frac{1-e^{-{\tau }_{\mathrm{ff}}(\nu )}}{{\tau }_{\mathrm{ff}}(\nu )}\right).$$

Full text

I. INTRODUCTION

Observing at low radio frequencies ($\lesssim 0.5$ GHz) is of a particular value; low energy and old plasma can be revealed, tracing and constraining physical conditions in galaxies. In star forming galaxies the low surface brightness plasma forms e.g. extended halos associated with winds and large scale magnetic fields, or diffuse emission from galactic disks. Furthermore, measurements at low radio frequencies can help to distinguish for instance between thermal and non-thermal plasma, and their absorbing mechanisms, responsible for the level of observed radio emission. It is expected that the Square Kilometre Array (SKA) will unravel a large star-forming galaxy population (e.g. Beswick et al., 2015; Jarvis et al., 2015), but before we can embark on a large scale study of star-forming and starburst galaxies and their evolution with continuum radio surveys, we need to understand the origin of the complex radio spectral energy distributions and morphologies of these galaxies. Nearby objects are ideal laboratories for this task.

NGC 253 is the dominant galaxy in the nearby Sculptor Group, at a distance of 3.94 Mpc from the Local Group Karachentsev et al. (2003) and velocity $cz=240$ km s*-1*. It is an almost edge-on SBc type galaxy de Vaucouleurs et al. (1991) observed at an inclination of $78.5^{\circ}$ Pence (1980) and is considered a prototype of nuclear starburst galaxies Rieke et al. (1980). Its estimated stellar mass is $\sim 4\times 10^{10}$ M*⊙, with a prominent stellar halo of $2.5\times 10^{9}$ M⊙* extending up to 30 kpc above the disk Bailin et al. (2011). As one of the closest and most prominent galaxies, NGC 253 has been extensively studied in all wavelengths, including broadband radio continuum, polarization and Hi observations (Beck et al., 1979; Klein et al., 1983; Hummel et al., 1984; Carilli et al., 1992; Beck et al., 1994; Ulvestad & Antonucci, 1997; Lenc & Tingay, 2006; Heesen et al., 2009a, b, 2011; Lucero et al., 2015, among others).

Radio emission from starburst galaxies originates from two principal components: the non-thermal synchrotron emission from relativistic electrons spiralling in the interstellar magnetic field, and the thermal emission from electrons colliding with ions in the ionized interstellar medium (ISM) around hot stars. The sources of the non-thermal emission are predominantly cosmic rays accelerated by supernova remnants (SNR) that in NGC 253 ultimately create a prominent synchrotron radio halo Carilli et al. (1992). Studies of the NGC 253 magnetic field suggest the disk wind model and large-scale dynamo action to be shaping the vertical structure of the field, which in turn enhances the cosmic ray transport through a collimation of strong, starburst driven superwind Beck et al. (1994); Pietsch et al. (2000); Strickland et al. (2002); Heesen et al. (2009a, b, 2011).

The starburst region of NGC 253 is violently active; the supernova rate of the inner 300 pc of the galaxy is estimated to be between $0.14$ and $2.4$ yr*-1*, and the star formation rate is $\sim 5$ M*⊙* yr*-1* Lenc & Tingay (2006); Rampadarath et al. (2014); Bendo et al. (2015); Lucero et al. (2015). It has been suggested that up to half of the radio sources in the central starburst region are dominated by thermal emission: i.e. Hii regions characterized by a flat radio spectral index111Radio spectral index $\alpha$ is defined such that the flux density $S_{\nu}$ at a frequency $\nu$ is $S_{\nu}\propto\nu^{-\alpha}$. $\alpha\simeq 0.1$ and including at least one large supercluster of stars Ulvestad & Antonucci (1997); Keto et al. (1999). Outside the central starburst region the radio emission at GHz frequencies is dominated by steep spectrum diffuse emission and SNRs, but several strong thermal sources are detected Ulvestad (2000). Based on integrated radio continuum spectra, Niklas et al. (1997) estimated 10% of the NGC 253 flux density at 1 GHz to be of thermal origin, increasing to 35% at 10 GHz.

At low radio frequencies both of these principal components become pronounced. Synchrotron emission has steep spectrum, becoming dominant at sub-GHz frequencies due to the population of old, low energy electrons. However, such emission may be also subject to self-absorption in the case of compact objects. Thermal emission also becomes increasingly more absorbed with decreasing frequency. The free-free absorption in the central starburst of NGC 253 has previously been measured Carilli (1996); Tingay (2004).

Here, we present extensive low radio frequency ($<230$ MHz) imaging of NGC 253 obtained with the Murchison Widefield Array (MWA; Bowman et al., 2013; Tingay et al., 2013). Our images are some of the deepest yet at these frequencies, and at low angular resolution they are especially sensitive to large-scale diffuse structure, allowing us to investigate the extent and frequency dependence of the radio halo. The paper is structured as follows. Our radio data and methods, including assumed models of radio spectra and model fitting, are described in §II and §III respectively. Results are presented in §IV. The synchrotron radio halo of NGC 253 is discussed in §V.1. We discuss low frequency radio emission from NGC 253, its radio spectral energy distribution and radio spectral maps in §V.2. Conclusions are given in §VI.

II. Observations and data reduction

We use radio continuum data from the Galactic and Extragalactic All-Sky MWA Survey (GLEAM; Wayth et al., 2015) and the MWA Epoch of Reionization experiment (MWA/EoR; Bowman et al., 2013; Pober et al., 2016). The GLEAM survey provides unprecedented spectral coverage between 72 and 231 MHz, while the MWA/EoR image at 169 MHz is almost twice as deep as the most sensitive GLEAM image at 200 MHz (rms noise 4.1 mJy beam*-1* and 7.3 mJy beam*-1* respectively). In addition, the data have been observed and processed independently, providing a verification of our flux density calibration.

II.1. The GaLactic and Extragalactic All-Sky MWA Survey (GLEAM)

The GLEAM survey observed the entire radio sky south of declination $+30^{\circ}$ at an angular resolution of approximately 1.7 arcmin (227 MHz) to 5 arcmin (76 MHz). At 154 MHz the GLEAM survey is sensitive to structures up to $10$ deg in angular scale, and has an instantaneous field of view of $25\times 25$ deg2. The observations were made in a meridian drift scan mode covering frequencies between 72 and 231 MHz with bandwidths of 7.68 MHz grouped in five 30.72 MHz-wide bands. These bands, centred on 87.7, 118.4, 154.2, 185.0 and 215.7 MHz (hereafter 88, 118, 154, 185 and 216 MHz), were observed sequentially as 112 sec snapshots; each frequency was observed every 10 min. During a night typically 8–10 h in hour angle were observed. Frequencies between 134 and 137 MHz were avoided due to satellite interference. For more details on the survey parameters and strategy see Wayth et al. (2015).

Here we use GLEAM data from the first year of observing (Data Release 1 from 2013 August – 2014 June; Hurley-Walker et al., 2017). The sky area covering NGC 253 was observed on 2013 August 10 and 2013 November 25. The full data reduction process is described in detail in Hurley-Walker et al. (2017); here we summarize only the main calibration and imaging steps.

The correlated data were first pre-processed with the cotter pipeline which performs flagging of data affected by radio frequency interference (RFI) and averaging of the data to 1s time and 40 kHz frequency resolution Offringa et al. (2010, 2015). Standard calibration (phase and amplitude bandpass calibration) was done with CASA (Common Astronomy Software Applications package; McMullin et al., 2007). Imaging and self-calibration were then performed using WSClean imager Offringa et al. (2014) that corrects for wide field $w$-term effects. Images of a 7.68 MHz bandwidth at 20 frequencies continuously distributed between 72 and 231 MHz (avoiding 134–137 MHz) and using a robust weighting $r=-1.0$ Briggs (1995) were then created. Deconvolution has been performed at this stage, and details are provided in Hurley-Walker et al. (2017).

The primary beam correction of our GLEAM observations was done with the Sutinjo et al. (2015) model down to the 10% level of the beam response. An additional calibration stage was necessary to correct for residual declination dependence of the flux density scale in the final mosaics arising from the limited accuracy of the adopted primary beam model. This was done by comparing flux density measurements of all unresolved sources extracted from GLEAM images above $8\sigma$ rms noise level to their radio spectra as predicted by three catalogues: VLA Low-Frequency Sky Survey redux (VLSSr; Lane et al., 2014), MRC and NRAO VLA Sky Survey (NVSS; Condon et al., 1998). The absolute flux density scale of the GLEAM images is accurate to 8%, which is included in the quoted uncertainties of the measurements (for details see Hurley-Walker et al., 2017).

Images for each of five central frequencies centered on frequencies of 88, 118, 154, 185 and 216 MHz and of bandwidth 30.72 MHz were made. The two highest frequency images are further combined to create a ‘deep’ image at 200 MHz with a 61.4 MHz bandwidth. We also use the 7.68-MHz images for construction of the high resolution radio spectrum of NGC 253. The final synthesised beam sizes, rms and background noise levels in the deep 200 MHz image are $2.22\times 2.12$ arcmin2, ${\rm PA}=-78^{\circ}$, 11 mJy beam*-1* and 7.3 mJy beam*-1* respectively, and their range between the lowest and highest GLEAM frequencies is listed in Table 1.

II.2. MWA Epoch of Reionization (EoR) data

The observed MWA/EoR field that contains the Sculptor Group (EoR0 field) is centered on $\text{RA}=0^{\rm h}$, $\text{Dec}=-27^{\circ}$, and was observed for a total of 30 hours between August and October 2013 in a combination of a tracking and drifting modes. In this hybrid mode the telescope tracks a set of discrete pointing centers through which the field of interest is drifting. The observations cover frequencies between $138.9-197.7$ MHz observed as two bands (low and high) each with an instantaneous bandwidth of $\Delta\nu=30.72$ MHz.

The correlated MWA/EoR0 data were pre-processed with the cotter pipeline Offringa et al. (2015) and averaged to 4s time and 40 kHz frequency resolution. Calibration of the data was performed as a direction-independent self-calibration using the mitchcal tool Mitchell et al. (2008) and was based on a bootstrapped sky model. The initial sky model was generated from the MWA Commissioning Survey Hurley-Walker et al. (2014), the MRC catalogue and the Sydney University Molonglo Sky Survey (SUMSS; Mauch et al., 2003). Imaging was performed with the WSClean software that corrects for the non-zero $w$-term effects. During the imaging process 2,500 sources were peeled, and the images were created with a uniform weighting. The primary beam, and so the flux density scale, was corrected by applying the Sutinjo et al. (2015) model. As shown by Hurley-Walker et al. (2014), this model is accurate to 10%, hence we add this error in quadrature to the quoted uncertainties of our measurements.

The final image used in this paper is centred at 169.6 MHz (thereafter 169 MHz) with a total bandwidth of $\Delta\nu=58.8$ MHz, synthesised beam size $2.3\times 2.3$ arcmin2 and rms noise $4.1$ mJy beam*-1*. The calibration and imaging process of the EoR0 data is presented and discussed in detail in Offringa et al. (2016).

II.3. Other low frequency radio surveys

There are additional two all-sky low frequency radio surveys that include NGC 253: the 74 MHz VLSSr Lane et al. (2014) and the 150 MHz Tata Institute of Fundamental Research (TIFR) Giant Metrewave Radio Telescope (GMRT) Sky Survey (TGSS) Alternative Data Release 1 (ADR1; Intema et al., 2016).

The TGSS survey observed the whole radio sky north of declination $-53^{\circ}$ at a frequency 150 MHz (bandwidth $\Delta\nu=16.7$ MHz) at an angular resolution $25\times 25/\text{cos}(\delta-19^{\circ})$ arcsec2 at declinations south of $+19^{\circ}$. The instantenous field of view of the survey at half power at 150 MHz is $3.1\times 3.1$ deg2, with sensitivity to structures up to 68 arcmin in angular scale Intema et al. (2016). Since the absolute flux density calibration of the TGSS ADR1 may be uncertain up to 50% in some sky regions222http://tgssadr.strw.leidenuniv.nl/doku.php, we independently verified the calibration in the area of NGC 253. We selected unresolved sources with flux density $>1$ Jy from the $5\times 5$ deg2 mosaic that included NGC 253 (R03_D17). We compared the TGSS ADR1 flux densities of these sources with the predicted values based on the spectral modeling in which we used the VLSSr, GLEAM (deep 200 MHz), MRC and NVSS surveys. We found that the TGSS mosaic required scaling by a factor 1.02 in flux density, and the absolute flux density calibration was accurate to 7%; we further added this error in quadrature to the quoted uncertainties of our measurements.

The VLSS survey Cohen et al. (2007) observed the radio sky north of $-30^{\circ}$ at a frequency 74 MHz. Here we use the recent re-reduction of the survey data, the VLSSr Lane et al. (2014). VLSSr images have an angular resolution of $75\times 75$ arcsec2 and a theoretical sensitivity to structures of 13–37 arcmin in angular scale.

We find that neither TGSS nor VLSSr are sensitive to the extended emission of NGC 253 (Figure 1). For this reason we use the TGSS data for the flux density measurement of the central starburst region, and both TGSS and VLSSr for measurements of the flux densities of background sources only.

III. Methods

III.1. Flux density measurements

Measurements of the total flux density of NGC 253 were performed with CASA task imstat that provides a summed flux density within a specified regions of the image corrected for the synthesised beam. We masked all pixels below $2.6\sigma$ local rms noise level Hales et al. (2012). For point sources the flux density was measured with the AIPS task jmfit; for each unresolved source we fit for two components, a Gaussian and a zero-level with a slope. The absolute flux density scale is set to the Baars et al. (1977) scale.

III.2. Models of radio spectra

Radio sources often show simple spectra that can be approximated by a power-law. However, at low radio frequencies (a few hundred MHz and below) radio spectra are consistently more curved Laing & Peacock (1980); Hummel (1991); Deeg et al. (1993); Blundell et al. (1999); Duffy & Blundell (2012); Marvil et al. (2015) until a turnover frequency below which the spectra become inverted. The spectral turnover is typically caused by either synchrotron self-absorption and/or thermal free-free absorption (Callingham et al., 2015, and references therein). If there is no evidence for a spectral turnover in the radio spectra we construct here (see §III.4 for model selection method), we proceed with fitting a polynomial. The curved radio spectra are then modeled with an $n$th-order polynomial, that in the logarithmic scale takes a general form of

[TABLE]

where $A_{0}$ is an offset parameter (equivalent to log($S_{0}$) in the simple power-law case), $A_{1}$ is the spectral index $-\alpha$, and $A_{n}$ are curvature parameters ($c_{n}$). In the linear space the model takes the following form

[TABLE]

which we use in our modeling to preserve Gaussian noise characteristics of the measurements.

Where the data suggest or show a spectral turnover, the following models are tested: synchrotron self-absorption, free-free absorption or a combination of these and power-law components.

III.2.1 Synchrotron self-absorption (SSA)

At low radio frequencies the intensity of the synchrotron radiation may become sufficiently high (optically thick regime) for re-absorption, termed synchrotron self-absorption, to take over. The process may be important, or even dominant, for compact sources Slysh (1990); Chevalier (1998). We model the synchrotron radio spectra that may turnover due to self-absorption at low radio frequencies as (e.g. Tingay & de Kool, 2003)

[TABLE]

[TABLE]

where $\nu_{\tau=1}$ is a frequency at which the optical depth ($\tau$) reaches unity.

III.2.2 Free-free absorption (FFA)

The self-absorbed bremstrahlung (i.e. free-free absorbed) radio spectrum can be expressed as (e.g. McDonald et al., 2002)

[TABLE]

where the opacity coefficient is given by

[TABLE]

As discussed in §I radio emission from NGC 253 is a mixture of synchrotron emitting cosmic rays from SNR and thermal emission Hii regions. The free-free absorption is expected to start dominating at low radio frequencies, where the intensity of the electrons in the ionized gas becomes high (optically thick regime). For NGC 253 it is a natural assumption that the thermal plasma co-exists with the synchrotron emitting electrons, hence the radio spectrum can be modeled as a synchrotron power-law with an internal free-free absorbing screen (SFA; Tingay & de Kool, 2003),

[TABLE]

III.3. Weighted non-linear least squares fitting

All measurements in this paper are considered independent of each other (in the GLEAM survey valid for flux densities $\gtrsim 5$ Jy; see Hurley-Walker et al. (2017), thus a simple form of $\chi^{2}$ statistic is used for the fitting of the radio spectra, which at the same time is the goodness-of-fit of the fitted model ($M_{i}$),

[TABLE]

for $i=1,..,n$ data points. In minimization of Eqn. 8 we use the Levenberg-Marquardt algorithm Levenberg (1944); Marquardt (1963) implemented in the Python333http://www.python.org module lmfit Newville et al. (2014).

III.4. Model selection

For the formal model selection we use the Bayesian inference method. We follow the prescription outlined in Callingham et al. (2015), with the log-likelihood function (the probability of observing the data given model parameters $\theta$) in the form of

[TABLE]

Under the hypothesis that the models being compared ($M_{1}$, $M_{2}$) are equally likely, the model selection can be performed based solely on the Bayesian evidence ($Z$), where

[TABLE]

and

[TABLE]

The dimensionality of the integration depends on the number of model parameters. If $\Delta\text{ln}(Z)\geqslant 3$ model $M_{2}$ is strongly favoured over $M_{1}$. If $1<\Delta\text{ln}(Z)<3$ model $M_{2}$ is only moderately favoured over $M_{1}$, and if $\Delta\text{ln}(Z)<1$ the preference of one model over the other is inconclusive. For more discussion on the theoretical background of the method used see Callingham et al. (2015). We use the MultiNest tool Feroz & Hobson (2008); Feroz et al. (2009, 2013); Buchner et al. (2014) for our calculations of the Bayesian evidence.

IV. Results

In what follows we refer to the ‘halo’ as the radio emission beyond the boundary of the optical disk of the galaxy. Given the low angular resolution of the MWA observations we distinguish only between the galaxy disk and the extended synchrotron halo.

IV.1. Background sources

There are three discrete radio sources located within the extended emission of NGC 253; the sources are marked in Figure 1, and their positions are based on the NVSS and Heesen et al. (2009a) measurements.

The discrete radio source no. 1 is located at RA(J2000)=$00^{h}47^{m}59^{s}.10$, Dec(J2000)=$-25^{\circ}18^{\prime}22^{\prime\prime}.45$ at 200 MHz, and is most likely a background AGN Carilli et al. (1992). The radio source is detected in the GLEAM images, but at low frequencies becomes increasingly confused with the NGC 253 halo emission. We measure the flux density of the source only in the deep GLEAM image (Table 4). The discrete radio source no. 2 is located at RA(J2000)=$00^{h}47^{m}12^{s}.01$, Dec(J2000)=$-25^{\circ}17^{\prime}43^{\prime\prime}.9$ and is most likely a faint background AGN Heesen et al. (2009a). This source is heavily embeded in the NGC 253 extended emission in the MWA images. The discrete radio source no. 3 is located at RA(J2000)=$00^{h}47^{m}44^{s}.91$, Dec(J2000)=$-25^{\circ}13^{\prime}38^{\prime\prime}.4$. This source is embeded in the extended emission in our MWA images, but is clearly detected in the TGSS ADR1 image (Figure 1). The flux density measurements of the background sources are listed in Table 4, and the spectral modeling results are given in Table 3. We subtract the estimated flux density contribution of these sources from the total flux density measurements of NGC 253. In addition, we model the background source no. 1 as a point source with a peak flux density $240$ mJy at 169 MHz and $223$ mJy at 200 MHz, and for pictorial purposes we subtract it directly from the radio image plane. The resulting radio contours are overlaid on H$\alpha$ and X-ray images in Figure 2 and discussed in §V.1.

IV.2. NGC 253

Radio images of NGC 253 at six chosen radio frequencies are presented in Figure 4. At 200 MHz, the deepest image from the presented here GLEAM observations, the size of NGC 253 is 1310 arcsec (major axis) and 535 arcsec (minor axis) measured at a $\text{PA}=52^{\circ}$, with a total radio luminosity density $2.4(\pm 0.1)\times 10^{22}$ W Hz*-1*. At 169 MHz (MWA/EoR0 image), the size increases by 3–6 percent, to 1440 arcsec (major axis) and 615 arcsec (minor axis), which may be a combination of the intrinsic increase in size and the uncertainty of the measurement.

IV.2.1 Total radio emission

Radio continuum spectra of NGC 253 between 76 MHz and 10.7 GHz are plotted in Figure 4 and the flux density measurements are tabulated in Table 4. Background radio sources (§IV.1) located within the diffuse emission of NGC 253 were subtracted from the total flux density measurements. In the construction of the radio spectrum we used archival data provided the measurements were of angular resolution comparable to GLEAM or were sensitive to low brightness emission on angular scales of at least 0.5 deg, the total flux density was integrated over the diffuse emission of NGC 253 and not fitted by Gaussian components, and the absolute flux density scale and the uncertainties of the measurements were quoted. We do not use measurements at angular resolution of $>20$ arcmin because of confusion of NGC 253 with nearby sources.

We find the best fitting model to be a 2nd-order polynomial with $S_{0}=7.30\pm 0.04$ Jy, $\alpha=0.56\pm 0.01$ and a curvature $c_{1}=-0.12\pm 0.01$ at a reference frequency of 1 GHz ($\chi^{2}=140$, with degrees of freedom: dof $=45$; Figure 4), which is significantly preferred to a simple power-law ($\Delta\text{ln}(Z)=100.5\pm 0.3$).

IV.2.2 Central starburst region

The angular resolution of the MWA data is too low to resolve the central starburst region of NGC 253; the highest angular resolution achieved is 102 arcsec at 227 MHz (GLEAM) and 138 arcsec at 169 MHz (EoR), which is over three times the size of the NGC 253 starburst region Condon et al. (1982); Antonucci & Ulvestad (1988); Carilli (1996). We construct radio spectra of the central starburst region using data from the literature and new measurements from the TGSS ADR1 survey (Figure 4, Table 4). We limit the measurements to those that are at an angular resolution comparable to the size of the central starburst region (approximately 20–30 arcsec).

We find the best fitting model to be a 2nd order polynomial with $S_{0}=2.28\pm 0.02$ Jy, $\alpha=0.20\pm 0.01$ and a curvature $c_{1}=-0.24\pm 0.01$ at a reference frequency of 1 GHz ($\chi^{2}=12.8$, dof $=13$; Figure 4). We further attempt to model the spectral turnover, and we find SFA (Eqn. 7) to be the best fitting model with $S_{\tau=1}=4.43\pm 0.14$ Jy, $\alpha=0.43\pm 0.01$ and $\nu_{\tau=1}=238\pm 15$ MHz ($\chi^{2}=42.9$, dof $=13$). Based on the Bayesian evidence the SFA model (synchrotron plasma absorbed by an internal free-free absorbing screen) is preferred to the pure synchrotron self-absorption model, SSA ($\Delta\text{ln}(Z)=8.3\pm 0.3$).

IV.2.3 Spectral index maps

Using the total flux density images at 200 MHz (GLEAM), 169 MHz (MWA/EoR) and 1.46 GHz Carilli et al. (1992), we created spectral index distribution maps and their corresponding uncertainty and signal-to-noise ratio maps (Figure 6). To create the spectral index maps we convolved the 1.46 GHz image to the resolution of the GLEAM image (for the 200 MHz–1.46 GHz spectral index map), and the MWA/EoR image (169 MHz–1.46 GHz spectral index map). The spectral index maps suggest an apparent variation in $\alpha$ across the NGC 253 disc and halo. The central regions of the galaxy are dominated by a flat component ($\alpha=0.31-0.34$), coinciding with the central starburst. The gradual steepening along the minor axis seen by Heesen et al. (2009a) is seen only on the northern side of the galaxy in our maps. Further, the region extending SW from the central starburst region seems to be flatter ($\alpha\lesssim 0.53$) than in the other parts of the galaxy outside the central starburst ($\alpha\sim 0.60-0.65$).

To verify the significance of the spectral index variation across the galaxy, we used the T-T method Turtle et al. (1962). The method allows one to estimate a spectral index within defined regions of a source between two frequencies. We define nine regions within NGC 253 (Figure 6). Due to the low angular resolution of our observations and to avoid oversampling, only one data point is associated with each region. The results are shown in Table 6 and Figure 6.

We find that between 200 MHz and 1.465 GHz the apparently flat regions (Region 8 and 9) are statistically different from the other regions within NGC 253 apart from Region 1 (Figure 6, Table 6). This spectral flattening is not present in the radio spectral index distribution map between 330 MHz and 1.46 GHz of Heesen et al. (2009a), even though the same high frequency map is used. This clearly indicates that the flattening occurs at $<300$ MHz. There is also a slight flattening of the spectral index in the NE region perpendicular to the major axis (Region 1), although we find this flattening to be statistically different only from Region 2 (eastern NW halo), 4 (radio spur) and 5 (including central starburst) in the $\alpha^{\rm 200MHz}_{\rm 1.4GHz}$ map. All regions further flatten at 169 MHz, reducing the differences between spectral indices of the regions.

V. Discussion

V.1. Low-frequency synchrotron radio halo

A large-scale radio halo in NGC 253 was discovered and confirmed by Carilli et al. (1992) and extensively studied by Beck et al. (1994), Heesen et al. (2009a, b) and Heesen et al. (2011). This synchrotron halo is most pronounced at low radio frequencies, with the estimated scale heights of $1.7\pm 0.1$ kpc at 1.4 GHz and $2.5\pm 0.2$ kpc at 330 MHz Heesen et al. (2009a). Both the deep 200 MHz GLEAM image and the MWA/EoR image at 169 MHz reveal the extended synchrotron halo, which is at least as extensive as one detected in the Carilli et al. (1992) 330 MHz map (Figure 1).

V.1.1 Maximum vertical extent

We measured the observed, projected maximum vertical extent of the NGC 253 disk and halo at 169 MHz as a function of the distance from the nuclear region along the major axis (Figure 7, Table 7). The extent is measured perpendicular to the major axis (at $\text{PA}=-38^{\circ}$, z-direction) in steps of 132 arcsec (2.4 kpc) separately for the North (filled circles) and South (empty circles) side of the disk and halo as divided by the major axis ($\text{PA}=52^{\circ}$). We find the projected radio halo to extend up to 4.75 kpc above the optical ($B$ band including 90 per cent total light, Lauberts & Valentijn, 1989) and 6.3 kpc above the infrared (total $K_{\rm s}$ band, Jarrett et al., 2003) edge of the galaxy, reaching up to a total 7.9 kpc in $z$-direction. This is consistent with previous radio measurements at higher radio frequencies (Reynolds & Harnett, 1983; Hummel et al., 1984; Carilli et al., 1992, Figure 1; but cf. Heesen et al. 2009a), as well as broadband X-ray observations of the galaxy’s extended extraplanar emission Pietsch et al. (2000). Heesen et al. (2009a, b) attributed decrease of scale height they measured and modeled to the increased synchrotron losses in the central regions, where the magnetic field is highest.

V.1.2 Halo morphology

The shape of the radio halo in our MWA radio images (Figure 2) resembles the ‘horn-like’ or ‘X-shaped’ structure seen at GHz radio frequencies Heesen et al. (2009a), in Hi Boomsma et al. (2005); Lucero et al. (2015), X-rays Fabbiano (1988); Pietsch et al. (2000); Bauer et al. (2008), H$\alpha$ (G. Meurer, priv.comm.; Figure 2), UV Hoopes et al. (2005) and far-IR Kaneda et al. (2009).

The radio halo was investigated in detail by Heesen et al. (2009b) who, through modeling of the large-scale magnetic field, attributed its origin to disk wind, confirming previous suggestions Carilli et al. (1992); Beck et al. (1994). This is also in line with the H$\alpha$ and optical analyses of the inner starburst-driven superwind Westmoquette et al. (2011). In our MWA maps both the NE and NW halo regions are pronounced. The extended soft X-ray emission ($<1$ keV) of the halo, detected in the north-western direction from the NGC 253 disk, is interpreted as bubbles of hot low density gas (Figure 2; Pietsch et al., 2000; Strickland et al., 2002). Heesen et al. (2009b) postulates that the large-scale magnetic field of the halo follows the walls of these bubbles, where it may be compressed, producing the X-shaped synchrotron radio halo as well as heating up pre-existing cold gas to X-ray energies. The northern halo can also be easily seen, in projection, in our Figure 2 where we overlay MWA/EoR intensity contours on the soft X-ray emission XMM-Newton image.

The SE region of the extended halo, the ‘spur’ Carilli et al. (1992), is contaminated by a background source. We modeled the background source as unresolved (at MWA angular resolution) and subtracted it from the 169 MHz EoR and deep 200 MHz GLEAM images as described in §IV.1. The residual emission, which we consider intrinsic to the ‘spur’ is shown in Figure 2 overplotted on the H$\alpha$ and X-ray images. Although slightly offset (as expected, cf. Fig 19 in Heesen et al., 2009b), the feature broadly coincides with the extended outflows at both frequencies as clearly seen in our figure. The spur has been previously interpreted as originating from the active star formation region in the NE end of the bar Waller et al. (1988); Carilli et al. (1992). We confirm the association of the radio ‘spur’ with an extended H$\alpha$ outflow, which can be clearly seen in our Figure 2 (G. Meurer, priv.comm.; Kennicutt et al., 2003). In our radio spectral map between 200 MHz and 1.46 GHz the spur displays the steepest spectral index within the galaxy (Figure 6), gradually steepening from $\alpha\sim 0.7$ to $\alpha\sim 0.9$ outwards from the galaxy disk (though it still contains the background source no.1 with $\alpha=0.57$), which is indicative of ageing unabsorbed synchrotron plasma.

V.2. Spectral properties of NGC 253

V.2.1 Broadband spectrum of total radio emission

The broadband spectrum of the total radio emission from NGC 253 is steep, although flattening at MHz radio frequencies (Figure 4). The total radio emission originates from SNRs, Hii regions (predominantly central starburst region; e.g. Ulvestad & Antonucci 1997, Ulvestad 2000, but cf. Waller et al. 1988, Hoopes et al. 1996) and electrons (cosmic rays) freely spiralling in the large scale magnetic field (radio halo; Heesen et al., 2009a). The steep spectrum is understood to be of a synchrotron origin. The flattening of the spectrum, however, may be due to a number of reasons, including some degree of absorption of the synchrotron emission (Eqn. 3) and a low energy cut-off of the electron population, where in general it is assumed the electron energy ($E$) spectrum can be described as a power-law $N(E)\propto E^{-p}$ with the index $p$ related to the radio spectral index as $\alpha=(p-1)/2$. The low radio frequency flattening of the spectra of starburst galaxies is not unusual and it has been observed previously (e.g. Marvil et al., 2015).

We attempt to separate the contribution of the extended and central starburst emission to the total flux densities to see if the spectral flattening can be attributed mostly to the absorption occurring in the central starburst region. To do this we simultaneously fitted the total flux density and the central starburst region spectra, assuming an underlying 2- or 3-component spectrum composed of the central starburst region (component C) and extended emission (component E, composed of one or two components). As an example of our fitting method, in the case where the component E is modeled as a power-law and component C as a free-free absorbed synchrotron emission, the model equation ($S_{\nu}^{\rm mod})$ takes the following form

[TABLE]

where

[TABLE]

[TABLE]

$\tau_{\rm ff}(\nu)$ is given by Eqn. 6 and the indices indicate: $\rm tot$ – total, $\rm cor$ – core, $\rm C$ – component C, $\rm E$ – component E, $\rm ext$ – extended. This model equation is then compared to our observed data ($S_{\nu}^{\rm obs}$), where $S_{\nu}^{\rm obs}=S_{\nu}^{\rm obs,tot}+S_{\nu}^{\rm obs,cor}$. In the 2-component model, we fit the component E with either simple power-law or a curved spectrum (2nd degree polynomial; Eqn. 1). In the 3-component model, we fit the component E with a combination of power-law and either a curved spectrum, synchrotron self-absorbed (Eqn. 3) or synchrotron free-free absorbed component (Eqn. 7). The component C is modeled with either 2nd-degree polynomial, self-absorbed synchrotron or synchrotron power-law emission with a free-free absorbing screen.

We find the best fitting model to be the 3-component model, with (1) the component E modeled as a combination of a simple power-law and 2nd order polynomial, with a total flux density $S_{0}=5.08\pm 0.50$ Jy and $\alpha=0.71\pm 0.01$ at reference frequency 1 GHz, and 34% of $S_{0}$ becoming absorbed at low frequencies as described by 2nd order polynomial with $c_{1}=-0.76\pm 0.17$, and (2) the component C modeled as a synchrotron plasma with an internal free-free absorbing screen, with $S_{\tau=1}=4.34\pm 0.11$ Jy, $\nu_{\tau=1}=231\pm 14$ MHz, $\alpha_{\rm SSA}=0.41\pm 0.01$. The 3-component model ($\chi^{2}=173$, dof $=58$) is favoured over any 2-component model, even if the component E is modeled as a 2nd order polynomial ($\Delta\text{ln}(Z)>4.3\pm 0.3$).

Preference of the 3-component model indicates that flattening of the component E spectrum at the lower frequencies is non-negligible. In principle, this flattening could be attributed to synchrotron self-absorption caused by shock re-acceleration of the halo/disk plasma, an external free-free absorbing screen, or an intrinsic low-energy cut-off of the electron distribution. Thermal free-free absorption can be largely excluded based on the limited evidence for high thermal content in the NGC 253 halo, especially in the SW region (see §V.2.2 for details). Although our Bayesian inference tests indicate that the flattening caused by the synchrotron self-absorption is moderately favoured over the low-energy cut-off in the electron distribution which could be inferred from the 2nd order polynomial fit ($\Delta\text{ln}(Z)=2.5\pm 0.3$), we find that the former fit is associated with very high uncertainties. We also find that any model invoking multiple internal components that we tested is strongly favoured over an external free-free absorbing screen. Our results do not change in the absence of the data points between 300 and 600 MHz that may seem unusually high, which further strengthens our result that the radio emitting plasma in the disk and halo of NGC 253 is composed of at least two spectral components that behave differently. This result is also in line with our findings on the radio spectral index variation across the galaxy (discussed further in the next section).

Furthermore, another important result of our spectral modeling is that the central starburst region is best modeled by the SFA model. Although, in principle, the 2nd order polynomial is statistically favoured, the SFA model is more realistic. Curved radio spectra can be explained by low energy cut-off of electron population, SSA, SFA or FFA models. As we have shown the SSA model is statistically ruled out (§IV.2.2). Given the overwhelming evidence of significant thermal component in the central starburst (e.g. Figure 2; Ulvestad & Antonucci, 1997; Keto et al., 1999; Kepley et al., 2011) coexisting with synchrotron plasma, the synchrotron free-free absorption model is more likely than the low energy cut-off in electron population. The plasma becomes optically thick around frequency $230$ MHz. Given this result, and under a simplified assumption that a uniform optical depth holds across the region, we estimate a typical emission measure Oster (1961); Spitzer (1978) of the absorbing gas towards the central starburst region to be very high, of the order $4\times 10^{5}$ pc cm*-6* (assuming electron temperature $T_{e}=4000$ K; cf. discussion in Carilli, 1996).

V.2.2 Radio spectral index distribution maps

We now consider the origin of the $\alpha$ variation across the NGC 253 disk and halo. The southern flattening occurs beyond the SW spiral arm, in the halo region. In Figure 8 we overplot the $\alpha^{\rm 200MHz}_{\rm 1.4GHz}$ and extraplanar Hi contours on the XMM-Newton soft X-ray image of NGC 253. The diffuse X-ray emission indicates ionized hot gas. As pointed out by Lucero et al. (2015), the neutral cold Hi gas seems to surround the X-ray emitting regions. The radio spectral index spatial variations seem to follow the distribution of the X-ray emission, with the $\alpha$ steepening occurring in the regions of intense soft X-ray emission (radio spur and NW halo) and the flattening around the voids of diffuse X-ray plasma (western SE halo, eastern NW halo). This distribution seems to also match the extraplanar Hi emission, especially in the western SE halo region. In H$\alpha$ we detect faint diffuse emission in the NE halo and the southern ‘spur’ (Figure 2), with the line fluxes measured down to $3\times 10^{-18}$ erg s*-2* s*-1* arcsec*-2*. In these regions the spectral index seems to steepen (Figure 8B). There is, however, almost no H$\alpha$ emission present in the western SE halo, while most of the flattening of the component E in the modeling of §V.2.1 is due to this region (Figure 8C, regions 8 and 9 in Figure 6).

It has previously been suggested that the halo gas originates from both galactic ‘fountains’ from the star-forming disk and a galactic superwind Pietsch et al. (2000); Bauer et al. (2008). This strong superwind may be pushing and collimating the neutral cold gas in the halo Lucero et al. (2015). In the case of strong collimation shocks one may observe flattening of radio spectra due to synchrotron self-absorption in transverse shocks. Our results seem to favour such a scenario. It is also worth noting that the spectrum flattening of the SW halo corresponds to an extended loop, or arch, seen in optical images (blue filter; Beck et al., 1982). However, if the SW halo is predominantly diffuse, and of low density, the flattening may be rather due to an intrinsic low-energy cut off of the electron distribution.

Another important note is that, based on our broadband radio spectrum modeling, the SW halo region cannot be fully responsible for the total spectrum flattening. Radio emission that becomes absorbed at lower frequencies constitutes more than 30% of the total extended radio emission at 1 GHz ($1.73\pm 0.36$ Jy), while the SW region is only 0.77 Jy at that frequency, which means that the flattening must be also occurring, although in a smaller degree, in other regions across the disk and halo.

Although the SW flattening of the radio spectral index is most likely of a synchrotron origin, an external free-free absorbing screen was also previously suggested. The foreground absorption model was favoured by Bauer et al. (2008) based on their X-ray data modeling and apparent differences of radio and X-ray halo morphologies. As proved later (e.g. Heesen et al., 2009a, b, this publication, §V.1), deep continuum and polarization radio observations at both GHz and MHz frequencies reveal the horn-like structure of the radio halo, which directly resembles the X-ray diffuse emission. Based on the equipartition assumptions, Heesen et al. (2009a, 2011) find the magnetic field within the halo to be very high, $7-12\mu$G, reaching as much as $160\pm 20\mu$G in the central regions and $46\pm 10\mu$G in the starburst outflow. The magnetic field in the central regions is strong enough for synchrotron emission to contribute a few per cent to the total X-ray emission Lacki & Thompson (2013). As discussed in the previous section, we also find that an external free-free absorbing screen is not a statistically preferred model. These new findings support models in which the total X-ray emission may indeed come from a combination of thermal and synchrotron plasma rather than multi-temperature pure thermal plasma with an externally caused absorption (cf. Bauer et al., 2008).

VI. Conclusions

We present deep, low-frequency radio continuum images and flux density measurements of a nearby, archetypal starburst galaxy, NGC 253. Our data are part of the Galactic and Extragalactic All-Sky MWA Survey and the MWA EoR observations. The images span frequencies between 76 and 231 MHz at angular resolution of 1.7 – 5 arcmin and rms noise levels of 4 – 75 mJy (depending on frequency), and present the deepest measurements of NGC 253 at these low radio frequencies yet.

Our main findings are summarized as follows.

We detect a large-scale synchrotron radio halo that at 154–231 MHz displays the X-shaped/horn-like structure seen at GHz radio frequencies, and is broadly consistent with other multiwavelength observations of NGC 253. 2. 2.

The projected maximum vertical extent of the synchrotron emission at 169 MHz extends up to 7.5 kpc NW (7.9 kpc SE) from the major axis of NGC 253, consistent with large-scale soft X-ray emission (extending 9 kpc NW) and X-ray outflow (6.3 kpc SE). 3. 3.

The radio spectrum of the central starburst region of NGC 253 is significantly curved at low radio frequencies, with the spectral turnover occurring around 230–240 MHz, which is for the first time statistically constrained. 4. 4.

The radio spectral index maps show significant spectral variations in the structure of NGC 253 between 200 MHz and 1.465 GHz. In particular, we isolate a region of statistically significant spectral flattening to the western side of the SE halo. However, as the SW region is rather faint at 1.46 GHz it cannot be fully responsible for the total spectrum flattening, which indicates that the flattening must be also occurring, likely in a smaller degree, in other regions across the disk and halo. 5. 5.

The broadband spectrum of integrated total radio emission of NGC 253 is best described as a sum of central starburst and extended emission, where the central starburst component is best modeled as an internally free-free absorbed synchrotron plasma, and the extended emission as synchrotron emission flattening at low radio frequencies. We also find that an external free-free absorbing screen is not a statistically preferred model when compared to models including multiple internal components. 6. 6.

We find that the extended emission of NGC 253 is best modeled by a combination of two synchrotron components, one of which becomes significantly absorbed at low radio frequencies. The flattening occurs at frequencies below $300$ MHz, and may be attributed to synchrotron self-absorption of shock re-accelerated electrons or an intrinsic low-energy cut off of the electron distribution.

Acknowledgments

ADK thanks P. A. Curran for valuable discussions on data modeling and constant encouragement in achieving the goals. The authors thank the anonymous referee for careful reading of the manuscript and suggestions that improved this paper. The authors thank W. Pietsch and D. Lucero for providing, respectively, X-ray and Hi fits images of NGC 253, and O.I. Wong and X. Sun for helpful comments. The authors thank V. Heesen for 1.465 GHz image of NGC 253 and for helpful discussions. SB acknowledges funding for the ICRAR Summer Scholarship.

This research was conducted under financial support of the Australian Research Council Centre of Excellence for All-sky Astrophysics (CAASTRO), through project number CE110001020. This scientific work makes use of the Murchison Radio-astronomy Observatory, operated by CSIRO. We acknowledge the Wajarri Yamatji people as the traditional owners of the Observatory site. Support for the operation of the MWA is provided by the Australian Government (NCRIS), under a contract to Curtin University administered by Astronomy Australia Limited. We acknowledge the Pawsey Supercomputing Centre which is supported by the Western Australian and Australian Governments. This research has made use of the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. This publication uses the following radio data reduction software: the Multichanel Image Reconstruction, Image Analysis and Display software (Miriad; Sault et al., 1995), the Common Astronomy Software Applications package (CASA; McMullin et al., 2007) and the Astronomical Image Processing System AIPS. AIPS is produced and maintained by the National Radio Astronomy Observatory, a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc.