HI contribution of XQ-100 subDLAs]Sub-damped Lyman $\alpha$ systems in the XQ-100 survey I – Identification and contribution to the cosmological HI budget

Berg et al.]
Trystyn A. M. Berg$^{1,2}$,  Sara L. Ellison$^2$, Rubén Sánchez-Ramírez$^{3,4}$, Sebastián López$^{5}$, Valentina D'Odorico$^{6,7}$, George D. Becker$^{8}$, Lise Christensen$^{9}$, Guido Cupani$^{6}$, Kelly D. Denney$^{10}$, Gábor Worseck$^{11,12}$

$^1$European Southern Observatory, Alonso de Cordova 3107, Casilla 19001, Santiago, Chile.

$^2$Department of Physics and Astronomy, University of Victoria, Victoria, British Columbia, V8P 1A1, Canada.

$^{3}$INAF, Istituto di Astrofisica e Planetologia Spaziali, Via Fosso del Cavaliere 100, I-00133 Roma, Italy.

$^{4}$Instituto de Astrofísica de Andalucía (IAA-CSIC), Glorieta de la Astronomía s/n, E-18008, Granada, Spain.

$^{5}$Departamento de Astronomía, Universidad de Chile, Casilla 36-D, Santiago, Chile.

$^{6}$INAF-Osservatorio Astronomico di Trieste, Via Tiepolo 11, I-34143 Trieste, Italy.

$^{7}$Scuola Normale Superiore Piazza dei Cavalieri, 7 I-56126 Pisa, Italy.

$^{8}$Department of Physics & Astronomy, University of California, Riverside, California, 92521, USA.

$^{9}$DARK, Niels Bohr Institute, University of Copenhagen, Lyngbyvej 2, 2100 Copenhagen, Denmark

$^{10}$Department of Astronomy, The Ohio State University, 140 West 18th Avenue, Columbus, OH 43210, USA.

$^{11}$Max-Planck-Institut für Astronomie, Königstuhl 17, D-69117 Heidelberg, Germany.

$^{12}$Institut für Physik und Astronomie, Universität Potsdam, Karl-Liebknecht-Str. 24/25, D-14476 Potsdam, Germany.

Sub-damped Lyman $\alpha$ systems (subDLAs;  column densities of $19.0\leq$ logN() $<20.3$) are rarely included in the cosmic  census performed at redshifts  $\gtrsim1.5$, yet are expected to contribute significantly to the overall  mass budget of the Universe. In this paper, we present a blindly selected sample of 155 subDLAs  found along 100 quasar sightlines (with a redshift path length $\Delta X=475$) in the XQ-100 legacy survey to investigate the contribution of subDLAs to the  mass density of the Universe. The impact of X-Shooter's spectral resolution on Ly$\alpha$ absorber identification is evaluated, and found to be sufficient for reliably finding absorbers down to a column density of logN()$\geq18.9$. We compared the implications of searching for subDLAs solely using  absorption versus the use of metal lines to confirm the identification, and found that metal-selection techniques would have missed 75 subDLAs. Using a bootstrap-Monte Carlo simulation, we computed the column density distribution function () and the cosmological  mass density () of subDLAs and compared with our previous work based on the XQ-100 damped Lyman $\alpha$ systems. We do not find any significant redshift evolution in  or  for subDLAs. However, subDLAs contribute 10–20 per cent of the total  measured at redshifts $2<z<5$, and thus have a small but significant contribution to the  budget of the Universe.

galaxies: high redshift – galaxies: ISM – quasars: absorption lines

§ INTRODUCTION

Studying the evolution of neutral gas reservoirs over cosmic time provides key constraints on aspects of galaxy evolution. The neutral atomic gas, primarily traced by , is thought to be an indicator of the gas reservoirs that eventually form stars in galaxies [Bird et al, 2015, Somerville & Davé, 2015]. Furthermore, the rate at which the  reservoirs change places direct constraints on the processes that can nourish or prevent future star formation [Davé et al, 2013, Lilly et al, 2013].

The comoving  mass density (; measured with respect to the critical mass density of the Universe) has typically been used to quantify the amount of neutral gas in galaxies. At $z\approx 0$,  is quantified from 21-cm emission density maps [Zwaan et al, 2005, Martin et al, 2010, Braun, 2012]. However at redshifts $z \gtrsim 0.2$ where current radio facilities cannot easily detect the faint 21-cm emission, quasars (QSO) are used to probe sightlines through intervening  clouds in the form of Ly$\alpha$ absorption <cit.>.

Of the many classes of Ly$\alpha$ absorption line systems seen towards QSOs, damped Lyman $\alpha$ systems (DLAs), defined as having  column densities of logN()$\geq 20.3$ [Wolfe et al, 1986, Wolfe et al, 2005], are typically used to measure . Although few in number compared to lower  column density counterparts, DLAs dominate the  column density distribution from $\sim5$ to the present epoch <cit.>. However, subDLAs ($19.0\leq$logN()$<20.3$) are thought to host a modest contribution of $\approx10$–20 per cent to the  budget of the Universe [O'Meara et al, 2007, Zafar et al, 2013]. The difficulty in using subDLAs to probe the neutral gas reservoirs of the Universe is that these systems are often not fully self-shielded from the cosmic UV background, and thus do not completely trace neutral gas reservoirs <cit.>. Nevertheless, there has been a substantial body of work calculating the contribution of  from subDLAs to the column density distribution function () and  [Péroux et al, 2003, Péroux et al, 2003, Péroux et al, 2005, O'Meara et al, 2007, Guimarães et al, 2009, Zafar et al, 2013], as well as understanding whether subDLAs probe different types of gaseous systems than DLAs [Péroux et al, 2003, Kulkarni et al, 2007, Dessauges-Zavadsky et al, 2009, Quiret et al, 2016].

It is expected that the observed  can place strong constraints on feedback and ionization prescriptions in cosmological simulations of neutral gas reservoirs in the Universe <cit.>. Both observational [Péroux et al, 2003, Péroux et al, 2003, O'Meara et al, 2007, Guimarães et al, 2009, Zafar et al, 2013] and simulation-based studies [Rahmati et al, 2013, Rahmati et al, 2015, Villaescusa-Navarro et al, 2018] of the nature of  for subDLA column densities have suggested that subDLAs are more common at redshifts $>3.5$ and contribute more than $\approx10 - 20$% of   traced solely by DLAs [Péroux et al, 2005, Guimarães et al, 2009, Noterdaeme et al, 2009, Zafar et al, 2013]. To date, samples of subDLAs have typically been selected using the presence of metal lines rather than using solely  absorption as is done with DLAs. Such a selection effect can potentially lead to sample incompleteness from missing metal-poor systems, and introduce biases in the comparison with DLAs. Additionally, most subDLA samples are too modest in size to detect any redshift evolution  at $>2\sigma$ confidence  <cit.>.

This paper aims to investigate the contribution of subDLAs to the  budget at high redshifts ($2.5\leq$  $\leq4.5$) using a completely blind sample of subDLAs identified in the XQ-100 survey [López et al, 2016] solely based on their Lyman series absorption. We searched for all  absorbers in the XQ-100 spectra down to a column density threshold of logN()$\geq18.8$ (i.e. the subDLA threshold minus the typical logN() error of $\pm0.2$ dex) in order to evaluate the completeness of our absorber identification and properly account for miscalssification of absorbers as subDLAs in our computation of  and . However, the standard subDLA  column density limits are still adopted in all our computations. The reader is referred to [Sánchez-Ramírez et al, 2016] for the HI statistics of the XQ-100 DLA sample. Through-out this paper we assume a flat $\Lambda$CDM Universe with ${\rm H_{0}=70.0~km~s^{-1}~Mpc^{-1}}$ and ${\rm \Omega_{M,0}=0.3}$.

§ DATA

§.§ XQ-100 Survey design and data reduction

The XQ-100 Legacy Survey (PI S. López) observed 100 QSO sightlines at redshift $\sim$ 3.5–4.5 with the X-Shooter Spectrograph [Vernet et al, 2011] on the Very Large Telescope (VLT) in Chile. The QSOs were purposefully chosen to be blind to any intervening absorption line systems, thus providing a random sample of sub-DLAs and DLAs to study the cosmological implications of these intervening absorbers  [Sánchez-Ramírez et al, 2016, Berg et al, 2016, Berg et al, 2017, Christensen et al, 2017]. For more details about the XQ-100 science cases and survey design, see [López et al, 2016].

For each QSO in the XQ-100 survey, the per-arm exposures were either $\sim0.5$ or $\sim1$ hour in length (depending whether the QSO was classified as `bright' or `faint'; respectively), providing median signal-noise ratios of $\approx30$ pixel$^{-1}$. The spectra for each arm were reduced using an internal IDL package developed by G. Becker and G. Cupani. For more details on the data reduction, see [López et al, 2016].

As the Ly$\alpha$ forest and intervening gas clouds along every sightline partially absorb all QSO photons below the Lyman limit, the QSO's flux is systematically absorbed for all wavelengths below the corresponding observed Lyman Limit at the redshift of the QSO ($\lambda_{LL}=(1+z_{\rm em})\times 911.7633~{\rm \AA{}}$), making both accurate continuum placement and subDLA identification difficult below this wavelength. For this reason, the wavelength regime below the Lyman limit is consequently not used for identifying the subDLAs in this work. The resulting redshift path length searched for subDLAs is $\Delta X=475$.

§.§ Identification of intervening absorbers

Starting with the manually placed QSO continuum fits outlined in [López et al, 2016], we visually identified candidate subDLAs by systematically scrolling through the XQ-100 spectra starting from the Ly$\alpha$ emission of the QSO to $\lambda_{LL}$.  Candidate absorbers are first identified by any saturated absorption broad enough to fit a simulated Ly$\alpha$ Voigt profile of column density logN()$\geq18.8$ <cit.>. Note that our adopted logN()$\geq18.8$ search threshold is below the minimum subDLA threshold (logN() $=19.0$) by the typical  column density error (i.e. $\pm0.2$) in order to enable a complete bootstrap error analysis (Section <ref>). Simultaneously to the Ly$\alpha$ fitting, a simulated higher order Lyman series line profile (usually Ly$\beta$) was checked at the corresponding wavelength. For each of these saturated regions identified, one or more Voigt profiles were fitted by adjusting simultaneously the original continuum fits, the logN(), and  until the modelled Ly$\alpha$ and Ly$\beta$ Voigt profiles matched the continuum-normalized spectra. Continuum adjustments were common for large column density absorption systems on top or near QSO emission lines. In many of these adjustments, the best-fit column logN() was constrained by simultaneously fitting both the Ly$\alpha$ and Ly$\beta$ profiles. Errors on the logN() column density were estimated as the minimum and maximum possible logN() that would provide a visually reasonable fit to the data. For cases where the continuum was uncertain, the assigned logN() error was increased to incorporate this uncertainty in the continuum placement.

To ensure higher redshift Lyman series lines were not confused with lower redshift Ly$\alpha$ absorption, a full model spectrum was generated to keep track of all fitted absorption systems along each sightline. Upon identifying a new absorber, the associated Voigt profiles for several key Lyman series lines (from Ly$\alpha$ to Ly$\delta$) were added to the full model spectrum. The identification of systems was terminated when the observed-frame wavelength of additional Ly$\alpha$ absorbers was shorter than $\lambda_{LL}$.

Once the absorber identification process was complete, the initial list of absorption line systems was cleaned by adjusting or removing any absorbing systems whose modeled Lyman series lines (Ly$\alpha$–Ly$\delta$, and the highest order Lyman series lines between 911Å and 920Å) did not fit the continuum normalized spectra. The final cleaned sample of identified candidate absorbers consists of 87 $18.8\leq$logN()$<19.0$ Lyman limit systems, 155 sub-DLAs and 34 DLAs. Figure <ref> shows the best-fitting Ly$\alpha$ and Ly$\beta$ profiles  for a representative sample of 14 absorbing systems, while the remaining profiles are displayed in Appendix <ref> (Figure <ref>). Table <ref> contains the full list of 276 candidate absorbers, including their  and logN().

Given the moderate resolution of the XQ-100 spectra, it is sometimes difficult to model regions of extended saturated absorption where the damping wings of the Voigt profile are not visible into individual absorption features. The troughs could be fit reasonably by one or more subDLAs, but also by several lower column density Lyman Limit systems ($17<$ logN() $<19$) that are clustered together.  Unless a higher-order Lyman series line could disentangle this degeneracy (such as the absorbers seen towards J0034+1639 at $\approx4.2$; see Figure <ref>), no candidate absorbers were fitted to the data in these degenerate regions. In many cases, we note that the edges of the Lyman series absorption appear to be better fit by Lyman Limit systems with high Doppler parameters rather than subDLAs. We therefore caution that our identification of subDLAs is potentially incomplete.

[]Each pair of panels displays the fitted Ly$\alpha$(first and third columns) and Ly$\beta$ (second and fourth columns) profiles (dark red lines) and associated errors (dark red shaded region) for 14 of the identified candidate system. The light blue line and shaded region denote the full model spectrum and associated error combined from all absorbers within the entire sightline. The red text denotes the properties of the absorber [QSO sightline name and redshift in the corresponding Ly$\alpha$ panel, and logN() in Ly$\beta$ panel], while the vertical black text labels the redshift and Lyman series lines associated with additional absorption included in the model. The x-axis of each panel is scaled such that integrated region of the Ly series line's Voigt profile displayed encompasses 98 per cent of the absorbed flux. The remaining 262 systems are displayed in Figure <ref>.

absorbers in XQ-100

QSO            logN()    Metal absorption

J0003-2603    4.010    2.78140    19.00$\pm$0.20    N

J0003-2603    4.010    3.05490    20.00$\pm$0.15    Y

J0003-2603    4.010    3.16080    19.10$\pm$0.15    Y

J0003-2603    4.010    3.39020    21.45$\pm$0.05    Y

J0003-2603    4.010    3.43470    18.90$\pm$0.10    N

J0006-6208    4.460    3.16890    18.90$\pm$0.20    N

J0006-6208    4.460    3.18430    18.90$\pm$0.15    N

J0006-6208    4.460    3.20250    20.80$\pm$0.20    Y

J0006-6208    4.460    3.45890    19.30$\pm$0.20    N

J0006-6208    4.460    3.77300    20.90$\pm$0.15    Y

§.§.§ Comparison to [Sánchez-Ramírez
et al, 2016]
In [Sánchez-Ramírez et al, 2016], we identified 41 DLAs in the XQ-100 sample using a similar absorber identification and fitting method. However, the identification of DLAs was completed over a slightly longer redshift path in <cit.> in that it was not truncated at $\lambda_{LL}$ as was done here. As a result, the initial search from [Sánchez-Ramírez et al, 2016] includes eight more DLAs, which have been removed from the sample in this paper. We also identify two XQ-100 DLAs missed in the [Sánchez-Ramírez et al, 2016] catalogue; one system is towards J0214$-$0517 ( $=3.7211$, logN() $=20.6\pm0.2$), and the other is along the J1111$-$0804 sightline ( $=3.60690$, logN() $=20.4\pm0.15$).

In order to facilitate a fair comparison between the  contributions of subDLAs and DLAs in XQ-100, it is important that we use the same continuum fits for both the DLA and subDLA identification; thus we re-measured the N() values for all the previously identified XQ-100 DLAs. All of the [Sánchez-Ramírez et al, 2016] DLAs within the common redshift range of the two studies were recovered, and their column densities are consistent (within the  column density errors) with those that we derived in our subDLA search. The largest discrepancies occurred in regions where the continuum is difficult to constrain, such as systems that are nearby the background QSO's Ly$\alpha$ emission line. For this paper, we use the newly measured logN() values of the XQ-100 DLAs for consistency.

§.§ Completeness analysis

§.§.§ Recovery test

To assess the completeness of the visual identification process and ensure that the modest XQ-100 spectral resolution (FWHM resolution R$\sim$5400–8900) did not affect the determination of column densities, an absorption line recovery test was necessary. To ensure the Ly$\alpha$ forest statistics and X-Shooter instrumental effects are reproduced while maintaining the same continuum fits used in the original identification process, artificial absorption systems were randomly inserted into the XQ-100 spectra.

Prior to insertion, a list of artificial systems was first generated to ensure that an approximately uniform distribution of column densities was used. In bins of column density from logN()$=18.0$ to 20.3[The recovery rate for DLAs in the XQ-100 dataset was previously measured at 100% in [Sánchez-Ramírez] et al, 2016]. Given that we recovered all the previously identified DLAs in an identical recovery test in [Sánchez-Ramírez et al, 2016], we did not repeat the experiment for DLAs. (of bin width 0.1 dex), $15\pm3$ systems were generated with column densities uniformly selected within the bin's limits. Therefore the total number of systems (between 276 and 414 systems) was randomly generated, preventing any subconscious counting during the recovery process. For each of these fake systems, a model Lyman series (down to Ly$\delta$) spectrum was generated, and randomly inserted into one of the XQ-100 QSO spectra at a random redshift. To simulate noise in the model spectrum, the flux of the model spectra was then varied randomly using the XQ-100 error spectrum flux as the standard deviation of a Gaussian distribution. Between 0 and 7 fake absorption systems were assigned to each XQ-100 sightline.

The identification process described in Section <ref> was repeated on the spectra with fake absorbers inserted, starting with the locations of the previously-identified real absorbers flagged in the model spectra. The top panel of Figure <ref> shows the column density distribution of fake systems inserted into the XQ-100 spectra (with observed-frame Ly$\alpha$ absorption at wavelengths larger than $\lambda_{LL}$; black line). The solid red bars indicate that we missed seven, logN()$\leq19.2$ systems in our recovery analysis.  However, the three missed systems with the higher column densities (logN()$\geq18.9$) happened to be randomly inserted on top of already identified, higher column density Ly$\alpha$ absorbers. Since the  contribution from real overlapping systems would have already be accounted in the original fit to the larger  systems, we do not include these three missed systems in the recovery rate derivation. Therefore only $\approx 27$ percent of systems with logN()$=18.8$ were missed in the recovery analysis.

For all the recovered systems, the bottom panel of Figure <ref> shows the difference in logN() between the measured column density of the fake absorber (N$_{\rm Meas}$) and the randomly generated value (N$_{\rm Fake}$). All but ten of the recovered systems ($\approx 3$% of the total number of systems) have consistent column densities with N$_{\rm Fake}$. For these ten systems, the column densities are only inconsistent by  $\sim0.1$ dex, which is statistically consistent with the random flux errors inserted into the fake model spectra.

In summary, the  recovery test suggests that our identification of XQ-100 absorbers is 100% complete for logN()$\geq18.9$, and 73% complete for logN()$=18.8$.

[XQ-100 N() recovery statistics ]Top panel: Distribution of the number of fake systems of a given  column density (logN$_{\rm Fake}$) used for the recovery test (black line) and those systems that were missed in the recovery process (solid red histogram). Note that only systems whose Ly$\alpha$ is beyond $\geq \lambda_{LL}$ are included in this plot. Bottom panel: The difference in measured and fake logN() column densities (i.e. logN$_{\rm Meas}$-logN$_{\rm Fake}$) as a function of logN$_{\rm Fake}$.

§.§.§ Effects from continuum placement

The location of the continuum becomes more uncertain towards bluer wavelengths as there is more absorption from additional Lyman series lines (i.e. Ly$\beta$ and Ly$\gamma$) removing regions of clean continuum. To test whether or not our continuum placements introduce significant errors in our column density determinations, we artificially inserted a fake absorption line in the un-normalized spectrum and fitted the line using our continuum placements (without continuum adjustment). This experiment was done in three regions within the rest-frame of the QSO, near the wavelengths: $\lambda=$ 1120.7, 999.1, and 961.1  Å. These three regions correspond to locations in the QSO spectrum where there is only Ly$\alpha$, Ly$\alpha+$Ly$\beta$, and Ly$\alpha+$Ly$\beta+$Ly$\gamma$ absorption (respectively). In order to remove the effects of individual sightlines on the experiment, five QSO spectra were selected at random for each region. However, we required that each of the selected QSOs had no previously identified Ly$\alpha$ absorption within the region where the line was inserted to avoid any effects from previous continuum placement modifications during the Ly$\alpha$ fit. The column density of the fake absorption line started at logN() $=18.8$, and was increased in steps of 0.2 dex up to a column density of logN() $=20.0$ (inclusive) to test the effects as a function of column density. For each fake absorption line, the wavelength was held fixed. For all but one QSO, we were able to recover the column density to within 0.2 dex (the typical quoted error on the  columns). These errors were not systematically offset in one particular direction (especially in the Ly$\alpha+$Ly$\beta$ and Ly$\alpha+$Ly$\beta+$Ly$\gamma$ region), suggesting there is no significant effect from the increased absorption (and thus more uncertain continuum) on our column density determinations. We note that for the one QSO sightline which provided inconsistent fits, the identification of the line was difficult due to significant blending with Ly$\alpha$ forest for the inserted Ly$\alpha$ profiles with column densities logN() $=19.0$ and 19.2, and the offsets were measured to be $+0.3$ dex each (i.e. the column densities were over-estimated), and is likely not strongly influenced from the continuum placement.

§.§.§ N() comparison to higher resolution data

Although the resolution of the XQ-100 data (FWHM resolution R$\sim$5400–8900) is sufficient to resolve the absorbing systems, it may be inadequate to resolve close blends of multiple lines. In order to test the absorber recovery and fit accuracy, we compared the fits derived from our X-Shooter data with spectra obtained with high resolution (R$\approx 80 000$) observations with the UVES instrument.  Using all available archival UVES spectra[Data accessed from <http://archive.eso.org/wdb/wdb/eso/uves/form>.] of the XQ-100 sightlines (19 sightlines total; with signal-noise of $\approx5-15$ pixel$^{-1}$), we repeated the Ly$\alpha$ profile fitting at the already identified redshifts of 62  absorbers with UVES coverage.  Comparing the obtained logN() from both the UVES and XQ-100 fits, all column densities are consistent within the errors (at most 0.2 dex difference), with two exceptions. For these particular exceptions, the XQ-100 absorbers (J0134$+$0400, =3.999; J0137$-$4224, =3.665) each appear to be better fit by two, smaller component systems using the UVES data. However, the total  column measured contained within the UVES absorbers remains consistent with the single  profile fit measured in the XQ-100 data.

Despite the overall consistency between the UVES and XQ-100 data, we note a systematic difference in logN() exists between the two datasets. Figure <ref> shows the mean N() offset between the UVES and XQ-100 ($\Delta$logN(); measured relative to the XQ-100 value) as a function of logN(), and indicates that the UVES data tend to provide lower logN() fits than the XQ-100 data for absorbers with logN()$<19.8$ by up to 0.1 dex on average. All the offsets are tabulated in Table <ref>. Table <ref> also includes the fraction of system which have a non-zero  $\Delta$logN() (n$_{\rm off}$) relative to the total number of absorbers (n$_{\rm tot}$) within each column density bin. As both the VIS (R$\approx$8900) and UVB (R$\approx$5400) arms of X-Shooter have different spectral resolutions, we checked for a dependence of the measured $\Delta$logN() for absorbers identified in the respective arms. Absorbers detected in the VIS arm have a slightly smaller $\Delta$logN() than those found in the lower resolution UVB arm by at most $\approx0.05$ dex, so the effects of the differing arm resolutions are minimal.

[]The red points show the mean offset of logN() measured between UVES and XQ-100 ($\Delta$logN(); positive values imply a higher UVES-derived column density) in bins of column density. Error bars show the standard deviation in each bin, while the number in parenthesis above each data point indicate the number of absorbers within the bin. Offsets measured for individual absorbers are displayed as the grey circles, where a darker grey represent a higher density of systems with the same logN() and $\Delta$logN().

Mean logN() offsets between UVES and XQ-100 data

XQ-100 logN() range     $\Delta$logN()     n$_{\rm off}$/n$_{\rm tot}$

$18.8\leq$logN()$<19.0$     $-0.08$     16/23

$19.0\leq$logN()$<19.2$     $-0.10$     8/13

$19.2\leq$logN()$<19.4$     $-0.06$     4/6

$19.4\leq$logN()$<19.6$          0/0

$19.6\leq$logN()$<19.8$     $-0.03$     1/3

$19.8\leq$logN()$<20.3$     $0.0$     0/4

logN()$\geq20.3$     $0.0$     0/5

§.§.§ False positives and metal lines

In addition to ensuring that all the subDLAs and DLAs were identified in the XQ-100 dataset, it is possible that two closely-separated Lyman limit systems or other Lyman series absorption can mimic a low column density subDLA, leading to falsely identified subDLA systems.  For low-to-modest resolution spectra or low redshift absorption systems, the Lyman series lines are difficult to detect within the blanket Lyman limit absorption, or separate from the Ly$\beta$ forest. One method to determine the false positive rate is to use higher resolution data where it is easier to resolve distinct absorption lines. In Section <ref>, we found only two cases of the 63 systems where two Ly$\alpha$ lines were a better fit to the higher resolution UVES data, suggesting these false positives are not a common occurrence in the XQ-100 data. However this is not a complete assessment of the false positive rate as it is still possible that lines can be blended at even closer separations.  We caution that our comparison with the UVES data does not uniformly probe the range in wavelength (with respect to the rest frame of the quasar) to properly assess the frequency of mis-identification of systems due to blending.

[N() vs.  for all XQ-100 absorber]The  column density (logN()) of all candidate XQ-100 absorbers (light grey diamonds) as a function of their redshift (). The red circles indicate the XQ-100  metal-selected absorbers. The total number of absorbers (N) for each sample are given in the top left corner of the main panel. The distributions of  and logN() are given in the top and right-most panels, using the same colour scheme for the two subsamples. For reference, the dotted and dashed horizontal lines denote the column density thresholds for subDLAs and DLAs (respectively).

Another common method to remove these `false positive' absorbers from a sample is to check for the presence of higher order Lyman series lines or metal lines associated with the absorber, particularly using high resolution data [Péroux et al, 2003, O'Meara et al, 2007, Zafar et al, 2013]. The typical metal lines used are: Cii, Civ, Siii, Siiv, Mgii, or Feii. However, using metal lines for confirmation of absorbers may either: (i) artificially remove systems with weak or no metal lines that are not easily identified (i.e. metal-poor or systems with large dust or ionization corrections to neutral species), or (ii) select blends of lower N() absorbers which still have strong metal lines. Either of these two cases can introduce different biases into the sample.

From the full sample of candidate absorption systems with logN()$\geq18.8$, we identified a subsample of metal-selected absorbers, where at least one of the following ionic species is detected: Civ ($\lambda\lambda$ 1548 Å & 1550 Å), Siiv ($\lambda\lambda$ 1393 Å & 1402 Å), Mgii ($\lambda\lambda$ 2796 Å & 2803 Å), or Feii ($\lambda$ 2344 Å, 2382 Å, & 2600 Å). The systems with observed metal lines are identified by a `Y' in the `Metal absorption detected?' column of Table <ref>. Figure <ref> shows the distribution of logN() and  for all the absorbers with logN()$\geq18.8$ (grey diamonds) and metal-selected absorbers (red circles). We point out that all 34 DLAs have metal lines associated.  It is clear that the metal selection process preferentially removes low redshift, low -column density (logN()$\leq19.4$) candidate absorbers from the full list of absorbers; these removed systems tend to be found where higher-order Lyman series blending is more likely to produce false positive Ly$\alpha$ absorption.

To understand what is driving the large discrepancy between the number of absorbers between metal-selected sample and all identified absorbers at low  and , Figure <ref> shows the cumulative number of subDLAs as a function of redshift for all the absorbers (solid black lines) and the metal-selected sample (dashed black lines). The distributions are split up into three panels for three different column density ranges. For comparison, the shaded regions represent the number of absorbers whose commonly-detected metal diagnostic lines (i.e. Civ 1548Å and Siiv 1393 Å; light blue and magenta, respectively) are outside of the Ly$\alpha$ forest. Although only Civ and Siiv doublets have been displayed, other commonly-observed low-ionization species (e.g. Siii 1260Å, Cii 1334Å) are bluer than the Siiv doublet, or are often blended with telluric lines at these redshifts (e.g. Feii 2586Å, and the Mgii doublet at 2796Å and 2803Å); and thus these lines are used less frequently for metal line identification in this work. Eleven of the absorbers with confirmed metal lines ($19.3\leq$ logN() $\leq20.2$) have been confirmed using only low ionization metal lines; eight of these systems have blended Civ or Siiv, and one subDLA (logN() $=19.8$) has no associated Civ and Siiv absorption.

The dark red shaded regions of Figure <ref> show the number of absorbers whose Ly$\beta$ is redward of $\lambda_{\rm LL}$ to summarize the frequency in which can confirm the presence of  absorption through Ly$\beta$. The redshift at which the metal-selected dashed line is greater than the red shaded region ( $\lesssim3.5$) denotes the highest redshift where it is no longer possible to remove false-positives in a metal-blind fashion using the Lyman series. This method for removing false-positives appears to be successful for the XQ-100 subDLAs, as the solid line for all logN() $\geq 19$ absorbers are very similar above  $\gtrsim3.5$. However, the large red-shaded region between the solid and dashed lines in the leftmost panels suggests that using Ly$\beta$ absorption to remove false-positives becomes challenging for systems with low  column densities (logN() $\lesssim19.4$) up to  $\approx4$ in the XQ-100 data. Upon investigation, this difficulty in removing false-positive using Ly$\beta$ is mostly due to blending of the weaker Ly$\beta$ line of these absorbers with the Ly$\alpha$ forest. As a result it is more difficult to confirm the observed absorption is primarily from Ly$\beta$ for logN() $<19$ XQ-100 absorbers at all redshifts.

Of the 86 systems without detected metal lines, 77 of these systems have the corresponding Civ and Siiv  absorption in the Ly$\alpha$ forest; this is represented by the unshaded regions below the solid black line.  Thus for nearly 40 per cent of all identified subDLAs, we cannot determine if the commonly detected metal lines are present and estimate the false positive rate. The subDLAs which populate the unshaded region are typically found at  $\lesssim3$ as a result of the XQ-100 QSO redshifts. We point out that this percentage is much less for subDLAs with logN()$\geq19.8$ as neutral metal species typically outside the forest and not in telluric lines (e.g. Feii 1608Å and Alii 1670Å) are more often detected due to larger column densities of gas.

In Figure <ref>, signatures of false positives would manifest as shaded regions above the black dashed lines (i.e. no diagnostic absorption lines are detected). We caution that this shaded regime can also be driven by many different factors other than false-positives such as regions where the line is blended (particularly for Ly$\beta$; red shaded region), systems which are metal-poor (blue and magenta regions only), or low signal-noise data. For the 14 absorbers which have undetected Civ or Siiv in uncontaminated regions of the spectrum, we measured detection thresholds of logN(Civ)$\leq12.97$ and logN(Siiv) $\leq12.76$. For Civ (which should be easily identified in absorbers within the shaded regime), the column density detection threshold is higher than the expected column densities observed in subDLAs from the literature <cit.>, and thus the XQ-100 sample can be limited by insufficient signal-noise to detect Civ for absorbers with logN() $<19.4$ at lower redshift ($z<3$). However, we note that generally higher metallicities are observed in low redshift subDLAs compared to DLAs [Pettini et al, 1999, Rafelski et al, 2014, Quiret et al, 2016].

To better understand whether the XQ-100 spectra have the capability of detecting low N(), metal-poor absorbers, we computed the signal-noise ratio ($SNR$) required to detect a neutral metal line absorption line (with rest-frame wavelength $\lambda$ and oscillator strength $f$), at 3$\sigma$ significance, for an absorber with an observed metallicity [M/H] (i.e. no dust or ionization corrections) using

\begin{equation}
\label{eq:snr}
SNR = \frac{3 m_{e} c ~ FWHM}{\pi e_{c}^{2} f \lambda^{2}}\left(N(\HI)10^{[M/H]+log\left(\frac{n_{\rm M_{\odot}}}{n_{\rm H_{\odot}}}\right)}\right)^{-1};
\end{equation}

where $m_{e}$ and $e_{c}$ are the electron mass and charge, and $\frac{n_{\rm M_{\odot}}}{n_{\rm H_{\odot}}}$ is the observed solar abundance of the species from [Asplund et al, 2009]. Equation <ref> assumes the metal line is detected at the corresponding full width half-maximum ($FWHM$; measured in the observed frame) of the instrumental resolution of the NIR arm ($R\approx 5600$).

Figure <ref> shows the range of signal-noise required to detect the strong neutral metal lines Feii $\lambda$ 2382 Å (left panel) and Mgii $\lambda$ 2796 Å (right panel), as a function of logN() and observed metallicity ([M/H]; coloured lines). The typical range of $SNR$ in the XQ-100 spectra is shown by the grey shaded regions. Figure <ref> demonstrates that subDLAs should have detectable Feii 2382 Å in our X-Shooter data for [Fe/H] $>-2.0$. Few subDLAs have previously been detected with an observed [M/H]$<-2$ [Quiret et al, 2016, Fumagalli et al, 2016]. Therefore, the ability of our X-Shooter data to detect metal lines for typical subDLA metallicities suggests that up to $\approx46$ per cent (128 systems which do not have confirmed metal lines) of the candidate subDLAs might be lower column density or blended systems. However the lack of low-metallicity subDLAs in the literature could be a result of metal-line selection, particularly as DLAs can probe down to [M/H]$\approx -3.0$ [Cooke et al, 2015, Quiret et al, 2016]. A detailed analysis of the metal abundances of the XQ-100 subDLAs will be presented in a future paper (Berg et al., in prep).

In summary, the region between the solid (all identified absorbers with logN()$\geq18.8$) and dashed lines (metal-selected absorbers) in Figure <ref> represents a combination of the false positives in our full sample of absorbers, metal-poor absorbers ([M/H]$\lesssim-2$), and absorbers whose metal lines lie in regions suffering from blending or insufficient signal-noise to classify into the metal-selected sample. As these two curves encompass the two extremes in the number of subDLAs identified, we use all the identified absorbers and metal-selected sample as bounds to the true number of subDLAs within the XQ-100 sample for our analysis. It is worth emphasizing that the discrepancy in the number of absorbers identified that do not have confirmed metal lines is strongest at low column densities (logN() $\leq19.4$) and for systems at low redshifts where metal lines are shifted into the Ly$\alpha$ forest ( $\lesssim3.2$). Although the false-positive rate may be significant at low column densities ($19.0\leq$ logN() $\leq19.4$), the difference in their contribution to the total  mass budget of subDLAs and DLAs is only a couple of per cent (see Section <ref>).

The cumulative distribution of XQ-100 absorbers (with $18.8\leq$ logN() $<20.3$) above a given redshift $z$ ($n_{\rm abs}(z_{\rm abs}\geq z)$). The four panels split the absorbers into different logN() bins ([18.8,19.0), [19.0,19.4), [19.4,19.8), [19.8,20.3); from left to right). The solid and dashed black line denote the curves using all the identified absorbers and metal-selected systems (respectively). The shaded regions show the distribution of absorbers whose key diagnostic lines (Ly$\beta$, Civ, and Siiv; dark red, light blue, and magenta, respectively) are in uncontaminated regions of the spectrum (i.e. Ly$\beta$ redward of the QSO's Lyman Limit; metal lines outside of the Ly$\alpha$ forest). Unshaded regions below the solid line denote the number of absorbers that cannot be confirmed as metal-selected using either Civ or Siiv, and contains $\approx40$ per cent of all the absorbers identified (mostly below  $\lesssim3.2$ and logN() $\lesssim19.4$). The shaded region above the dashed line estimates the inability to detect the given metal line (due to low signal-noise, blending, or being metal-poor) combined with false-positives. This regions is significantly populated by logN()$<19.4$ absorbers.

The predicted spectral $SNR$ required to observe the metal lines Feii $\lambda$ 2382 Å (left panel) and Mgii $\lambda$ 2796 Å (right panel) as a function of absorber logN().  A range in metallicities is represented by the different coloured lines. The typical range in signal-noise of the XQ-100 spectra is denoted by the light grey region.

§ RESULTS

In [Sánchez-Ramírez et al, 2016], we computed both  and  for the XQ-100 DLA sample and analyzed the evolution of . In this section, we will present the added information obtained from calculating  and  from the XQ-100 subDLA sample. For consistency with [Sánchez-Ramírez et al, 2016], we restricted the XQ-100 sample to exclude systems with an absorption redshift within 5000  of the QSO emission redshift. These proximate systems likely trace a different population of systems compared to their intervening counterparts, [Ellison et al, 2002, Ellison et al, 2010, Berg et al, 2016, Perrotta et al, 2016] and are typically ignored when computing .  To properly measure  for DLAs, a sample size much larger than the 33 DLAs identified in the XQ-100 sample must be used in order to statistically probe the high column density end of  (logN()$>21.5$). We thus combined the XQ-100 DLA sample with literature sample compiled in [Sánchez-Ramírez et al, 2016], which consists of 1019 DLAs measured in various literature surveys across a redshift range of $1.6<$$<5.1$.  It is important to note that this supplementary DLAs sample does not include any subDLAs, and probes a very different redshift path than the XQ-100 subDLAs.  For this reason, we divide the entire XQ-100 absorber sample into DLAs and subDLAs in the rest of this analysis.

§.§ Bootstrap Monte Carlo re-sampling

To estimate the errors in deriving  and  from the XQ-100 dataset, the bootstrap Monte Carlo (further abbreviated B-MC) resampling technique presented in [Sánchez-Ramírez et al, 2016] was used. This technique provides a robust method of combining the random errors quantified for measurements (such as logN()) with the sampling errors inherent to a relatively small sample size. In essence, the B-MC method takes the original absorber N() distribution and duplicates the sample a large number of times to create N$_{\rm samp}$ simulated mock surveys. For each mock survey, the logN() of each absorbers is randomly varied within its error-bar (assuming a Gaussian error distribution) to create a `new' set of logN() measurements. From this new set of measurements, absorbers are randomly drawn (with replacements) until a new mock survey is generated with the same number of absorbers. As the drawing of values is done with replacement, this mimics the effect of sampling errors by randomly excluding or duplicating absorbers in the sample. The choice of N$_{\rm samp}$ is arbitrary, but needs to be very large ($\gtrsim 10000$) in order to provide a robust statistical sample.

We also experimented with statistically incorporating the effects of using the modest resolution of X-Shooter to measure  column densities over high resolution datasets (e.g. UVES; see Section <ref>). Using the frequency and average measured offset seen between the X-Shooter and UVES column densities (as presented in Table <ref>), we varied a random fraction of the column densities in the B-MC simulation by the average measured N() offset. After implementing these offsets, we did not detect a significant difference in the results of this paper. In addition to not having a sufficient sample size of absorbers with UVES coverage to provide robust statistics for different redshifts and column densities (see Section <ref>), we did not implement this resolution experiment into the final B-MC simulation.

As the XQ-100 absorber sample is split into two subsamples based on column density (subDLAs and DLAs), the B-MC technique also needs to take into account the effects of removing/including systems that could be misclassified into the wrong  absorber type (e.g. a subDLA with a measured logN()$=20.1\pm0.2$ could be a DLA in reality, or conversely a logN()$=18.8\pm0.2$ absorber could be a subDLA). By design, this misclassification effect is naturally taken into account by the random variation in logN() during each resampling of the identified absorbers down to a column density threshold of logN()$=18.8$ (Section <ref>).  However, as the supplementary literature DLA sample from [Sánchez-Ramírez et al, 2016] does not contain subDLAs and covers a different redshift path; only absorbers along the associated redshift path of the XQ-100 sightlines are used in the subDLA computation, even if the supplementary literature DLA would be classified as a subDLA. We emphasize that even though all absorbers identified with logN()$\geq18.8$ are used in the B-MC simulation, absorbers with a randomly-varied logN() are still classified as subDLAs and DLAs based on the classical logN() cuts of $19.0\leq$ logN() $<20.3$ and logN() $\geq20.3$, respectively, within each mock generated survey. As already stated, the B-MC simulations were repeated using all identified absorbers with logN()$\geq18.8$ and those that were metal-selected to provide bounds on the quantities computed in this section. The subDLAs generated from these two sets of B-MC simulations will be further referred to as the full and metal-selected samples (abbreviated FS and MS), respectively.

§.§ The frequency of subDLAs

Using the B-MC resampling of the XQ-100 and supplementary DLA samples,  was computed using:

\begin{equation}
\label{eq:fNX}
f(N,X) dN dX=  \frac{m_{\rm abs}}{\Delta N \times \sum_{i}^{n_{\rm QSO}} \Delta X_{i}} dN dX,
\end{equation}

where $m_{\rm abs}$ is the total number of absorbers with column densities between $N-\Delta N/2$ and $N+\Delta N/2$ along the observed absorption distance $\Delta X_{i}$ of the $i^{th}$ QSO sightline, and $n_{\rm QSO}$ is the number of  QSO sightlines observed. The total redshift path covered by the QSO sightline ($\Delta X_{i}$) is computed as

\begin{equation}
\label{eq:zpath}
\Delta X(z)=\int_{z_{min}}^{z_{max}} (1+z)^2 [\Omega_M(1+z)^3 + \Omega_\Lambda]^{-0.5} dz,
\end{equation}

where $z_{min}$ and $z_{max}$ are the minimum and maximum redshift observed along the QSO path, while $\Omega_M$ and $\Omega_\Lambda$ are the matter and dark energy densities of the Universe observed at the current epoch. $z_{min}$ is set by the Lyman limit of the QSO, while $z_{max}$ is set as  $-5000$  to remove proximate systems that can be associated with the QSO [Ellison et al, 2002, Ellison et al, 2010, Berg et al, 2016, Perrotta et al, 2016].

Rather than using the standard bins of column density, which can hide or amplify features depending on the choice of bins, we elected to use a sliding logN() bin to quantify the median  in each column density bin. The sliding bin keeps a fixed bin width, but the centre of the bin shifts by a small amount in each step to capture the effects of individual data being included and excluded from the bin.  Using a logN() bin width of 0.4 dex, we computed the median, 1$\sigma$ and 2$\sigma$ percentiles of the  underlying  distribution of the 10 000 B-MC samples. For the subDLAs, the centre of the bin starts at logN()$=19.2$ (such that bin edge includes the minimum subDLA threshold of logN()$=19.0$), and is then shifted by 0.05 dex until the the centre of the bin reaches logN()$=20.1$ (such that bin edge includes logN()$=20.3$). The same is done for the DLAs, but using logN() bin centres from 20.5 to 21.8. Note that the sliding bins of both the DLA and subDLA samples do not pass the logN()$=20.3$ threshold due to the different redshift path coverage of the two samples. The B-MC computation of  was repeated using both the FS and MS absorbers from XQ-100.

[Column density distribution of subDLAs and DLAs]The column density distribution of subDLAs in the XQ-100 sample. The median   for the full sample of subDLAs (FS) is given by the solid line, while the dark and light shaded regions represent the $1\sigma$ and $2\sigma$ percentiles (respectively). The dashed line shows the median  for the metal-selected sample of subDLAs (MS). The confidence intervals for the MS selected sample are similar to those of the FS, but are not shown for clarity. Note that the dashed lines from both of the MS redshift bins overlap below logN()$\approx19.8$.  Both subDLA samples were split by the median redshift of the XQ-100 absorbers ($=3.2$), such that those in the lower redshift bin ($2.3\leq$  $<3.2$) are in dark purple while those at higher redshifts ($3.2\leq$  $<4.5$) are in light orange. The greyed region below logN()$=19.4$ represents the regime where incompleteness and false-positive effects could  significantly influence the determined  (Section <ref>). The square points show the data for <cit.> subDLAs in two similar redshifts ($1.5\leq$  $<3.1$ and $3.1\leq$  $<5$; light purple and dark orange squares, respectively). The inset panel show the median  curves of both the subDLA and DLA samples for reference.

Figure <ref> shows the resulting median column density distribution functions of the FS and MS subDLA samples (solid and dashed lines; respectively). The $1\sigma$ and $2\sigma$ errors are shown as the dark and light shaded regions (respectively) for the FS absorbers. Each sample was split into two redshift bins by the median redshift of the XQ-100 absorbers ($=3.2$);  for the lower and higher redshift bins are denoted by the purple and orange colours (respectively). The numerical values of the FS and MS  curves are tabulated in the appendix (Tables <ref> and <ref>; respectively). For reference, the full  range for subDLAs and DLAs is shown in the inset panel of Figure <ref> with the same notation.

There is a clear difference in the shape of  between all candidate subDLAs (FS) and those absorbers which have metal lines (MS), resulting in up to a $\approx0.4-0.6$ dex difference in   at low column densities. It is interesting to note that the difference between the FS and MS curves changes with logN(), and is strongest at the lowest column densities. This discrepancy could result from either incompleteness in the MS (i.e. lack of ability to detect metal lines at low column densities) or a high false positive rate in the FS [i.e. from blending of low logN()]. We emphasize again that using both the FS and MS likely bounds the true nature of  at these column densities. The derived  curves were compared with the values computed by <cit.> for subDLAs over approximately the same redshift range. The [Zafar et al, 2013] sample consisted of a search for subDLAs and DLAs in archival high resolution spectra taken with UVES in addition to a literature compilation. The identification of absorbers presented in [Zafar et al, 2013] started with an automated search for candidate subDLAs and DLAs characterized by saturated absorption over a large range of wavelength. The absorber was confirmed as a (sub)DLA if higher order Lyman series or metal lines were present; which is analogous to the methods used in determining the XQ-100 MS absorbers. It is therefore both encouraging and not surprising that  for the MS absorbers is roughly consistent with [Zafar et al, 2013] over the similar high redshift bin (3.1 $\leq z<$ 5.0; orange squares in Figure <ref>). The discrepancy seen between the low redshift data is likely a result of different redshift coverage, as nearly two thirds of the absorbers in the low redshift bin of <cit.> are at a redshift  $<2.5$. We note that requiring coverage Ly$\beta$ redward of $\lambda_{LL}$ to remove potential false-positives (see Section <ref>) does not significantly affect the resulting  curve generated for the 3.2 $\leq z<$ 4.5 redshift bin. Given the redshift range of the QSOs of the XQ-100 sample, there is insufficient spectral coverage of Ly$\beta$ to perform the same experiment within the lower redshift bin (2.3 $\leq z<$ 3.2).

To better quantify the significance of the difference in  between the MS and FS subDLAs, as well as search for any global evolution in  with redshift, we computed the best-fit power law to  (i.e. $=K$N$^\alpha$; where $K$ and $\alpha$ are the parameters of the fit). For each mock sample generated in the B-MC simulation, an unweighted least-squares fit to  was computed using the sliding bin data. The median and 1$\sigma$ percentile errors of $\alpha$ and $K$ of all 10 000 samples least squares fits are presented in Table <ref>. We point out that a single power law does not appear to be a good fit for  for the $3.2\leq$  $<4.5$ FS subDLAs. Comparing the power law fit parameters $\alpha$ and $K$ between the MS and FS datasets demonstrates that the power law fits between the two subDLA ($19.0\leq$ logN() $<20.3$) samples are inconsistent at the $\approx2\sigma$ level, demonstrating that using metal selection could have a profound effect on the interpreted .

Comparing the power law indices from each redshift bin, an insignificant evolution of $\alpha$ with redshift is detected at $\approx1.3\sigma$ and $\approx0.6\sigma$ significance for FS and MS subDLAs (respectively), and slight evolution for DLAs at $\approx1.5\sigma$. The increased number of logN() $\gtrsim19.7$ subDLAs at higher redshifts  is consistent with previous observational studies <cit.>. However, we note that the simulations predict a consistently higher  across all column densities at higher redshifts [Rahmati et al, 2013, Villaescusa-Navarro et al, 2018], whereas both XQ-100 datasets suggest this is not the case at logN() $\lesssim 19.7$.  In the FS subDLAs, the increased number of low redshift, low column density subDLAs could partly be explained by a high false-positive rate due to blending of low redshift Ly$\alpha$ with the Ly$\beta$ forest. However the MS also shows no redshift evolution in  for logN() $\lesssim 19.7$, which is approximately consistent with the lack of  redshift evolution seen in simulations when local ionizing sources are included [Rahmati et al, 2013].

Previous analyses comparing the power law indices of DLAs and subDLAs have indicated a flattening in  in the subDLA regime [Péroux et al, 2005, O'Meara et al, 2007, Guimarães et al, 2009, Zafar et al, 2013], which has been attributed to the loss of Ly$\alpha$ self-shielding in low N() subDLA systems [Zheng & Miralda-Escudé, 2002]. [O'Meara et al, 2007] noted that this flattening in  at low column densities is dominated by logN()$\leq19.3$ systems; an effect which was tentatively seen in their power law fits to  when only logN()$\geq19.3$ absorbers were included. We attempted to search for this flattening as well by restricting the subDLA sample only to absorbers with $19.3\leq$logN()$<20.3$ and recomputing the power law fit to limited column density range (see Table <ref> for fit parameters). The power law indices between the two logN() regimes are inconsistent at $1.3\sigma$ and $0.6\sigma$ significance for the FS and MS (respectively), and thus we do not find significant evidence of this flattening.

Power law fit parameters to $=KN^{\alpha}$

logN()          2c|$\alpha$     2clog($K$)

FS absorbers     MS absorbers    FS absorbers     MS absorbers

0pt3ex$\in$[19.0,20.3)     $\in$[2.3,3.2)     $-2.00_{-0.176}^{+0.176}$     $-1.46_{-0.196}^{+0.276}$     $18.64_{-3.369}^{+3.381}$     $7.74_{-5.293}^{+3.847}$

0pt4ex     $\in$[3.2,4.5)     $-1.37_{-0.119}^{+0.083}$     $-1.12_{-0.198}^{+0.186}$     $6.18_{-1.602}^{+2.294}$     $1.15_{-3.688}^{+3.752}$

0pt3ex$\in$[19.3,20.3)     $\in$[2.3,3.2)     $-1.94_{-0.399}^{+0.230}$     $-1.41_{-0.305}^{+0.337}$     $17.46_{-4.457}^{+7.831}$     $6.72_{-6.660}^{+6.050}$

0pt4ex     $\in$[3.2,4.5)     $-1.14_{-0.160}^{+0.260}$     $-0.90_{-0.293}^{+0.373}$     $1.68_{-5.255}^{+3.077}$     $-3.30_{-7.461}^{+5.769}$

0pt3ex$\in$[20.3,22.0)     $\in$[2.3,3.2)     $-2.10_{-0.098}^{+0.114}$     $-2.08_{-0.084}^{+0.062}$     $21.23_{-2.392}^{+2.058}$     $20.88_{-1.319}^{+1.713}$

0pt4ex     $\in$[3.2,4.5)     $-1.96_{-0.082}^{+0.131}$     $-1.95_{-0.116}^{+0.111}$     $18.52_{-2.729}^{+1.690}$     $18.20_{-2.308}^{+2.354}$

§.§.§ The line density of subDLAs

In addition to , the sightline density of absorbers () provides additional constraints on the nature of the absorbers.  is defined as:

\begin{equation}
\ell(X) = \int_{N_{min}}^{N_{max}} f(N,X) dN = \frac{m_{\rm abs}}{\sum_{i}^{n_{\rm QSO}} \Delta X_{i}}.
\end{equation}

Using the same sliding bin technique described above, we computed  in bins of $\Delta X=2.5$ with redshift steps of 0.01. Figure <ref> shows the resulting redshift evolution curves of  for the subDLAs (both FS and MS, blue solid and dashed lines) and DLAs <cit.>. The curves presented in Figure <ref> are available in Table <ref> (online only). Previous studies have disagreed on whether there is indeed any redshift evolution in  for subDLAs [O'Meara et al, 2007, Zafar et al, 2013]. However, any redshift evolution present in   can physically manifests as a change in either (or both) the  cross-sections of the absorbers or their number density <cit.>.

After applying a unweighted least squares linear fit to the redshift evolution of  for each of the 10 000 iterations in the MC-B simulation, we find the median slope of the redshift evolution is almost negligible; $-0.02_{-0.013}^{+0.011}$ for the FS subDLAs, $0.01_{-0.009}^{+0.006}$ for the MS subDLAs, and $0.01_{-0.009}^{+0.009}$ for both DLA samples. We point out that the evolution of  in the MS is shallower than the evolution detected by <cit.>. The difference in the evolution of the  curves of the MS and FS (both in terms of the slope and relative offset) results from the increased number of low column density FS absorbers seen in Figure <ref>.

The lack of evolution in  seen for DLAs has been suggested to be due to a balance between the decreasing density of the Universe and a decrease in the strength of ionizing UV background in order to maintain an approximately constant  [Prochaska et al, 2005, O'Meara et al, 2007]. Despite subDLAs being about two to four times more common than DLAs along a QSO sightline, the redshift evolution of both the subDLAs frequency and cross-sections appear to follow the evolution of the respective DLA quantities over these redshifts.

The redshift evolution of  as probed by DLAs (red) and subDLAs (blue) in the XQ-100 survey.  The solid line shows the median  evolution for the full sample of subDLAs (FS) while the dashed line represents when only the metal-selected subDLAs (MS) is used. The shaded regions show the $1\sigma$ and $2\sigma$ confidence intervals for the FS; the confidence intervals for the MS are approximately the same size, but are not shown for clarity. The values of  computed for the <cit.> subDLAs are shown as the grey squares.

§.§ The HI mass contribution of subDLAs

was estimated using:

\begin{equation}
\label{eq:OmegaInt}
\Omega_{\rm H \textsc{i}} (X) = \frac{H_{0} m_{H}}{c \rho_{0}} \int_{N_{min}}^{N_{max}} N f(N,X) dN,
\end{equation}

where $H_{0}$ is the Hubble constant at the current epoch, $c$ is the speed of light, $m_{H}$ is the atomic mass of hydrogen, and $\rho_{0}$ is the critical density of the Universe.

Following [Sánchez-Ramírez et al, 2016],  is computed in sliding bins of redshift. Rather than holding the bin width fixed in redshift space, [Sánchez-Ramírez et al, 2016] elected to use a fixed redshift path length $\Delta$X$=2.5$, with redshift bin steps of 0.01. For consistency, we elected to use the same sliding bin setup to reproduce the observed  for the XQ-100 DLAs. For the derivation of  (Equation <ref>), we start with the same $z_{min}$ and $z_{max}$ used in [Sánchez-Ramírez et al, 2016] for all the XQ-100 and supplementary QSO sightlines. However, due to the incompleteness of subDLA absorbers below $\lambda_{LL}$ (Section <ref>), we increased $z_{min}$ when necessary to remove all XQ-100 absorbers whose Ly$\alpha$ fell below $\lambda_{LL}$.

Similar to Figure <ref>, the top panel of Figure <ref> shows the redshift evolution of  for both the DLAs (red) and subDLAs (blue), using both the FS (solid line) and MS absorbers (dashed line; subDLAs only). The  curves presented in Figure <ref> are available in the Appendix (Table <ref>). Although there is a slight but significant difference in the magnitude of  for the subDLAs between the FS and MS (detected at the $\approx1\sigma$ level), no significant evolution in  with redshift is detected for subDLAs in the range $2.5\lesssim$  $\lesssim4.2$. We remind the reader that beyond  $\gtrsim4.2$  the redshift path coverage of the XQ-100 survey significantly drops <cit.>, so it is unclear if the drop in  beyond  $\gtrsim4.2$ for subDLAs is real.

The relatively constant  as a function of redshifts as traced by subDLAs is remarkably similar to the gradual decrease in  measured in DLAs over the same redshift range ($2.5\lesssim$  $\lesssim4.2$). The canonical picture that has emerged to explain the gradual evolution in  in DLAs is a constant replenishment of gas in star-forming galaxies to counteract the consumption of gas from star formation and ejection from galactic outflows [Noterdaeme et al, 2009, Noterdaeme et al, 2012, Zafar et al, 2013]. The similar lack of significant redshift evolution in the  statistics of both DLAs and subDLAs (,  and ) is suggestive that the same astrophysical processes are responsible in the host galaxies of these absorbers.

Previous cosmological simulations have demonstrated that only higher mass galaxies contribute significantly to  below  $\lesssim 3$ while lower mass galaxies dominate the contribution at  $\gtrsim 3$, which leads to a relatively constant  with redshift <cit.>. We speculate that the lack of significant evolution in  for both subDLAs and DLAs for $2.5\lesssim$  $\lesssim4.2$ in Figure <ref> could potentially suggest that the average subDLA and DLA absorbers are probing the same range of host galaxy halo masses. As such, this potential equality in halo mass range for subDLA and DLA host galaxies would be in tension with previous works which use the observed higher metallicity of subDLAs (relative to DLAs) as a proxy for a higher host galaxy mass <cit.>.

[ of XQ-100 subDLAs and DLAs]Top panel: The redshift evolution of  as probed by DLAs (red) and subDLAs (blue) in the XQ-100 survey.  The notation for the DLA and subDLA samples are the same as Figure <ref>. The values of  computed for subDLAs in <cit.> are shown as the grey squares. Bottom panels: The cumulative contribution of absorbers to the total  for a given maximum column density ().  is defined to be zero when only including absorbers with logN()$<19.0$, and has the value one when including all absorbers with logN()$\leq22.0$. The left and right panels in the bottom row show  for the respective full (FS) and metal-selected (MS) subDLAs, while the colours show the two redshift ranges $2.3\leq$  $<3.2$ (purple) and $3.2\leq$  $<4.5$ (orange). The vertical dashed line indicates the logN()$=20.3$ threshold that separates DLAs and subDLAs, while the greyed region below logN()$=19.4$ represents the regime where incompleteness and false positives could significantly affect the calculation (Section <ref>).

To estimate the relative mass contribution of subDLAs to the  budget, we computed the relative ratio of  from subDLAs to the total  (i.e. DLAs and subDLAs). To reduce the effects of the vastly different redshift sampling of the DLA and subDLA samples, we elected to use the median  curves (Figure <ref>) for the two redshift bins using the discrete limit of Equation <ref>, namely

\begin{equation}
\label{eq:OmegaFrac}
\Omega_{\rm H \textsc{i}} (X) = \frac{H_{0}  m_{H}}{c \rho_{0}} \frac{\sum_{j}^{m_{\rm abs}} N_{i}}{\sum_{j}^{n_{\rm QSO}} \Delta X_{i}}.
\end{equation}

The bottom two panels of Figure <ref> shows the cumulative fraction contribution to  from absorbers of a given  column density or lower (). Since we are only focussing on a range of N() column densities, we restrict  to have the value 0 for all absorbers that have logN()$<19.0$, and is 1 for all absorbers with logN() $<22.0$.

The bottom left and right panels of Figure <ref> show  for both the FS and MS (respectively), with the different colours representing the two redshift bins: $2.3\leq$  $<3.2$ (purple) and $3.2\leq$  $<4.5$ (orange). We note that the discontinuity in  at N() $\approx20.3$ is a result from the switch between the subDLA and DLA  curves. We remind the reader that this discontinuity is a result of the different redshift paths covered by the DLA and subDLA samples, resulting in a truncated  at logN() $=20.3$. Overall, the relative contribution of mass from subDLAs appears to be between $\approx12$–22% (for the MS and FS, respectively), which is consistent with estimates in the literature [Péroux et al, 2005, Guimarães et al, 2009, Noterdaeme et al, 2009, Zafar et al, 2013].  Despite the potentially high false-positive rate discussed in Section <ref> for subDLAs with logN()$\leq19.4$ at  $\lesssim3$, there is at most a 3 per cent contribution to the  mass contribution (i.e. the difference between the FS and MS curves in the left panel) from the potential false psoitives identified. We also note the lack of a redshift evolution in  for the MS, as previously seen by [Zafar et al, 2013]. Despite the increase in the number of logN()$\gtrsim19.7$ absorbers at higher redshifts (Figure <ref>), there is little effect on the  contribution from these systems.

§ CONCLUSION

We have presented a metal-blind catalogue of 242 candidate absorbers with  column densities of $18.8\leq$ logN() $<20.3$, identified solely using their  absorption, in the XQ-100 legacy survey. By comparing to archival UVES data, we demonstrated that X-Shooter has sufficient spectral resolution (FWHM resolution R$\sim$5400–8900) to identify potential subDLAs and measure their  column density accurately to within the typically quoted errors ($\pm0.2$ dex), which are primarily due to continuum placement. Using associated metal absorption lines as a method of removing false positive systems, we note that the 128 absorbers without confirmed metal lines were preferentially found to have logN() $\lesssim19.4$ and are at lower redshifts ( $\lesssim 3.0$). At these low column densities and redshifts, the corresponding prominent metal lines (Civ and Siiv) are shifted into the Ly$\alpha$ forest, Ly$\beta$ shifts beyond the Lyman limit of the QSO, and blending of the Ly$\alpha$ and Ly$\beta$ forests can mimic low column density subDLAs. Therefore there is a competing effect between incompleteness from using a metal-selected subsample and the true false-positive identification rate of $\lesssim49$ per cent of the sample. The false positive rate in XQ-100 is likely significant at low redshifts ( $\lesssim 3.5$) and low  column densities (logN() $\lesssim19.4$).

Using a Monte Carlo-Bootstrap simulation applied to the candidate absorber sample, we computed the column density distribution function (; Figure <ref>), line density (, Figure <ref>) and cosmic  mass density (; Figure <ref>) for subDLAs selected using the standard column density range $19.0\leq$ logN() $<20.3$. We compared the  computed for two different redshift bins ($2.3\leq$ $<3.2$ and $3.2\leq$ $<4.5$; split by the median absorber redshift of XQ-100) for both the metal-selected and full catalogue of subDLAs (dashed and solid lines in Figure <ref>, respectively). Irrespective of which subDLA catalogue is used, we find that solely high column density subDLAs (logN()$\gtrsim19.7$) are more frequently found at higher redshifts ($\geq3.2$), whereas all column densities are expected to increase in frequency in simulations [Rahmati et al, 2013, Villaescusa-Navarro et al, 2018]. However, the power law fits to the entire subDLA  curves suggest there is no  significant redshift evolution in  ($\lesssim1.3\sigma$ significance).

Using the full and metal-selected samples of subDLAs to bound the true value of  and ,  we demonstrated that there is little evolution in  and  as a function of redshift for DLAs; a result that is also seen for DLAs [Prochaska et al, 2005, Noterdaeme et al, 2012, Sánchez-Ramírez et al, 2016]. In comparison with simulations [Davé et al, 2013, Lagos et al, 2014] where the contribution from varying galaxy mass influences the nature of the redshift evolution  of , the lack of significant redshift evolution in  and  is suggestive that the subDLAs and DLAs could potentially trace host galaxies with roughly the same range of halo masses and that the evolution of their absorption cross-sections trace each other. Despite subDLAs being about two to four times more common than DLAs (Figure <ref>), subDLAs contribute a modest amount (between $\approx$10–20 per cent) to the total  budget over redshifts $2.3\leq$ $<4.5$ (Figure <ref>).

§ ACKNOWLEDGMENTS

We thank the anonymous referee for their comments on improving the clarity of this paper, and Céline Péroux for providing reduced 1D UVES data for comparison. SL was partially funded by UCh/VID project ENL18/18.

[Asplund et al, 2009]
Asplund M., Grevesse N., Sauval A. J., Scott P., 2009, , 47, 481

[Berg et al, 2017]
Berg T. A. M. et al., 2017, , 464, L56

[Berg et al, 2016]
Berg T. A. M. et al., 2016, , 463, 3021

[Bird et al, 2015]
Bird S., Haehnelt M., Neeleman M., Genel S., Vogelsberger M.,
Hernquist L., 2015, , 447, 1834

[Bird et al, 2014]
Bird S., Vogelsberger M., Haehnelt M., Sijacki D., Genel S., Torrey
P., Springel V., Hernquist L., 2014, , 445, 2313

[Braun, 2012]
Braun R., 2012, , 749, 87

[Christensen et al, 2017]
Christensen L. et al., 2017, , 608, A84

[Cooke et al, 2015]
Cooke R. J., Pettini M., Jorgenson R. A., 2015, , 800, 12

[Crighton et al, 2015]
Crighton N. H. M., Hennawi J. F., Simcoe R. A., Cooksey K. L., Murphy
M. T., Fumagalli M., Prochaska J. X., Shanks T., 2015, , 446, 18

[Davé et al, 2013]
Davé R., Katz N., Oppenheimer B. D., Kollmeier J. A., Weinberg
D. H., 2013, , 434, 2645

[Dessauges-Zavadsky et al, 2009]
Dessauges-Zavadsky M., Ellison S. L., Murphy M. T., 2009, , 396,

[Dessauges-Zavadsky et al, 2003]
Dessauges-Zavadsky M., Péroux C., Kim T.-S., D'Odorico S.,
McMahon R. G., 2003, , 345, 447

[Ellison et al, 2010]
Ellison S. L., Prochaska J. X., Hennawi J., Lopez S., Usher C.,
Wolfe A. M., Russell D. M., Benn C. R., 2010, , 406, 1435

[Ellison et al, 2001]
Ellison S. L., Ryan S. G., Prochaska J. X., 2001, , 326, 628

[Ellison et al, 2002]
Ellison S. L., Yan L., Hook I. M., Pettini M., Wall J. V., Shaver
P., 2002, , 383, 91

[Fox et al, 2007]
Fox A. J., Ledoux C., Petitjean P., Srianand R., 2007, , 473, 791

[Fumagalli et al, 2016]
Fumagalli M., O'Meara J. M., Prochaska J. X., 2016, , 455, 4100

[Guimarães et al, 2009]
Guimarães R., Petitjean P., de Carvalho R. R., Djorgovski S. G.,
Noterdaeme P., Castro S., Poppe P. C. D. R., Aghaee A., 2009, ,
508, 133

[Kulkarni et al, 2007]
Kulkarni V. P., Khare P., Péroux C., York D. G., Lauroesch J. T.,
Meiring J. D., 2007, , 661, 88

[Kulkarni et al, 2010]
Kulkarni V. P., Khare P., Som D., Meiring J., York D. G.,
Péroux C., Lauroesch J. T., 2010, , 15, 735

[Lagos et al, 2014]
Lagos C. d. P., Davis T. A., Lacey C. G., Zwaan M. A., Baugh C. M.,
Gonzalez- Perez V., Padilla N. D., 2014, , 443, 1002

[Lanzetta et al, 1995]
Lanzetta K. M., Wolfe A. M., Turnshek D. A., 1995, , 440, 435

[Lilly et al, 2013]
Lilly S. J., Carollo C. M., Pipino A., Renzini A., Peng Y., 2013,
, 772, 119

[López et al, 2016]
López S. et al., 2016, , 594, A91

[Martin et al, 2010]
Martin A. M., Papastergis E., Giovanelli R., Haynes M. P., Springob
C. M., Stierwalt S., 2010, , 723, 1359

[Neeleman et al, 2016]
Neeleman M., Prochaska J. X., Ribaudo J., Lehner N., Howk J. C.,
Rafelski M., Kanekar N., 2016, , 818, 113

[Noterdaeme et al, 2012]
Noterdaeme P. et al., 2012, , 547, L1

[Noterdaeme et al, 2009]
Noterdaeme P., Petitjean P., Ledoux C., Srianand R., 2009, , 505,

[O'Meara et al, 2007]
O'Meara J. M., Prochaska J. X., Burles S., Prochter G., Bernstein
R. A., Burgess K. M., 2007, , 656, 666

[Péroux et al, 2003]
Péroux C., Dessauges-Zavadsky M., D'Odorico S., Kim T.-S.,
McMahon R. G., 2003a, , 345, 480

[Péroux et al, 2005]
Péroux C., Dessauges-Zavadsky M., D'Odorico S., Sun Kim T.,
McMahon R. G., 2005, , 363, 479

[Péroux et al, 2003]
Péroux C., McMahon R. G., Storrie-Lombardi L. J., Irwin M. J.,
2003b, , 346, 1103

[Perrotta et al, 2016]
Perrotta S. et al., 2016, , 462, 3285

[Pettini et al, 1999]
Pettini M., Ellison S. L., Steidel C. C., Bowen D. V., 1999, , 510,

[Prochaska et al, 2005]
Prochaska J. X., Herbert-Fort S., Wolfe A. M., 2005, , 635, 123

[Prochaska & Wolfe, 2009]
Prochaska J. X., Wolfe A. M., 2009, , 696, 1543

[Quiret et al, 2016]
Quiret S. et al., 2016, , 458, 4074

[Rafelski et al, 2014]
Rafelski M., Neeleman M., Fumagalli M., Wolfe A. M., Prochaska J. X.,
2014, , 782, L29

[Rahmati et al, 2013]
Rahmati A., Pawlik A. H., Raičević M., Schaye J.,
2013a, , 430, 2427

[Rahmati et al, 2015]
Rahmati A., Schaye J., Bower R. G., Crain R. A., Furlong M.,
Schaller M., Theuns T., 2015, , 452, 2034

[Rahmati et al, 2013]
Rahmati A., Schaye J., Pawlik A. H., Raičević M.,
2013b, , 431, 2261

[Rao & Turnshek, 2000]
Rao S. M., Turnshek D. A., 2000, , 130, 1

[Rao et al, 2006]
Rao S. M., Turnshek D. A., Nestor D. B., 2006, , 636, 610

[Sánchez-Ramírez et al, 2016]
Sánchez-Ramírez R. et al., 2016, , 456, 4488

[Somerville & Davé, 2015]
Somerville R. S., Davé R., 2015, Annual Review of Astronomy and
Astrophysics, 53, 51

[Storrie-Lombardi & Wolfe, 2000]
Storrie-Lombardi L. J., Wolfe A. M., 2000, , 543, 552

[Vernet et al, 2011]
Vernet J. et al., 2011, , 536, A105

[Villaescusa-Navarro et al, 2018]
Villaescusa-Navarro F. et al., 2018, , 866, 135

[Wolfe et al, 2005]
Wolfe A. M., Gawiser E., Prochaska J. X., 2005, , 43, 861

[Wolfe et al, 1986]
Wolfe A. M., Turnshek D. A., Smith H. E., Cohen R. D., 1986, , 61,

[Zafar et al, 2013]
Zafar T., Péroux C., Popping A., Milliard B., Deharveng J.-M.,
Frank S., 2013, , 556, A141

[Zheng & Miralda-Escudé, 2002]
Zheng Z., Miralda-Escudé J., 2002, , 568, L71

[Zwaan et al, 2005]
Zwaan M. A., Meyer M. J., Staveley-Smith L., Webster R. L., 2005,
, 359, L30

§ LY$\ALPHA$ FITS TO XQ-100 ABSORBERS

[]Each pair of panels displays the fitted Ly$\alpha$(first and third columns) and Ly$\beta$ (second and fourth columns) profiles (dark red lines) and associated errors (dark red shaded region) for 14 of the identified candidate system. The light blue line and shaded region denote the full model spectrum and associated error combined from all absorbers within the entire sightline. The red text denotes the properties of the absorber [QSO sightline name and redshift in the corresponding Ly$\alpha$ panel, and logN() in Ly$\beta$ panel], while the vertical black text labels the redshift and Lyman series lines associated with additional absorption included in the model. The x-axis of each panel is scaled such that integrated region of the Ly series line's Voigt profile displayed encompasses 98 per cent of the absorbed flux.

§  AND  CURVES