Smooth HI Low Column Density Outskirts In Nearby Galaxies
Roger Ianjamasimanana Fabian Walter W. J. G. de Blok George H. Heald Elias Brinks


 
  1Max-Planck Institut für Astronomie, 
  Königstuhl 17, 69117, Heidelberg, Germany
    2Netherlands Institute for
  Radio Astronomy (ASTRON), Postbus 2, 7990 AA Dwingeloo, The  Netherlands
    3Astrophysics, Cosmology and Gravity Centre, 
  Department of Astronomy, University of Cape
  Town, Private Bag X3, Rondebosch 7701, South Africa  
    4Kapteyn Astronomical Institute, 
    University of Groningen, PO Box 800, 9700 AV, Groningen, The Netherlands
    5CSIRO Astronomy and Space Science, 26 Dick Perry Avenue, Kensington WA 6151, Australia
    6Centre for Astrophysics Research, University of Hertfordshire, College Lane, Hatfield, AL10 9AB, UK
   ianjamasimanana26@gmail.com 
   walter@mpia.de
   blok@astron.nl
   George.Heald@csiro.au
   e.brinks@herts.ac.uk
for Publication in the Astronomical Journal

    The low column density gas at the outskirts of galaxies as
    traced by the 21 cm hydrogen line emission (H i) 
    represents the interface between galaxies and the intergalactic medium, i.e., 
    where galaxies are believed to get their supply of gas to
    fuel future episodes of star formation. Photoionization models predict a 
    break in the radial profiles of H i at a column density of ∼ 5 × 10^19 cm^-2 
    due to the lack of self-shielding against extragalactic ionizing photons. To investigate the 
    prevalence of such breaks in galactic disks and to characterize what
    determines the potential edge of the H i disks, we study the 
    azimuthally-averaged H i column density profiles of 17 nearby galaxies from The H i
    Nearby Galaxy Survey (THINGS) and supplemented in two cases with published Hydrogen Accretion in
    LOcal GAlaxieS (HALOGAS) data. To detect potential faint H i emission that would
    otherwise be undetected using conventional moment map analysis, we line up individual profiles to the same 
    reference velocity and average them azimuthally to derive stacked radial profiles. 
    To do so, we use model velocity fields created from a simple 
    extrapolation of the rotation curves to align the profiles in velocity at radii beyond 
    the extent probed with the sensitivity of traditional 
    integrated H i maps. With this method, 
    we improve our sensitivity to outer-disk H i emission by up to an order of magnitude. 
    Except for a few disturbed galaxies, none show evidence for a sudden change in the slope of the H i radial
    profiles, the alleged signature of ionization by the extragalactic background.







§ INTRODUCTION
 
The 21 cm atomic hydrogen (H i) line can trace the interstellar medium (ISM) out 
to large radii in galaxies <cit.>, 
typically well beyond the optical disk. 
Detailed studies of the shapes of the radial H i  
distribution are important to understand the physical conditions of the 
outer disk of galaxies, e.g., constrain models of the 
H i to H ii transition, to trace the influence of 
environment on H i, and to use the H i gas as a tracer of the (dark) 
matter distribution.
In general, the radial H i distribution 
tends to be approximately flat out to about the edge of the 
optical disk (with a central depression for many spiral galaxies), 
followed by a gradual decline 
<cit.> at a column density of 
about a few times  10^20 cm^-2. A second break in the radial H i 
profiles has been reported using high sensitivity observations 
of M33 <cit.> and NGC 3198 <cit.>. 
In these studies, a break was reported 
at a column density of a few times 10^19 cm^-2. 
Such a behavior has been interpreted in the context of the 
existence of a critical column density below which H i is ionized by the 
extragalactic radiation field <cit.>. 
However, as mentioned by <cit.>, 
the degree of steepness of the radial H i fall-off predicted by 
ionization models depends on many factors such 
as the vertical structure of the H i disk, the H i 
mass distribution, and more importantly the flux of the 
ionizing photons, which are not well 
constrained by observations <cit.>.  

Although there have been many discussions in the literature 
explaining the edge of the H i disk from theoretical 
viewpoints <cit.>, 
detailed observations of the radial extent of the outer H i disk in 
a sample of galaxies are lacking. 

Using single-dish observations with the Arecibo telescope, 
<cit.> examined the outer radial H i 
distribution of the northern part of the major axis of M33 
down to a (3σ detection limit) column 
density of ∼ 2 × 10^18 cm^-2. The column density is 
∼ 2 × 10^21 cm^-2 in the central disk and the value drops 
below ∼ 2 × 10^20 cm^-2 just outside the optical disk, which then 
gradually declined to ∼ 3 × 10^19 cm^-2. 
A second break in the radial 
profile was observed at ∼ 3 × 10^19 cm^-2, below which the value 
decreased to ∼ 2 × 10^18 cm^-2 within ∼ 1 kpc. 
The NRAO[The National Radio Astronomy 
Observatory is a facility of the National Science 
Foundation operated under cooperative agreement by Associated Universities, Inc.] 
Very Large Array (VLA) observation of NGC 3198 by <cit.>, with a 
3σ column density limit of 4 × 10^18 cm^-2, 
also revealed the existence of 
a radial break of the H i disk below ∼ 2 × 10^19 cm^-2. 
The VLA data showed that the H i surface density of NGC 3198 decreased 
by an order of magnitude from ∼ 30 to ∼ 33 kpc. 
The abrupt decline in NGC 3198 was, however, not resolved within 
the synthesized beam of the observation (∼ 2.7 kpc).
Recently <cit.> studied the H i edge in M83 using interferometric 
observations with the Karoo Array Telescope <cit.>. 
These observations reached 
a 3σ column density sensitivity 
of 5.6×10^18 cm^-2. Based on the 
analysis of the H i column density map 
of M83, they concluded that this galaxy also 
has an H i 
edge similar to what was found for 
NGC 3198 and M33. <cit.> considered 
ram pressure or ionization by the intergalactic 
medium as a possible cause of the H i edge. 
There are also other observations targeting the H i outskirts  
of nearby galaxies <cit.> but 
we restrict ourselves to those that specifically address 
the edge of the H i disk.                          

In this paper, we study the shapes of the radial H i 
column density profiles for a large 
sample of nearby galaxies and investigate whether the presence of 
an edge is a common feature of the H i disk. 
We do this by employing a new method to extract faint H i emission at 
large galactocentric radii from H i data cubes. 
Our paper is organized as follows. 
In Section <ref>, we present our data and 
sample galaxies. In Section <ref>, 
we describe our new stacking approach to derive the 
radial H i column density profiles. 
In Section <ref>, we discuss our results. 
In Section <ref>, we present our summary and conclusions. 


§ DATA AND SAMPLE

Analyzing the distribution of the H i emission at 
large radii requires high sensitivity observations. 
With the limitation of current state-of-the-art telescopes, 
only observations of nearby galaxies 
provide enough spatial resolution and sensitivity to resolve the H i edge 
and investigate the H i distribution at large radii. 
Based on the model by <cit.> and the previous study of NGC 3198 by <cit.>, 
a column density sensitivity limit of ∼10^19 cm^-2 and a spatial resolution better than 2.7 kpc 
are required to detect and resolve the H i edge. 
Using a stacking method that will be described later, we satisfy this criteria 
using data from The H i Nearby Galaxy Survey <cit.>. 
THINGS provides VLA 21-cm line maps of 34 nearby 
spiral and dwarf galaxies at an average angular 
resolution of ∼11 (linear resolution is given in Table <ref>). 
l c c c c c

0pt
The Sample galaxies 


1cGalaxy     incl     Dist.     
1cB_min   B_maj     Lin. res.
 
    []   [Mpc]    []   []   [pc]

1c1   1c2   1c3
    4     5     1c6

IC 2574     53     4.0     11.9     12.8     240 

NGC 3621     65     6.6     10.2     15.9     419 

NGC 2366     64     3.4     11.8     13.1     206 

DDO 154     66     4.3     12.6     14.1     278 

NGC 6946     33     5.9     5.6     6.0     167 

NGC 4736     41     4.7     9.1     10.2     220 

NGC 3521     73     10.7     11.2     14.1     656 

NGC 4826     65     7.5     9.3     12.2     391 

NGC 3627     62     9.3     8.8     10.6     438 

NGC 925     66     9.2     5.7     5.9     260 

NGC 3198     72     9.4     11.6     13.0     560 

NGC 7331     76     14.7     5.6     6.1     418 

NGC 2841     74     14.1     9.4     11.1     698 

NGC 2403     63     3.2     7.7     8.8     127 

NGC 2903     65     8.9     13.3     15.3     617 

NGC 7793     50     3.9     10.8     15.6     250 

NGC 2976     65     3.6     6.4     7.4     121 
Column 1: name of galaxy; 
Column 2: inclination; Column 3: distance; Column 4&5: 
minor and major axes of synthesized beam in arcsec using a natural weighting; 
Column 6: average linear resolution.




In this work, we use the high resolution rotation curves derived by <cit.> 
(see Section <ref>). Thus we require that our galaxies are part of the 
<cit.> sample. However, we exclude NGC 5055 as this galaxy is strongly warped. We also exclude 
M81 as its H i emission is heavily disturbed, due to its interaction with M82 and NGC 3077. 
Thus, our sample is reduced to 17 spiral and dwarf galaxies. 
Their observational properties and names are 
shown for reference in Table <ref> <cit.>. 

In this paper, we use the THINGS natural-weighted data and do not 
apply any residual flux correction nor blanking so that 
we preserve the noise properties of the data.  
However, to remove the flux discrepancy due to 
the mismatch between the shapes of the dirty beam and the restoring clean beam 
<cit.>, 
we clean the data close to the noise level 
(1.5σ, as opposed to the 2.5σ of the standard THINGS data release).
When cleaning to this depth, we find better agreement between 
the flux derived from the natural non-residual scaled cubes and the flux from the 
standard THINGS moment zero maps, which were derived 
from residual-scaled cubes <cit.>. 
For completeness, we also present in this paper the 
radial H i column density profiles 
derived from the standard THINGS moment zero maps.

A few of the THINGS galaxies have also been observed as 
part of the Hydrogen Accretion in LOcal GAlaxieS Survey 
<cit.>. For two galaxies (NGC 3198 and NGC 2403) 
which overlap with the HALOGAS sample and 
fulfill our sample selection criterion, we compare our results from the THINGS data 
with those from the HALOGAS data. HALOGAS is a deep H i survey of 
24 nearby galaxies, with 23 galaxies observed by the 
WSRT[Westerbork Synthesis Radio Telescope] 
<cit.> and one 
by the VLA <cit.>. The WSRT data by <cit.> have a 
typical (5σ per 4.1 km s^-1 channel) 
column density sensitivity limit of ∼5 
× 10^18 cm^-2 at 30 resolution using 
a Brigg's <cit.> robustness parameter of 0 with a Gaussian taper. 
The VLA data (NGC 2403) by <cit.> reached a 
(5σ per 5.2 km s^-1 channel) column density 
sensitivity of ∼5× 10^19 cm^-2 
using a Robust (robust = 0.2) or Briggs weighting scheme 
at 15 (230 pc) angular resolution.  


§ METHODOLOGY


We aim to derive H i column density values as a function of 
radius, beyond the extent probed by conventional pixel-by-pixel moment map 
analysis. To do so, we divide each galaxy into concentric elliptical annuli, 
each of them with a width along the major axis 
similar to the size of the synthesized beam of the observations. 
For THINGS and the VLA data from <cit.>, 
we adopt a common width of 15; for the 
HALOGAS data of NGC 3198, we use a 30 width. Using a stacking technique that 
will be described in the next section, we derive the 
average H i column density within each annulus to derive the radial profiles. 
Since we aim to recover 
H i flux out to large radii, correction for
sensitivity bias due to the response of the primary beam of the telescope is
important. We therefore apply a primary beam correction in the analysis that follows.



 §.§ Stacking of H i spectra


Detecting H i emission in the extreme outskirts of galaxies is 
challenging due to the low signal-to-noise (S/N) of the 
individual profiles.   
We here aim to push the sensitivity limit of the available H i data 
to much deeper levels. To do so, we stack and sum individual profiles within the 
elliptical annuli mentioned in the previous section to get high 
S/N azimuthally averaged H i profiles. 
Individual profiles at different positions within a data cube peak at different 
velocities predominantly due to the rotation 
of the gas within the galaxy disk. Consequently, 
they cannot be simply summed or averaged but individual 
spectra need to be aligned to a common 
reference velocity prior to the averaging. 
Faint H i 
emission that would otherwise be undetected 
(or “blanked”) using conventional pixel-by-pixel moment 
map analysis can thus be recovered. 
We therefore shift individual profiles 
by the amount suggested by the galaxy's velocity 
field to attempt to align them to the same velocity 
<cit.>. This technique is implemented in the GIPSY[The Groningen Image
Processing System, <https://www.astro.rug.nl/ gipsy/articles/gipsypaper.html>]
task SHUFFLE. The stacking is  
straightforward in the bright inner disk of galaxies, 
where the velocity fields are well determined 
<cit.>. However, in the faint outskirts of the 
disks, where the peak emission of the individual profiles falls below 3 
times the rms noise, direct 
measurements of the velocity field cannot be derived. 
At such large radii, we will use extrapolations 
of the rotation curves derived in the inner disk 
by <cit.> to create model velocity fields. 
As most galaxies reach the flat part of the rotation 
curve at large radii, we here simply assume that the rotation curve will 
continue to stay flat beyond the last point where the 
rotation curves can be directly measured. 
We also assume that the position angle and inclination remain constant as well, i.e., 
there is no warp beyond the outermost radius directly detected in H i emission. 
These assumptions are the key elements of our stacking procedure. 
As we will see below, 
stacked H i signal is recovered beyond the measured 
rotation curves, which implies that the assumption of flat, 
extended rotation curves holds. 




 §.§ Effects of noisy profiles
 
 
Since we are stacking profiles azimuthally, we are summing both low and high S/N profiles.
We thus assign a weight, w_i, to each individual 
profile or spectrum, S_i, before the stacking. 
We assume 
that the noise in the data cubes is uniform before the primary beam correction. 
This is a good assumption for our data cubes. 
After the primary beam correction, the noise is no longer uniform and we 
define our stacked profiles as

    S=∑_i=1^N w_iS_i∑_i=1^N w_i,

 where w_i=p_i^2σ_i^2 is the weight assigned 
 to spectrum S_i, with p_i being the associated primary beam 
correction factor (we adopt a Gaussian shape, with p_i = 1 at the center and 
0.5 at 15.4 and 17.0 from the pointing center of the VLA and WSRT, respectively) 
and σ_i is the rms noise before the correction. 
Since σ_i is assumed to be equal to σ for all 
profiles in a data cube, Equation <ref> becomes 

    S=∑_i=1^Np_i^2S_i∑_i=1^N p_i^2.
 
We use a simple propagation of 
errors to define the noise of the stacked profiles 
in each elliptical annulus. Some stacked profiles show 
negative bowls due to unsampled (u, v) data (i.e., “missing short spacings”). 
We correct for this by fitting a single or double Gaussian function with a 
polynomial background to the stacked profiles. 
We then convert the integrated flux of the fitted Gaussian 
components with respect to the fitted baseline to a column density. Finally, to correct for the orientation of the galaxies, 
we multiply the results by the cosine of the inclination of the galaxies to derive the actual surface density 
profiles. The adopted inclination values are listed in 
Table <ref>. 

   

 
 
 



§ RESULTS
 


 §.§ Examples of azimuthally stacked spectra
  
In Figure <ref> we show examples of the 
azimuthally averaged stacked profiles of NGC 7331, derived using the model velocity field 
obtained from the model rotation curve shown as the blue solid line in the figure. 
We show the profiles for all sample galaxies in Figure <ref> of the Appendix. 
In the region where the rotation curve is defined, we fit the rotation curve 
with a spline function. 
Beyond this, we assume a flat rotation curve. 
We derive extended model velocity fields from 
the model rotation curves and use them as inputs 
for the stacking method implemented in the 
GIPSY task SHUFFLE. 
Note that if our estimate of 
the rotation curves significantly deviates from the true rotation curves, 
we would not be able to detect a signal.  
Figure <ref> demonstrates that we are 
indeed able to detect signal beyond what is reached by 
conventional pixel-by-pixel analysis. The outermost point of the observed rotation curve 
measured by <cit.> is at ∼350 (using a 3 σ cut-off) 
but we are able to still 
detect stacked H i emission out to a radius of ∼530. 
The fairly constant width of the profiles indicates 
that our assumption of a flat rotation curve is reasonable. 
The profiles remain narrow and 
we do not see multiple peaked or broad profiles as would be expected from an incorrect 
shift of the individual velocity profiles. 
There is a shift away from `zero' velocities, though. This shift is the consequence of 
the northern side of NGC 7331 having slightly more H i emission 
than the southern side. 
 


 §.§ Radial profiles for all galaxies

In Figure <ref>, we show the radial column density profiles of the 
sample galaxies as derived using our stacking method. 
Out of the 17 galaxies studied, we reach a column density sensitivity limit of ∼ 10^19 cm^-2 
for 9 galaxies, and below ∼ 10^19 cm^-2 for 8 galaxies. 
As a check, we also show the radial column density profiles 
derived from the standard THINGS 
moment-zero maps in Figure <ref>. There is still a difference between the column densities 
derived from the moment-zero maps and from the stacking method despite the deeper cleaning of our data cubes. 
This small difference is due to the stacking of residual uncleaned emission in the data cubes, leading to 
somewhat higher fluxes. This effect can be corrected for in the moment maps, but less easily so for stacked profiles. 
We refer to <cit.> for an extensive discussion. The small offset, 
however, does not affect the conclusions of this paper. 
The moment-zero values can only be reliably determined above 3 times the local 
rms noise. Note that the rms noise in a moment zero map 
is not the same as the rms noise in 
a data cube channel map since a variable number of channels, N, contribute 
to any given pixel in the moment map. 
We therefore calculate the rms noise in the moment maps as the rms noise in one channel 
map multiplied by the square root of the number of the contributing channels. 
This assumes channels are independent, which is a valid 
assumption as the THINGS data are Hanning-smoothed <cit.>.
We indicate the (moment-zero) column density sensitivity corresponding to 3 times the 
rms noise in Figure <ref>. We also show in the same figure the 3 times 
rms noise column density sensitivity of the azimuthally stacked profiles. On average, 
the stacking method offers an order of magnitude 
sensitivity improvement compared to the traditional moment-zero analysis. 
Close inspection of Figure <ref> shows that the stacked 
profiles do not always reach the S/N = 3 limit. In these cases 
we are limited by the presence of low-level artefacts in the 
data (mostly due to deconvolution residuals) that prevent us 
from reaching lower limits. Nevertheless, the fact that with the 
current data set stacking already yields limits much lower than 
traditional moment-zero limits, clearly shows the promise of the method. 

The most surprising finding from Figure <ref> is that most galaxies show no evidence of a break. 
The profiles are mostly flat in the inner disk and then decline smoothly down to the sensitivity limit 
of our (stacked) data. 
Our assumption of a flat rotation curve 
 cannot be the reason for this absence of an H i break. 
 Any incorrect estimates of the rotation curves would 
 result in stacked profiles that are systematically 
 lower in amplitude but broader in width and this 
 would not remove the signature of a break. 
 A major velocity offset between real and assumed velocities 
 would result in non-detections, which would mimic an 
 artificially strong break in the profiles. 
 
 A possible cause of the absence of a break is that the H i 
 surface density distribution is azimuthally asymmetric. Our sample galaxies 
 exhibit this asymmetry. These include, among others, NGC 7331, NGC 3198, NGC 2903, NGC 3621. 
 We also note that the break in M33 observed by <cit.> is based on a single radial slice 
 along the major axis of the galaxy. In addition, the break in NGC 3198 observed by 
 <cit.> <cit.> is seen in only one side of the galaxy.   
 Thus the stacking method may smooth out any potential break as we are averaging 
 over the full azimuthal angle. To explore this possibility, we also divide the galaxies into 
 30 degree wide “pie-wedge-shaped” sectors and do the stacking inside annular sectors instead of the 
 full azimuthal annuli. The results are shown in Figure <ref> of the Appendix. Although the 
 radial profiles are steeper in some sectors than in others, the profiles are smooth without a clear break. 
 We thus conclude that the absence of a break in the (full) azimuthally averaged profiles are not caused 
 by the presence of an asymmetric H i distribution. 


 §.§ Exceptions

 A close inspection of the individual galaxies reveals that there are a 
 few exceptions to the conclusion that radial breaks 
 are not present in this sample. These are  discussed in the following.
 The radial profiles of NGC 4736 and NGC 3621 decline steeply 
 around 10^20 cm^-2 and for NGC 4826 at ∼5×10^19 cm^-2. 
 These three galaxies, however, have some structural peculiarities. 
 NGC 4736 has a ring-like H i morphology which manifests 
 itself as bumps in the H i column density 
 profiles shown in Figure <ref>. Earlier 
 lower resolution H i observations 
 by <cit.> also revealed this pattern. The authors 
 classified the galaxies as 
 consisting of 5 different zones in which a faint outer 
 ring surrounds the main body of the galaxy. 
 The H i edge around 7 kpc corresponds to the transition from the main disk of the 
 galaxy and the faint outer rings identified by <cit.>. 
 The origin of the ring 
 is not clear but could be related to a past accretion event. 
 Due to the high column-density 
 at which it occurs (∼10^20 cm^-2), 
 it is unlikely that the H i break is related to the breaks driven by ionization. 
 NGC 3621, on the other hand, has an extended H i distribution. 
 The integrated H i intensity map shows a tail or stream of 
 gas extending along the major axis, 
 which could be the result of past interactions or signatures of gas accretion. 
 The break around 10^20 cm^-2 for this galaxy therefore corresponds to the
 transition between the main disk of NGC 3621 and the stream or tail (See also Figure <ref>). 
 The profile of NGC 4826 declines steeply around 5×10^19 cm^-2.
  However, it is unlikely that we are tracing effects due to ionization 
 by extragalactic photons. Kinematical studies 
 by <cit.> and by <cit.> 
 revealed that the inner disk is counter-rotating with respect to the outer disk. 
 As reported by <cit.>, this counter-rotation 
 may indicate that the outer disk of NGC 4826 has been built from the merging of 
 gas rich dwarf galaxies or from the infall of H i gas with a retrograde orbit. 
 Regardless of the exact cause, it is conceivable that we are tracing the 
 relics of a merging event or episodic gas 
 infall rather than the effects of photoionization. In the 
 discussion that follows we do not 
 include the three galaxies discussed in this section.          


 §.§ Comparison with HALOGAS

For two galaxies we have both THINGS and HALOGAS data available. 
This provides us with the opportunity to check for possible systematics. 
In Figure <ref>, we overplot the azimuthally stacked THINGS 
and HALOGAS H i 
column density profiles for NGC 2403 and NGC 3198. 
In both cases, we find excellent agreement in the measurements 
(where they overlap). Also, neither galaxies column density 
show any indication for a break in H i, but feature rather smooth profiles.


 §.§ Identification of any possible breaks and comparison with models

To better identify possible breaks in the azimuthally stacked profiles, we measure 
the change in slope of the profiles by calculating their first logarithmic derivatives. 
We also use the first logarithmic derivatives to gauge whether any break seen 
in the data resembles that of the photoionization models of <cit.> 
and <cit.>. The derivative will be zero for a constant profile, 
negative if it is decreasing and positive if the 
profile is increasing. 
As illustrated in Figure <ref>, 
the change in slope is close 
to zero for many galaxies and does not resemble that of the photoionization models. 
Note that the models 
by <cit.> and by <cit.> only show 
N_H  I(R) for 
R ≥ 20 kpc. Purely for illustrative purposes, we thus extrapolated 
the models to R = 0 kpc. 
The clear difference between the observed and model profiles shows that 
with our present data, we do not find convincing evidence that would suggest the 
presence of a break caused by photoionization 
in our column density profiles.        

To quantify the amount of change in slope, 
we fit the profiles using Sérsic profiles <cit.> 
with different slope indexes, n. 
We here use the formalism introduced by <cit.> who, for n < 1, 
defines the Sérsic profile in terms of a break radius, R_b, instead of the 
half-light radius that is mostly used in the literature. 
Thus, following <cit.>, we define the Sérsic profile as: 
    I(R) = I(0) exp{-(1-n)(RR_b)^1/n},
 
where I(0) is the central intensity. We vary the slope 
indexes, n, but fix I(0) and R_b and 
derive the first logarithmic derivatives of the resulting profiles. 
As seen in Figure <ref>, a Sersic index of 0.46 roughly corresponds to a 
derivative value of -0.2 dex cm^-2 kpc^-1. 
We create a second set of synthetic profiles consisting of a 
straight line without a break and straight lines that break at a fixed radius, R_b, 
with varying slopes. 
We illustrate those results in Figure <ref>. 
Here an order of magnitude decrease within 3 kpc 
(corresponding to f=10 in the Figure) corresponds to a 
derivative value change of -0.3 dex cm^-2 kpc^-1.
Note that the values of the derivatives of the real 
data are typically ≳ -0.2 dex cm^-2 kpc^-1. 
The photoionization 
models by <cit.> and <cit.>, on the other hand, have derivative 
values ≲ -0.3 dex cm^-2 kpc^-1. 
We conclude that, down to our sensitivity limit, 
a radial column density break is not a common feature of disk galaxies. 
More high sensitivity H i observations are needed to investigate 
the H i distribution at even lower column density values.  




 §.§ NGC 3198: Comparison with <cit.>

NGC 3198 is a galaxy for which we do not find a global 
break in neither the THINGS nor the HALOGAS data. However, <cit.> reported a 
break toward the approaching 
side of this galaxy. She showed that the 
column density of NGC 3198 decreased by an order of 
magnitude within a synthesized beam of 2.7 kpc. 
When we redo the analysis for the approaching 
side of NGC 3198, we did not find a noticeable break as reported by 
<cit.>. This prompted us to redo the analysis 
for the approaching and receding sides of all galaxies. 
To allow a 
better comparison of our data with that of <cit.>, 
we divide the radial profiles of the approaching and receding sides of 
the galaxies in 2.7 kpc bins, thus using the same definition for the slope as 
<cit.>. We then calculate the slopes of 
the outer N_H  I(R)≤ 10^20 cm^-2 as 
N_H  I(R)/N_H  I(R+2.7) (see Figure <ref>). 
From this analysis, we do not find evidence that an H i break associated with ionization is a common feature in galaxies. 











                    


§ SUMMARY AND CONCLUSION
    
We studied the radial H i column density profiles of 17 nearby galaxies. 
We used data from the THINGS and, for two galaxies, the HALOGAS survey.
To push the sensitivity limit of our data to deeper levels and trace the H i 
distribution to previously unexplored radii, we stacked individual profiles azimuthally 
using a stacking technique. 
Model velocity fields derived from a simple (flat) 
extrapolation of the rotation curves were used 
to stack individual profiles at large radii. In general the radial column density profiles 
are flat in the inner disk and then smoothly decline down 
to the sensitivity limit of the data. 
Thus, our present 
data do not confirm the prediction of photoionization models that an H i 
break at a column density of ∼ 5× 10^19 cm^-2 is a 
common feature of H i disks. The absence of an H i 
break may indicate that ionization by extragalactic photons is not the limiting 
factor of the extent of the H i disk. Thus, the outskirts of the H i 
disk may instead correspond to the transition to a low column density gas accreted from 
the cosmic web at later evolutionary stages of disk formation. 
However, in order to map the morphology and kinematics of 
the low column density outskirts, 
we require the next generation of radio telescopes <cit.>. 
Until then we can only start to trace the average properties at the outermost radii, 
as done in this paper.
The authors thank the anonymous referee for excellent comments that helped improve the
presentation of this paper. 
R.I. acknowledges funding from the Alexander von Humboldt foundation 
through the Georg Forster Research Fellowship 
programme. 
EB acknowledges support from the UK Science and 
Technology Facilities Council [grant number ST/M001008/1].  
We thank the HALOGAS team for providing us with the NGC 3198 and NGC 2403 data.  


§ AZIMUTHALLY STACKED SPECTRA

In this Section, we show the model rotation curves used to derive the radial surface density profiles of our sample galaxies
(see Figure <ref>). 
We also show examples of the azimuthally averaged stacked spectra at different radius (see Section <ref> for details). 















§ RADIAL COLUMN DENSITY PROFILES AS A FUNCTION OF AZIMUTHAL POSITIONS

Here we present the radial profiles derived inside the annular sectors described in Section <ref>. 
In summary, the galaxies are divided into concentric rings with widths close to the synthesized beam of the observations. 
All rings are then split into 30 degree wide sectors, resulting in annular sectors as shown in Figure <ref> 
(only the outer-most annular sectors are shown). The stacking of individual profiles at a given radius are 
done inside the annular sectors to obtain the radial column density profiles shown in 
Figure <ref>.











[Bergeron & Gunn(1977)]bergerongunn77 Bergeron, J., & Gunn, J. E. 1977, , 217, 892
[Bosma et al.(1977)]bosmaetal77 Bosma, A., van der Hulst, J. M., 
& Sullivan, W. T., III 1977, , 57, 373
[Bosma(1981)]bosma81 Bosma, A. 1981, , 86, 1825 
[Braun et al.(1994)]braunetal94 Braun, R., 
Walterbos, R. A. M., Kennicutt, R. C., Jr., & Tacconi, L. J. 1994, , 420, 558 
[Briggs(1995)]brigs Briggs, D. S. 1995, High Fidelity Deconvolution of Moderately Resolved Sources, PhD Thesis, New Mexico Institute of Mining and Technology
[Broeils & van Woerden(1994)]broeilsetal94 Broeils, A. H., 
& van Woerden, H. 1994, , 107 
[Carignan et al.(2013)]carignanetal13 Carignan, C., Frank, B. S., 
Hess, K. M., et al. 2013, , 146, 48 
[Corbelli et al.(1989)]corbellietal89 Corbelli, E., 
Schneider, S. E., & Salpeter, E. E. 1989, , 97, 390 
[Corbelli & Salpeter(1993)]corbellisalpeter93 Corbelli, 
E., & Salpeter, E. E. 1993, , 419, 104 
[de Blok et al.(2008)]debloketal08 de Blok, W. J. G., Walter, F., 
Brinks, E., et al. 2008, , 136, 2648
[de Blok(2011)]meerblok de Blok, W. J. G. 2011, Tracing the Ancestry of Galaxies, 277, 96 
[Dove & Shull(1994)]doveshull94 Dove, J. B., & Shull, J. M. 1994, , 423, 196 
[Fraternali et al.(2002)]fraternalietal02 Fraternali, F., van Moorsel, G., 
Sancisi, R., et al. 2002, , 123, 3124
[Fumagalli et al.(2017)]fumagallietal2017 Fumagalli, M., 
Haardt, F., Theuns, T., et al. 2017, , 467, 4802 
[Heald et al.(2011)]healdetal11 Heald, G., Józsa, G., 
Serra, P., et al. 2011, , 526, 118
[Heald et al.(2016)]healdetal16 Heald, G., de Blok, W. J. G., 
Lucero, D., et al. 2016, , 462, 1238 
[Ianjamasimanana et al.(2012)]ianjamasimananaetal12 Ianjamasimanana, R., 
de Blok, W. J. G., Walter, F., & Heald, G. H. 2012, , 144, 96
[Ianjamasimanana et al.(2017)]Ianjamasimananaetal17 Ianjamasimanana, R., de Blok, W. J. G., & Heald, G. H. 2017, , 153, 213
[Irwin(1995)]irwin95 Irwin, J. 1995, , 107, 715 
[Irwin et al.(2009)]Irwinetal09 Irwin, J. A., Hoffman, G. L., Spekkens, K., et al. 2009, , 692, 1447 
[Kam et al.(2017)]kametal17 Kam, S. Z., Carignan, C., Chemin, L., et al. 2017, , 154, 41
[Maloney(1993)]maloney93 Maloney, P. 1993, , 414, 41
[Oh et al.(2008)]ohetal08 Oh, S.-H., de Blok, W. J. G., Walter, F., 
Brinks, E., & Kennicutt, R. C., Jr. 2008, , 136, 2761-2781 
[Oosterloo et al.(2007)]oosterlooetal07 Oosterloo, 
T., Fraternali, F., & Sancisi, R. 2007, , 134, 1019
[Portas(2010)]Portas Portas, A. M. P., From giants to dwarfs : 
probing the edges of galaxies,  2010, Ph.D. Thesis
[Rix et al.(1995)]rixetal95 Rix, H.-W. R., Kennicutt, R. C., Jr., 
Braun, R., & Walterbos, R. A. M. 1995, , 438, 155
[Sancisi(1983)]sancisi83 Sancisi, R. 1983, Internal Kinematics 
and Dynamics of Galaxies, 100, 55
[Sersic(1968)]sersic Sersic, J. L. 1968, Cordoba, 
Argentina: Observatorio Astronomico, 1968, 
[Sunyaev(1969)]sunyaev69 Sunyaev, R. A. 1969, , 3, 33
[van Gorkom(1993)]vangorkom93 van Gorkom, J. H. 1993, 
In The Environment and Evolution of Galaxies, ed. J. M. 
Shull and H. A. Thronson, Jr. (Dordrecht, Kluwer), p. 345
[Walter et al.(2008)]walteretal08 Walter, F., Brinks, E., 
de Blok, W. J. G., et al. 2008, , 136, 2563
      
