## Abstract

The first black hole seeds, formed when the Universe was younger than ∼500 Myr, are recognized to play an important role for the growth of early (z ∼ 7) supermassive black holes. While progresses have been made in understanding their formation and growth, their observational signatures remain largely unexplored. As a result, no detection of such sources has been confirmed so far. Supported by numerical simulations, we present a novel photometric method to identify black hole seed candidates in deep multiwavelength surveys. We predict that these highly obscured sources are characterized by a steep spectrum in the infrared (1.6-4.5 μm), i.e. by very red colours. The method selects the only two objects with a robust X-ray detection found in the CANDELS/GOODS-S survey with a photometric redshift z ≳ 6. Fitting their infrared spectra only with a stellar component would require unrealistic star formation rates (≳2000 M yr−1). To date, the selected objects represent the most promising black hole seed candidates, possibly formed via the direct collapse black hole scenario, with predicted mass >105 M. While this result is based on the best photometric observations of high-z sources available to date, additional progress is expected from spectroscopic and deeper X-ray data. Upcoming observatories, like the JWST, will greatly expand the scope of this work.

## INTRODUCTION

The most distant objects with a spectroscopic redshift measurement detected in the Universe to date are extremely powerful gamma-ray bursts (z ∼ 8.2; Tanvir et al. 2009), and Lyman-break galaxies (z ∼ 11.1; Oesch et al. 2016). In the same period of the cosmic time, the reionization epoch, the first supermassive black holes (SMBHs) were already in place, but the seeds out of which these extremely massive objects were born have not been observed yet. Several observations (e.g. Mortlock et al. 2011; Wu et al. 2015) have detected the presence of accretion-powered objects with masses in excess of 109–10 M already at z ∼ 7, when the Universe was less than one billion years old. These observations are in tension with the standard theory of black hole growth, which would require, assuming Eddington-limited accretion, a longer time to produce these massive objects (Fan et al. 2006; Haiman 2013) from stellar-mass seeds, born out of the first population of stars (Pop III). Assuming that low angular momentum gas is always available for feeding, a black hole grows in mass exponentially, with an e-folding time ∼0.045 Gyr, for radiatively efficient accretion models. Starting from a stellar-mass seed (≲100 M), this process would require a constant accretion at the Eddington rate to produce a ∼109 M SMBH by z ∼ 7, a physically unlikely condition.

Two main solutions to the mass growth enigma have been proposed (Volonteri 2010). First, the black hole seeds from which the accretion started could have been as massive as ∼104–6 M (Lodato & Natarajan 2006; Devecchi & Volonteri 2009; Davies, Miller & Bellovary 2011). One possibility to build up massive black hole seeds at z ∼ 10–15 is through the direct collapse black hole (DCBH) scenario (Bromm & Loeb 2003; Shang, Bryan & Haiman 2010; Johnson et al. 2012). The collapse of a primordial atomic-cooling halo (with a virial temperature Tvir ≳ 104 K) may lead, in the presence of a strong flux of Lyman–Werner photons (energy hν = 11.2 − 13.6 eV) dissociating the H2 and thus preventing gas fragmentation, to the formation of DCBHs with a typical mass around 105 M (Ferrara et al. 2014). Secondly, the growth rates might not be capped by the Eddington value (Volonteri & Rees 2005; Alexander & Natarajan 2014; Madau, Haardt & Dotti 2014; Pacucci, Volonteri & Ferrara 2015a; Volonteri, Silk & Dubus 2015) $$\dot{M}_{\mathrm{Edd}} \equiv L_{\mathrm{Edd}}/(\eta c^2)$$, where LEdd ≈ 1.2 × 1038M erg s−1, M is the black hole mass in solar masses, η is the efficiency factor for mass–energy conversion and c is the speed of light. In highly obscured environments, radiation trapping could allow super-Eddington accretion rates to take place, dramatically speeding up the black hole growth.

At the present time, no detection of early SMBH progenitors has been confirmed. This may be due to their extreme faintness, but possibly also to the current lack of clear theoretical indications about their mass, host halo properties and typical accretion rates, resulting in large uncertainties in the prediction of their observational signatures. We proposed (Pallottini et al. 2015; see also Agarwal et al. 2015; Hartwig et al. 2015; Dijkstra, Gronke & Sobral 2016a; Smidt, Wiggins & Johnson 2016; Smith, Bromm & Loeb 2016; Visbal, Haiman & Bryan 2016) that a z ≈ 6.6 object named CR7, the brightest Ly α emitter discovered so far (Sobral et al. 2015), could be powered by a DCBH with an initial mass ∼105 M, associated with a standard stellar emission. This was motivated by the peculiarities of its spectrum, namely: (i) very strong Ly α and He ii lines emission, and (ii) absence of metal lines, within the detection threshold. This proposition needs deeper X-ray observations to be confirmed or rejected.

In previous works (Pacucci & Ferrara 2015; Pacucci et al. 2015b), employing a combination of radiation-hydrodynamic simulations and spectral synthesis codes, we computed the time-evolving spectrum emerging from the halo hosting a black hole seed. The bulk of the emission occurs in the observed infrared-submm (1–1000 μm) and X-ray (0.1–100 keV) bands. Here, we present a method to select DCBH candidates in ultradeep fields purely based on infrared photometry.

The outline of this paper is as follows. In Section 2, we describe the numerical setup of our simulations, while Section 3 presents the photometric tool that we developed. In Section 4, we show our results, including the identification of two DCBH candidates in CANDELS/GOODS-S. Finally, in Section 5, we discuss the caveats of this work and in Section 6 we provide some further discussion and the conclusions. Throughout, we adopt recent Planck cosmological parameters (Planck Collaboration XIII 2015): (Ωm, ΩΛ, Ωb, h, ns, σ8) = (0.32, 0.68, 0.05, 0.67, 0.96, 0.83).

## NUMERICAL IMPLEMENTATION

In this Section, we introduce the physical and numerical implementation of our simulations. The interested reader is referred to (Pacucci & Ferrara 2015, , radiation-hydrodynamic module) and to (Pacucci et al. 2015b, , spectral module) for a much more detailed description.

### General physical framework

A high-z black hole seed, with initial mass M(t = 0) in the range 104–5 M, is placed at the centre of a dark matter halo with primordial composition gas (H and He, with helium fraction YP ≈ 0.247, Planck Collaboration XIII 2015), total mass (baryonic and dark matter) Mh and virial temperature Tvir = Tvir(Mh, z) ∼ 104 K (Barkana & Loeb 2001):

(1)
$$T_{\rm vir} \approx 1.98\times 10^4\ \left(\frac{\mu }{0.6}\right) \left(\frac{M_{\rm h}}{10^8\ h^{-1} \ \mathrm{M_{\odot }} }\right)^{2/3} \left(\frac{1+z}{10}\right)\ {\rm K} .$$
In this formula, μ is the mean molecular weight and h is the reduced Hubble constant. At z ∼ 10, Tvir ∼ 104 K corresponds to Mh ∼ 108 M. In this work, we do not assume any specific relation between the initial black hole mass M(t = 0) and its host halo mass Mh. The mass of a black hole seed is assigned following the initial mass function derived in Ferrara et al. (2014). The black hole accretes mass from the inner regions of the host halo, with an accretion rate, $$\dot{M}$$, self-regulated by the combined effects of gravity, gas pressure, and radiation pressure.

Our radiation-hydrodynamic code (Pacucci & Ferrara 2015) evolves self-consistently the standard system of ideal, non-relativistic Euler's equations for a gas accreting radially on to the central black hole, assumed at rest and already formed at t = 0. We do not simulate the formation process of the massive black hole, but our implicit assumption is that it forms via the DCBH scenario, since: (i) the host halo is metal-free, and (ii) the surrounding gas has not been photo-evaporated by Pop III stellar emission. Other formation channels would not comply with these two assumptions.

The spatial domain of the simulations spans from 0.002 to 200 pc, largely encompassing the characteristic spatial scale for accretion, the Bondi radius, which, for a 105 M black hole, is ∼3 pc. The gas initially follows the isothermal density profile derived from the simulations in Latif et al. (2013a), where ρ(r) ∝ (r/a)−2 and a ∼ 2 pc is the core radius of the baryonic matter distribution. This density profile is the one resulting after the formation of a massive seed of mass ∼105 M at the halo centre (Latif et al. 2013b; Latif, Niemeyer & Schleicher 2014).

The value of the black hole mass M(t) changes with time, due to the accretion, modelled with the accretion rate $$\dot{M}_{\bullet }=4 \pi r^2 \rho |v|$$ (r is the radial coordinate, ρ the gas density, and v is the radial velocity, all computed at the inner boundary of the spatial domain). The accretion rate generates an emitted bolometric luminosity Lbol computed as:

(2)
$$L_{\rm bol} \equiv \eta c^2 \dot{M} ,$$
where η ≈ 0.1. The radiation pressure accelerates the gas via arad(r) = κ(ρ, T)Lbol(r)/(4πr2c), where the gas opacity κ(ρ, T), function of the gas density and temperature, includes Thomson (Begelman, Rossi & Armitage 2008) and bound-free terms.

### Observed DCBH spectrum

Our radiation-hydrodynamic code computes the frequency-integrated radiative transfer through the gas, taking into consideration the appropriate: (i) cooling and heating terms, (ii) matter-to-radiation coupling, and (iii) energy propagation through a two-stream approximation method. Full details are provided in Pacucci & Ferrara (2015). The frequency-dependent radiative transfer through the host halo is then performed using the spectral synthesis code cloudy (Ferland et al. 2013). This code computes the time-evolving spectrum emerging from the host halo (Pacucci et al. 2015b) employing as input data: (i) the spatial profiles for hydrogen number density and temperature, (ii) the source spectrum of the central object, and (iii) the bolometric luminosity of the source (computed self-consistently with equation 2). The source spectral energy distribution (SED), extended from far-infrared to hard X-ray, is a standard active galactic nuclei (AGN) spectrum, computed for a metal-free gas (Yue et al. 2013). The source spectrum depends on the black hole mass, and therefore evolves with time as M increases.

## PHOTOMETRY AND SED FITTING

In this Section, we describe the main photometric tool used in this work, along with the optical/infrared and X-ray data employed.

### The colour–colour plot

The shape of the spectrum emerging from a high-z dark matter halo hosting a DCBH depends on the black hole mass and on the absorbing hydrogen column density of the surrounding gas. If the halo is Compton-thick (NH ≳ 1.5 × 1024 cm−2) high-energy photons (≲10–50 keV) are largely absorbed and the DCBH may be invisible in the X-ray. None the less, sufficiently massive black hole seeds, M(t = 0) ≳ 104–5 M, are always visible in the infrared (Pacucci et al. 2015b), due to the reprocessing of high-energy radiation into lower energy photons. Moreover, the infrared spectrum is less dependent on the specifics of the accretion flow.

For these reasons, we devised a method to select DCBH candidates through infrared photometry. The computation of the predicted photometry for a DCBH follows from: (i) the time evolution of the spectrum emerging from the host halo, and (ii) the shapes of the photometric filters employed. Denoting by $${\cal F}(F)$$ the flux measured through the filter F, the generic colour index for the couple of filters A and B is defined as follows: $$A-B \equiv -2.5 \, \mathrm{Log}_{10} \left[ {\cal F}(A)/{\cal F}(B) \right]$$.

We employ three photometric filters to probe the infrared spectrum: the Hubble Space Telescope (HST) H (1.6 μm), the Spitzer IRAC1 (3.6 μm) and IRAC2 (4.5 μm). We construct a colour–colour diagram (Fig. 1) in which we plot both data points for objects in the CANDELS/GOODS-S multiwavelength survey (see Section 3.2 for a description of the data) and our predictions for the colours of DCBHs. All colours are observed quantities. The data points roughly fall on a line since the observed infrared spectra are to a good approximation power laws, with a broad range of slopes. The position in the colour–colour plot of DCBH seeds with initial masses in the range 104–5 M, computed via radiation-hydrodynamic simulations, is shown with filled circles whose colour depends on M(t = 0), as shown in the colour bar.

Figure 1.

Colour–colour diagram for the infrared filters H, IRAC1, and IRAC2. GOODS-S objects, brighter than the 27th magnitude in the H band (H < 27) and with 3.5 ≲ z ≲ 10, are shown with green points. Numerical predictions for the colours of DCBHs are shown, at z ∼ 7, with filled circles, whose colour depends on the initial mass of the seed (see the colour bar). Larger black hole masses are associated with redder spectra (i.e. more negative colours). All colours are observed quantities. An example of a photometric track (see Section 4.1) for a DCBH of initial mass ∼8 × 104 M is shown in orange. Its position has been shifted vertically to avoid information overload.

Figure 1.

Colour–colour diagram for the infrared filters H, IRAC1, and IRAC2. GOODS-S objects, brighter than the 27th magnitude in the H band (H < 27) and with 3.5 ≲ z ≲ 10, are shown with green points. Numerical predictions for the colours of DCBHs are shown, at z ∼ 7, with filled circles, whose colour depends on the initial mass of the seed (see the colour bar). Larger black hole masses are associated with redder spectra (i.e. more negative colours). All colours are observed quantities. An example of a photometric track (see Section 4.1) for a DCBH of initial mass ∼8 × 104 M is shown in orange. Its position has been shifted vertically to avoid information overload.

To compute the photometric points for DCBHs as a function of their mass, a total number of ∼100 simulations have been performed, initializing the black hole seed with masses in the range 104–5 M. For each simulation, the mass of the host halo has been modified in order to set the virial temperature Tvir ∼ 104 K for the redshift of choice (equation 1). Different sets of simulations have been performed at redshifts z = 13, z = 10, and z = 7. Since the shape of the SED at the wavelengths of our interest is roughly a power law, the effect of a redshift change on the position of the photometric points is small. Each simulation has been run for a time which varies with Mh and with M(t = 0), between ∼10 and ∼80 Myr, in order to reach a steady state accretion without depleting the gas reservoir. Sections 4.1 and 4.2 describe the dependence of DCBH colours on column density and on black hole mass.

Our method is perfectly adaptable to the JWST photometric system. The photometric filters F150W (centre wavelength λ0 = 1.5 μm), F356W0 = 3.6 μm) and F444W0 = 4.4 μm) are equivalents to the H, IRAC1, and IRAC2 filters, respectively, and may be used to implement our selection method with the JWST.

### The CANDELS/GOODS-S survey

The photometric data employed in this work are based on the GOODS-S field of the official CANDELS catalogues (Guo et al. 2013). Object selection has been performed in the H band of the near-infrared Wide Field Camera-3. The covering area of the GOODS-S survey is ∼170 arcmin2, to a mean 5σ depth of 27.5 mag in the H band. These imaging data include photometry over a wide range of wavelengths, from the U band (0.36 μm) to the IRAC4 band (7.9 μm). Importantly for our purposes, the catalogue includes very deep imaging with the IRAC instrument from the Spitzer Extended Deep Survey (Ashby et al. 2013), covering the CANDELS fields to a 3σ depth of 26 AB magnitude at both 3.6 μm (IRAC1) and 4.5 μm (IRAC2). Overall, the detection of objects in the GOODS-S field is realistic up to 27.5–28.0 mag in the H band, at ∼90 per cent completeness limit (Grazian et al. 2015).

For the photometric data that we employed in this work, each source has been visually inspected. The total number of sources detected in the GOODS-S field is 34 930. In our redshift range of interest, 3.5 ≲ z ≲ 10 (see Fig. 1), the total number of objects is 2037, while at high redshifts (6 ≲ z ≲ 10) is 97. In our analysis we included only GOODS-S objects with a precise value (i.e. not upper limits) for the three magnitudes H, IRAC1, and IRAC2. In total, there are 2272 spectroscopic redshifts of good quality in the GOODS-S field. The number of spectroscopic redshifts available at z > 3.5 is 173 (25 at z > 6.5). For sources lacking spectroscopic information, photometric redshifts were computed by optimally combining six different photometric redshifts (Dahlen et al. 2013). For a broader discussion on the uncertainties in the redshift determination, see Section 5.2.

### X-ray detected objects in the GOODS-S survey

In a recent work (Giallongo et al. 2015), 22 faint AGN candidates were identified in the GOODS-S field as having an X-ray counterpart (but see Weigel et al. 2015 and Cappelluti et al. 2015). The selection method is the following one: high-z galaxies are selected in the H band down to very faint levels (H ≤ 27) using reliable photometric redshifts. At z > 4, this corresponds to a selection criterion based on the galaxy rest-frame UV flux. AGN candidates are picked up from this parent sample if they show X-ray fluxes above a threshold of FX ∼ 1.5 × 10−17 erg cm−2 s−1 in the soft energy band (0.5–2 keV).

For three of these sources the photometric redshift is z ≳ 6: object 33160 at z ≈ 6.06, object 29323 at z ≈ 9.73, and object 28476 at z ≈ 6.26. H-band images of the first two sources, with X-ray contours from the Chandra 7 Ms field, are shown in Fig. 2. As we will show in Section 4, these two sources are our DCBH candidates. The photometric redshift of 29323 is more uncertain, but its probability distribution function suggests that z ≳ 6 with high confidence (see Fig. 2). Objects 33160 and 29323, both at z ≳ 6, are to be considered robust detections, being also included in the Xue et al. (2011) catalogue, with identification numbers 85 and 156, respectively. Moreover, the source 29323 is also detected in Cappelluti et al. (2015), but these authors suggest that the determination of its photometric redshift could be affected by potential artefacts in the SED.

Figure 2.

H-band images with X-ray contours (1–4 keV) of the two GOODS-S sources selected as DCBH candidates. X-ray contours are from the Chandra 7 Ms field, in square root scale, from 2σ to 10σ. H-band counterparts are circled. On the right, probability distribution functions (taken from Giallongo et al. 2015) for the photometric redshifts of the same objects, along with their optical/infrared spectra, fitted from their photometry.

Figure 2.

H-band images with X-ray contours (1–4 keV) of the two GOODS-S sources selected as DCBH candidates. X-ray contours are from the Chandra 7 Ms field, in square root scale, from 2σ to 10σ. H-band counterparts are circled. On the right, probability distribution functions (taken from Giallongo et al. 2015) for the photometric redshifts of the same objects, along with their optical/infrared spectra, fitted from their photometry.

Instead, object 28476 is a possible X-ray detection: for this source there are 15 detected counts in the soft-band (0.5–2 keV), compared to ∼5 background counts within a radius of 3 pixel, corresponding to a probability of spurious detection of ∼10−4. Anyway, the association between this X-ray source and its optical counterpart is more problematic. A full description, including coordinates, of all the 22 X-ray detected sources is provided in Giallongo et al. (2015).

## RESULTS

In this Section, we describe how the photometry of DCBHs changes with the column density of the host halo and with the black hole mass. Then, we present the black hole seed candidates selected in the CANDELS/GOODS-S field.

### Column-density dependence of DCBH colours

The initial black hole mass M(t = 0) is set for any given simulation in the range 104–5 M. At some redshift z, the gas mass Mg of the host halo is also set from equation (1), which in turn translates into an initial column density NH, with values generally in the Compton-thick range. The column density controls how photons with different frequencies are absorbed or transmitted. As long as the column density is large (NH ≳ 1.5 × 1024 cm−2), the photometry of DCBHs does not vary largely. As the black hole grows with time, its bolometric luminosity increases as well, but this effect does not modify the colours, since the increase in flux is similar in each photometric band. Once the column density falls below the Compton-thick limit, the evolution is very rapid and short-lived (≲1–2 Myr; Pacucci & Ferrara 2015): the DCBH performs a short photometric track in the colour–colour plot (see Fig. 1). The filled circles shown in Fig. 1 are computed at z = 7 and are averages, weighted over time, of the photometric track performed in each simulation.

### Mass dependence of DCBH colours

The trend of DCBH colours with the black hole mass (Fig. 1) may be simply explained with the following argument. The column density NH is proportional to the gas mass Mg and to the radius of the baryonic matter distribution Rg in the following way: $$N_H \propto M_{\rm g}/R_{\rm g}^2$$. Defining the concentration parameter of the gas distribution as $${\cal C} \equiv M_{\rm g}/R_{\rm g} \propto M_{\bullet }$$ (larger black holes are more efficient and rapid in concentrating the gas around them), we derive that the column density of haloes hosting more massive black holes is larger: $$N_H \propto M_{\rm g}/R_{\rm g}^2 \propto {\cal C}/R_{\rm g} \propto M_{\bullet }$$. Here, we are assuming that, on sufficiently large spatial scales (≫RB), the characteristic radius of the baryonic matter distribution of the host halo is not modified by the mass of the black hole seed that forms at its centre. The gas around larger black holes is more concentrated and then more efficient in reprocessing the radiation to lower energies. Consequently, the steepness of the infrared spectrum of massive black hole seeds is directly proportional to the black hole mass.

### Identification of DCBH candidates in GOODS-S

DCBHs with larger masses are associated with steeper (i.e. redder) infrared spectra. These massive objects extend down to a region in the colour–colour diagram where there is a relative paucity of GOODS-S objects (see Fig. 1). We pick three GOODS-S objects characterized by different observational features and compare in Fig. 3 their SEDs with the one computed for a ∼5 × 106 M black hole, grown from a DCBH with initial mass around ∼105 M. The SED fitting has been executed with a standard stellar population model (Bruzual & Charlot 2003). As described in Castellano et al. (2014) and in Grazian et al. (2015), we have included the contribution from nebular emission lines following Schaerer & de Barros (2009).

Figure 3.

Comparison between the stellar SEDs of three GOODS-S objects with the computed SED of a ∼5 × 106 M black hole, born out of a DCBH with initial mass around ∼105 M. The three photometric bands employed in our work are shown as shaded regions. Filled circles (observations) and stars (numerical simulations) show the magnitudes in the three filters, with error bars. Objects 29323 and 14800 have X-ray counterparts (i.e. they are likely associated with a black hole), while 24021 has not (i.e. it is likely a normal galaxy). Moreover, object 29323 is characterized by very negative colours (i.e. its infrared SED is very steep, as we predict for DCBHs), while objects 14800 and 24021 are not. The steepness of the SED and the infrared magnitudes for the object 29323 are well fitted by the spectrum predicted for a ∼5 × 106 M black hole. In the computed SED for a DCBH, the He ii line (0.164 μm rest-frame) is visible and it is marginally inside the H-band at z ≈ 9.7.

Figure 3.

Comparison between the stellar SEDs of three GOODS-S objects with the computed SED of a ∼5 × 106 M black hole, born out of a DCBH with initial mass around ∼105 M. The three photometric bands employed in our work are shown as shaded regions. Filled circles (observations) and stars (numerical simulations) show the magnitudes in the three filters, with error bars. Objects 29323 and 14800 have X-ray counterparts (i.e. they are likely associated with a black hole), while 24021 has not (i.e. it is likely a normal galaxy). Moreover, object 29323 is characterized by very negative colours (i.e. its infrared SED is very steep, as we predict for DCBHs), while objects 14800 and 24021 are not. The steepness of the SED and the infrared magnitudes for the object 29323 are well fitted by the spectrum predicted for a ∼5 × 106 M black hole. In the computed SED for a DCBH, the He ii line (0.164 μm rest-frame) is visible and it is marginally inside the H-band at z ≈ 9.7.

Objects 29323 (red line) and 14800 (blue line) are objects with an X-ray counterpart whose colours are very negative (29323) or close to zero (14800). Object 24021 (green line) is instead an object without an observed X-ray emission, whose colours are close to zero. The predicted infrared SED for a DCBH is redder than the ones observed for normal galaxies (e.g. 24021, in which radiation comes predominantly from stars) and AGNs (e.g. 14800, in which accretion-powered radiation co-exists with the stellar one). Current theoretical models suggest that most of the high-z black hole growth occurred into heavily obscured hosts (see e.g. Comastri et al. 2015), providing a theoretical ground for the redness of the infrared SED of DCBHs.

Objects in the GOODS-S field may be detected up to H ≈ 28 with a completeness ∼90 per cent. From our simulations, we conclude that DCBH seeds with initial masses below ∼6 × 104 M cannot be detected in current surveys, assuming realistic growth rates. Then, we predict to observe DCBH candidates with initial masses ≳6 × 104 M, i.e. with colours IRAC1 − H ≲ −1.8 and IRAC2 − H ≲ −1.8 (see Fig. 1). The infrared SEDs of objects within this area of the colour–colour diagram are very red: they are promising DCBH candidates.

The two GOODS-S objects with a robust X-ray detection at z ≳ 6 (see Section 3.3) fall within our selection region for DCBH candidates, as shown in Fig. 4, which is the main result of this work. At z ≳ 6, the photometry of these sources is very similar to the one predicted for DCBHs, characterized by a very red infrared SED, caused by large absorbing column densities. At lower redshifts (z ≲ 6), these pristine spectra are likely modified by several events, like star formation and large-scale outflows, which decrease the column density of the host halo. Moreover, also the pollution of the gas by metals may contribute to the modification of the SED. In Section 4.4, we discuss other kind of sources possibly mimicking DCBH colours. Conservatively, we have not included in Fig. 4 the position of the z ≳ 6 object 28476, which also falls in our selection region, but whose X-ray detection is less robust.

Figure 4.

Colour–colour diagram for objects at z ≳ 6. The yellow-shaded area indicates the region in the diagram where we expect SEDs compatible with the ones predicted for DCBHs, and observable with current surveys. Green points are GOODS-S objects without an X-ray counterpart, while red stars are GOODS-S objects detected in the X-ray and likely powered by accretion on to a collapsed object. Filled circles are, as in Fig. 1, DCBH simulations. In addition, we report the position in the colour–colour plot of several z > 6 sources. ULAS J1120+0641 (Mortlock et al. 2011; Barnett et al. 2015) is the most distant-known QSO (z ≈ 7.1). SDSS J1148+5251 (Fan et al. 2003) and SDSS J1048+4637 (Maiolino et al. 2004) are two highly obscured and dusty QSOs at z ≈ 6.4 and z ≈ 6.3, respectively. GN-z11 (Oesch et al. 2016) is the highest redshift galaxy discovered to date (z ≈ 11.1).

Figure 4.

Colour–colour diagram for objects at z ≳ 6. The yellow-shaded area indicates the region in the diagram where we expect SEDs compatible with the ones predicted for DCBHs, and observable with current surveys. Green points are GOODS-S objects without an X-ray counterpart, while red stars are GOODS-S objects detected in the X-ray and likely powered by accretion on to a collapsed object. Filled circles are, as in Fig. 1, DCBH simulations. In addition, we report the position in the colour–colour plot of several z > 6 sources. ULAS J1120+0641 (Mortlock et al. 2011; Barnett et al. 2015) is the most distant-known QSO (z ≈ 7.1). SDSS J1148+5251 (Fan et al. 2003) and SDSS J1048+4637 (Maiolino et al. 2004) are two highly obscured and dusty QSOs at z ≈ 6.4 and z ≈ 6.3, respectively. GN-z11 (Oesch et al. 2016) is the highest redshift galaxy discovered to date (z ≈ 11.1).

In Fig. 4, we also show the position of other high-z objects, to compare it with the locations of our DCBH candidates. ULAS J1120+0641 (Mortlock et al. 2011; Barnett et al. 2015) is the most distant known Quasi-Stellar Object (QSO, z ≈ 7.1), powered by an SMBH with an estimated mass ∼109 M. SDSS J1148+5251 (Fan et al. 2003) and SDSS J1048+4637 (Maiolino et al. 2004) are two highly obscured and dusty QSOs at z ≈ 6.4 and z ≈ 6.3, respectively. GN-z11 (Oesch et al. 2016) is the highest redshift galaxy discovered to date (z ≈ 11.1). The position of these objects in the colour–colour diagram is outside the region of detectability for sources powered by a DCBH: their H-band magnitudes would be too faint to be currently observed. Consequently, our photometric selection method seems to separate the spectra of pristine DCBHs from those of other contaminant high-z sources.

For the couple of z ≳ 6 X-ray detected sources in the CANDELS/GOODS-S field, we estimate a black hole mass in excess of ∼106 M. For the source 33160, we are able to reproduce simultaneously with simulations both its H-band magnitude, H = 25.9, and its 2–10 keV X-ray luminosity emerging from the host halo, Log10(Lx) = 43.65. Since we assume the photometric redshift reported in Giallongo et al. (2015), also its X-ray flux is automatically reproduced. Source 29323 at z ≈ 9.73 has a X-ray luminosity slightly above our predicted range (corresponding to a maximum black hole mass of ∼5 × 106 M), while we successfully reproduce its H-band magnitude (see Fig. 3). If this object is instead at redshift z ≲ 6, its X-ray luminosity would be downgraded to Log10(Lx) ≲ 43.5, well inside our predicted range. More generally, the X-ray detected objects in CANDELS/GOODS-S that are selected by our method (at any redshift) are compatible with the X-ray luminosities (2–10 keV) predicted for a population of DCBHs with initial masses in the range 6 × 104–105 M. Within this ensemble of sources, more negative colours are associated with a larger X-ray luminosity (i.e. with a larger black hole mass, as in Fig. 1).

In the simulated DCBH spectrum of the object 29323 in Fig. 3, the He ii line (0.164 μm rest-frame) is visible and it is marginally inside the H band at z ≈ 9.7. This emission line is an important indicator of DCBH activity (see e.g. Pallottini et al. 2015) and its detection, with future spectroscopic observations, may be a discriminant in the determination of the real nature of this high-z source.

### Star formation rates in the X-ray detected sample

The photometry of all the X-ray detected objects in GOODS-S has been fitted by a stellar SED-fitting model (Bruzual & Charlot 2003). The SED fitting assumes that all the luminosity is produced by stars, without any contribution from accreting objects.

Fig. 5 shows the distribution of star formation rates (SFRs) of the X-ray detected objects in GOODS-S (excluding object 28476, whose X-ray detection is still debated).

Figure 5.

Distribution of star formation rates of the X-ray detected objects in GOODS-S (excluding object 28476, whose X-ray detection is still debated). Each column of the histogram is coloured according to the average (IRAC1–H) colour index of the objects falling within the SFR bin, as shown in the colour bar. SFRs for the Milky Way (Robitaille & Whitney 2010), a typical ULIRG (Le Floc'h et al. 2005) and for a massive maximum-starburst galaxy at z = 6.34 (Riechers et al. 2013) are shown for comparison. The fraction of objects in each bin, at any redshift, falling into our selection region is shown as a red line, on the right vertical axis. The two DCBH candidates selected by our method are shown with yellow stars.

Figure 5.

Distribution of star formation rates of the X-ray detected objects in GOODS-S (excluding object 28476, whose X-ray detection is still debated). Each column of the histogram is coloured according to the average (IRAC1–H) colour index of the objects falling within the SFR bin, as shown in the colour bar. SFRs for the Milky Way (Robitaille & Whitney 2010), a typical ULIRG (Le Floc'h et al. 2005) and for a massive maximum-starburst galaxy at z = 6.34 (Riechers et al. 2013) are shown for comparison. The fraction of objects in each bin, at any redshift, falling into our selection region is shown as a red line, on the right vertical axis. The two DCBH candidates selected by our method are shown with yellow stars.

An increase in the SFR corresponds to a more negative value of the average (IRAC1–H) colour index, as the colour bar shows. As expected, a redder photometry is in general due to a larger dust extinction, which requires star formation to be active in the host galaxy. None the less, fitting the photometry of extremely red spectra, as the ones of our two DCBH candidates (shown as yellow stars) requires SFRs comparable to or larger than the highest ever measured, for a massive maximum-starburst galaxy at z = 6.34 (Riechers et al. 2013, see also Barger et al. 2014 where the authors find a characteristic maximum SFR of ∼2000 M yr−1). The computed SFR for object 29323 is ∼2 times larger than the maximum SFR ever measured, ∼20 times larger than the typical SFR of Ultra-Luminous Infrared Galaxies (ULIRGs) (Le Floc'h et al. 2005) and ∼5000 times the SFR of our Galaxy (Robitaille & Whitney 2010). Instead of evoking extremely large and unrealistic SFRs, we propose that these objects host a central DCBH (see Fig. 3 where the photometry of object 29323 is very well fitted by our DCBH model). The redness of their spectra is explained by the large, absorbing gas column densities in the host halo implicit in the DCBH evolution scenario.

A final remark: at low redshifts and for mildly negative colour indexes, the colours of DCBHs may be mimicked by star-forming galaxies. Fig. 5 shows with a red line the fraction of X-ray detected objects at any redshift that falls in our selection region. It is likely that the redness of some of these spectra may be genuinely explained by a burst of star formation underway in the host galaxy. At high redshifts and for very negative colour indexes (IRAC1–H ≲ − 2.5), the fraction of objects selected by our method is close or equal to unity and the required SFRs become unrealistic. We suggest that these objects are most likely high-z black hole seeds.

## CAVEATS

In this Section, we revise some of the most important caveats of the current work, related to the assumptions of our radiation-hydrodynamic code and to the uncertainties in the photometric redshift estimate of the X-ray detected sources.

### Assumptions of the radiation-hydrodynamic code

As described in Section 2, our radiation-hydrodynamic code simulates the accretion on to an already-formed DCBH with some assumptions on the geometry of the flow and on the environment hosting the compact object.

On sufficiently large scales (comparable with the Bondi radius of the DCBH), the geometry is assumed to be spherical, with gas moving radially without angular momentum. As discussed in Pacucci et al. (2015a), this assumption is supported by several studies (e.g. Choi, Shlosman & Begelman 2013, 2015) on the triaxiality of primordial haloes, which should allow the gas to flow radially well beyond the centrifugal radius of the system.

The accretion rate is computed self-consistently from the density and velocity fields of the gas at the inner boundary of our simulation domain. While some accretion models allow for a fraction of the inflowing mass to be lost on small scales during the accretion process (see e.g. Begelman 2012; Coughlin & Begelman 2014), our simulations assume that the gas reaching the inner boundary of our spatial domain is entirely accreted on to the black hole. The bolometric luminosity emerging from the accretion disc is then computed consistently from the accretion rate (see equation 2).

Star formation is not active in our simulations. Consequently, the emerging spectrum is computed for a metal-free environment. The presence of a stellar component in the host halo (Pop II or Pop III) may modify the emerging spectrum in the optical/UV region, while we checked that the infrared part, upon which our photometric method is based, is totally unaffected. Note, however, that, in the classical DCBH scenario, fragmentation, and hence star formation, is heavily suppressed by the external UV background.

### Redshift estimates for DCBH candidates

The redshift estimates for sources 33160 and 29323, selected by our method as DCBH candidates (see Fig. 4), are photometric and subject to uncertainties (see Fig. 2, right-hand panel).

As already mentioned, the photometric redshifts of the GOODS-S galaxies have been computed with the unconventional technique described in Dahlen et al. (2013). In summary, it is an optimal Bayesian combination of six different probability distribution functions of the photometric redshift computed by six different groups within the CANDELS Team. For the object 29323, the individual probability distribution functions are centred at z ∼ 6.5 and at z ∼ 9.7, and only one solution out of six gives z ≈ 10. The 68 per cent (95 per cent) confidence level region is between zmin = 5.46 and zmax = 9.86 (zmin = 2.36, zmax = 9.97). From the final probability distribution function for this object, we can conclude that the redshift is unconstrained at z ≳ 5.5 at 1σ (or 2.4 at 2σ), due to the power-law SED without significant breaks at the bluer wavelengths. Given the extreme steepness of the infrared SED of object 29323 (see Fig. 3), however, it would still be selected as a DCBH candidate down to z ∼ 2.

Other photometric redshift catalogues on the same field have been published (e.g. by the 3D-HST collaboration, Brammer et al. 2012, and by Hsu et al. 2014). None the less, we decided to adopt the CANDELS photometric redshifts for several reasons. First, the 3D-HST photometric redshifts are very accurate for bright objects with strong emission lines. For faint galaxies such as the ones in this paper, we prefer to use the CANDELS photometric redshifts, which are based on point spread function deconvolved photometry, which is quite robust also for faint and contaminated objects (e.g. object 29323). Moreover, the CANDELS photometric redshifts have the unique feature of averaging the differences between different photometric redshift codes, which could have their own biases. To conclude, the accuracy in the determination of the photometric redshift at z ≳ 4 relies mainly on the absorption of the intergalactic medium, which is an accurate feature independently of the adopted templates. The photometric redshifts of our DCBH candidates are z(29323) ≈ 5.3 and z(33160) ≈ 4.9 in the 3D-HST catalogue and z(29323) ≈ 4.7 and z(33160) ≈ 3.4 in the Hsu et al. (2014) catalogue.

## DISCUSSION AND CONCLUSIONS

Supported by radiation-hydrodynamic simulations, we have devised a novel method to identify black hole seed candidates at z ≳ 6: at high-z the infrared SED of pristine DCBHs is predicted to be significantly red. At lower redshifts, intervening processes (e.g. star formation, metal pollution, outflows) make the selection more uncertain. By applying our method to state-of-the-art photometric observations of the CANDELS/GOODS-S field, we select the only two sources at z ≳ 6 with a robust X-ray detection. The SFRs required to mimic their extremely red spectra are unrealistic (≳2000 M yr−1), being comparable to or larger than the highest SFR ever measured, in a maximum-starburst galaxy. These objects represent the most stringent observational identification of black hole seed candidates, likely formed as DCBHs, so far obtained. The possible presence of markers of DCBH activity, such as the predicted He ii emission line, needs to be investigated with spectroscopic follow-up observations.

The extension of this work to more infrared photometric bands may prove to be useful in distinguishing the different classes of objects observed, especially at lower redshifts. However, our two-dimensional colour–colour plot aims at conveying the simple and clear idea that, in looking for pristine DCBH candidates, observational efforts should be focused on objects with very steep infrared spectra.

This work extends previous efforts in understanding the observational features of DCBHs. Our previous studies (Pacucci et al. 2015b; Pallottini et al. 2015) investigated the spectrum emitted from these sources, foreseeing the possibility that they could be observed by the JWST. We suggested that CR7, the brightest Ly α emitter discovered to date, could be powered by a typical DCBH of initial mass ∼105 M. Recently, CR7 has been the subject of great interest in the community, with a large wealth of studies investigating its nature (e.g. Agarwal et al. 2015; Hartwig et al. 2015; Dijkstra et al. 2016a; Smidt et al. 2016; Smith et al. 2016; Visbal et al. 2016), most of them confirming the possibility that this source may host a DCBH. The efforts to observe these objects have recently spread to other regions of the electromagnetic spectrum. For instance, Dijkstra et al. (2016b) suggested that the large neutral hydrogen column densities of primordial gas at ∼104 K with low molecular abundance (all key aspects of the DCBH scenario) give rise to stimulated fine structure emission at ∼3 cm in the rest-frame. Detecting this signal is challenging, but would provide direct evidence for the DCBH scenario.

Our current work has, of course, some limitations. First, high-z sources with a robust X-ray detection in the GOODS-S field are few: making a precise assessment on the robustness of our selection criterion is beyond current reach. Secondly, only deep high/medium-resolution spectroscopy may conclusively confirm the DCBH nature of these objects. A spectroscopic analysis of these sources is beyond the capabilities of current observatories (e.g. HST, VLT, Keck) because they are too faint. With the start of JWST operations, a large wealth of infrared spectroscopic data (extended above the ∼ 2.5 μm limit of ground telescopes) will be available, disclosing to our eyes the first glimpses of light in the Universe. Currently, we probe only the most massive and luminous black holes, the peak of their mass distribution. JWST will mark a breakthrough in this field, by detecting the light from the most distant stars and accreting black holes, probing the mass range (104–5 M) of the first black hole seeds, if they were formed via the DCBH scenario. Our work establishes a solid theoretical framework to interpret these data with the aim of finding the first black holes in the Universe.

We thank N. Cappelluti, C. Feruglio, S. Gallerani, A. Loeb, and A. Mesinger for helpful discussions and comments. This work is based on observations taken by the CANDELS Multi-Cycle Treasury Program with the NASA/ESA HST. Observations were also carried out using the VLT, the ESO Paranal Observatory, the Spitzer Space Telescope and the Chandra X-ray Observatory.

## REFERENCES

Agarwal
B.
Johnson
J. L.
Zackrisson
E.
Labbe
I.
van den Bosch
F. C.
Natarajan
P.
Khochfar
S.
2015
preprint (arXiv:1510.01733)
Alexander
T.
Natarajan
P.
2014
Science

345
1330
Ashby
M. L. N.
et al
2013
ApJ

769
80
Barger
A. J.
et al
2014
ApJ

784
9
Barkana
R.
Loeb
A.
2001
Phys. Rep.

349
125
Barnett
R.
et al
2015
A&A

575
A31
Begelman
M. C.
2012
MNRAS

420
2912
Begelman
M. C.
Rossi
E. M.
Armitage
P. J.
2008
MNRAS

387
1649
Bromm
V.
Loeb
A.
2003
ApJ

596
34
Brammer
G. B.
et al
2012
ApJS

200
13
Bruzual
G.
Charlot
S.
2003
MNRAS

344
1000
Cappelluti
N.
et al
2015
preprint (arXiv:1512.00510)
Castellano
M.
et al
2014
A&A

566
A19
Choi
J.-H.
Shlosman
I.
Begelman
M. C.
2013
ApJ

774
149
Choi
J.-H.
Shlosman
I.
Begelman
M. C.
2015
MNRAS

450
4411
Comastri
A.
Gilli
R.
Marconi
A.
Risaliti
G.
Salvati
M.
2015
A&A

574
L10
Coughlin
E. R.
Begelman
M. C.
2014
ApJ

781
82
Dahlen
T.
et al
2013
ApJ

775
93
Davies
M. B.
Miller
M. C.
Bellovary
J. M.
2011
ApJ

740
L42
Devecchi
B.
Volonteri
M.
2009
ApJ

694
302
Dijkstra
M.
Gronke
M.
Sobral
D.
2016a
preprint (arXiv:1602.07695)
Dijkstra
M.
Sethi
S.
Loeb
A.
2016b
ApJ

820
10
Fan
X.
et al
2003
AJ

125
1649
Fan
X.
et al
2006
AJ

131
1203
Ferland
G. J.
et al
2013
Rev. Mex. Astron. Astrofis.

49
137
Ferrara
A.
S.
Yue
B.
Schleicher
D.
2014
MNRAS

443
2410
Giallongo
E.
et al
2015
A&A

578
A83
Grazian
A.
et al
2015
A&A

575
A96
Guo
Y.
et al
2013
ApJS

207
24
Haiman
Z.
2013
Wiklind
T.
Mobasher
B.
Bromm
V.
Astrophysics and Space Science Library, Vol. 396, The First Galaxies

Springer-Verlag
Berlin
293
Hartwig
T.
et al
2015
preprint (arXiv:1512.01111)
Hsu
L.-T.
et al
2014
ApJ

796
60
Johnson
J. L.
Whalen
D. J.
Fryer
C. L.
Li
H.
2012
ApJ

750
66
Latif
M. A.
Schleicher
D. R. G.
Schmidt
W.
Niemeyer
J.
2013a
MNRAS

433
1607
Latif
M. A.
Schleicher
D. R. G.
Schmidt
W.
Niemeyer
J. C.
2013b
MNRAS

436
2989
Latif
M. A.
Niemeyer
J. C.
Schleicher
D. R. G.
2014
MNRAS

440
2969
Le Floc'h
E.
et al
2005
ApJ

632
169
Lodato
G.
Natarajan
P.
2006
MNRAS

371
1813
P.
Haardt
F.
Dotti
M.
2014
ApJ

784
L38
Maiolino
R.
Oliva
E.
Ghinassi
F.
Pedani
M.
Mannucci
F.
Mujica
R.
Juarez
Y.
2004
A&A

420
889
Mortlock
D. J.
et al
2011
Nature

474
616
Oesch
P. A.
et al
2016
ApJ

819
129
Pacucci
F.
Ferrara
A.
2015
MNRAS

448
104
Pacucci
F.
Volonteri
M.
Ferrara
A.
2015a
MNRAS

452
1922
Pacucci
F.
Ferrara
A.
Volonteri
M.
Dubus
G.
2015b
MNRAS

454
3771
Pallottini
A.
et al
2015
MNRAS

453
2465
Planck Collaboration XIII
2015
preprint (arXiv:1502.01589)
Riechers
D. A.
et al
2013
Nature

496
329
Robitaille
T. P.
Whitney
B. A.
2010
ApJ

710
L11
Schaerer
D.
de Barros
S.
2009
A&A

502
423
Shang
C.
Bryan
G. L.
Haiman
Z.
2010
MNRAS

402
1249
Smidt
J.
Wiggins
B. K.
Johnson
J. L.
2016
preprint (arXiv:1603.00888)
Smith
A.
Bromm
V.
Loeb
A.
2016
preprint (arXiv:1602.07639)
Sobral
D.
Matthee
J.
Darvish
B.
Schaerer
D.
Mobasher
B.
Röttgering
H.
Santos
S.
Hemmati
S.
2015
ApJ

808
139
Tanvir
N. R.
et al
2009
Nature

461
1254
Visbal
E.
Haiman
Z.
Bryan
G. L.
2016
preprint (arXiv:1602.04843)
Volonteri
M.
2010
A&AR

18
279
Volonteri
M.
Rees
M. J.
2005
ApJ

633
624
Volonteri
M.
Silk
J.
Dubus
G.
2015
ApJ

804
148
Weigel
A. K.
Schawinski
K.
Treister
E.
Urry
C. M.
Koss
M.
Trakhtenbrot
B.
2015
MNRAS

448
3167
Wu
X.-B.
et al
2015
Nature

518
512
Xue
Y. Q.
et al
2011
ApJS

195
10
Yue
B.
Ferrara
A.
Salvaterra
R.
Xu
Y.
Chen
X.
2013
MNRAS

433
1556