Cosmic evolution of radio-excess active galactic nuclei in quiescent and star-forming galaxies across 0 < z < 4

Yijun Wang; Tao Wang; Daizhong Liu; Mark T. Sargent; Fangyou Gao; David M. Alexander; Wiphu Rujopakarn; Luwenjia Zhou; Emanuele Daddi; Ke Xu; Kotaro Kohno; Shuowen Jin

doi:10.1051/0004-6361/202347787

Home

All issues

Volume 685 (May 2024)

A&A, 685 (2024) A79

Full HTML

Open Access

Issue		A&A Volume 685, May 2024


Article Number		A79
Number of page(s)		25
Section		Extragalactic astronomy
DOI		https://doi.org/10.1051/0004-6361/202347787
Published online		08 May 2024

A&A, 685, A79 (2024)

Cosmic evolution of radio-excess active galactic nuclei in quiescent and star-forming galaxies across 0 < z < 4

Yijun Wang¹^,2, Tao Wang¹^,2, Daizhong Liu³, Mark T. Sargent⁴, Fangyou Gao¹^,2, David M. Alexander⁵, Wiphu Rujopakarn⁶^,7, Luwenjia Zhou¹^,2, Emanuele Daddi⁸, Ke Xu¹^,2, Kotaro Kohno⁹^,10 and Shuowen Jin¹¹^,12

¹ School of Astronomy and Space Science, Nanjing University, Nanjing 210093, PR China
² Key Laboratory of Modern Astronomy and Astrophysics (Nanjing University), Ministry of Education, Nanjing 210093, PR China
e-mail: taowang@nju.edu.cn
³ Max-Planck-Institut für extraterrestrische Physik, Giessenbachstrasse 1, 85748 Garching, Germany
⁴ International Space Science Institute (ISSI), Hallerstrasse 6, 3012 Bern, Switzerland
⁵ Centre for Extragalactic Astronomy, Department of Physics, Durham University, South Road, Durham DH1 3LE, UK
⁶ National Astronomical Research Institute of Thailand, Don Kaeo, Mae Rim, Chiang Mai 50180, Thailand
⁷ Department of Physics, Faculty of Science, Chulalongkorn University, 254 Phayathai Road, Pathumwan, Bangkok 10330, Thailand
⁸ CEA, IRFU, DAp, AIM, Université Paris-Saclay, Université Paris Diderot, Sorbonne Paris Cité, CNRS, 91191 Gif-sur-Yvette, France
⁹ Institute of Astronomy, Graduate School of Science, The University of Tokyo, 2-21-1 Osawa, Mitaka, Tokyo 181-0015, Japan
¹⁰ Research Center for the Early Universe, Graduate School of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033, Japan
¹¹ Cosmic Dawn Center (DAWN), Copenhagen, Denmark
¹² DTU Space, Technical University of Denmark, Elektrovej 327, 2800 Kgs. Lyngby, Denmark

Received: 23 August 2023
Accepted: 9 February 2024

Abstract

Context. Radio-excess active galactic nuclei (radio-AGNs) are essential to our understanding of both the physics of black hole (BH) accretion and the interaction between BHs and host galaxies. Recent deep and wide radio continuum surveys have made it possible to study radio-AGNs down to lower luminosities and up to higher redshifts than previous studies, and are providing new insights into the abundance and physical origin of radio-AGNs.

Aims. Here we focus on the cosmic evolution, physical properties, and AGN-host galaxy connections of radio-AGNs selected from a total sample of ∼400 000 galaxies at 0 < z < 4 in the GOODS-N and COSMOS fields.

Methods. Combining the deep radio continuum data with multi-band, de-blended far-infrared, and submillimeter data, we were able to identify 983 radio-AGNs out of the entire galaxy sample through radio excess relative to the far-infrared–radio relation.

Results. We studied the cosmic evolution of 1.4 GHz radio luminosity functions (RLFs) for both star-forming galaxies (SFGs) and radio-AGNs, which can be well described by a pure luminosity evolution of L_⋆ ∝ (1 + z)^{−0.34 × z + 3.57} and a pure density evolution of Φ_⋆ ∝ (1 + z)^{−0.77 × z + 2.69}, respectively. We derived the turnover luminosity, above which the number density of radio-AGNs surpasses that of SFGs. We show that this crossover luminosity increases with increasing redshifts, from 10^22.9 W Hz⁻¹ at z ∼ 0 to 10^25.2 W Hz⁻¹ at z ∼ 4. At the full redshift range of 0 < z < 4, we further derive the probability (p_radio) of SFGs and quiescent galaxies (QGs) hosting a radio-AGN, as a function of stellar mass (M_⋆), radio luminosity (L_R), and redshift (z), which yields p_radio ∝ (1+z)^3.08 M_⋆^1.06 L_R^−0.77 for SFGs, and p_radio ∝ (1+z)^2.47 M_⋆^1.41 L_R^−0.60 for QGs, respectively.

Conclusions. The quantitative relation for the probabilities of galaxies hosting a radio-AGN indicates that radio-AGNs in QGs prefer to reside in more massive galaxies with higher L_R than those in SFGs. The fraction of radio-AGN increases toward higher redshift in both SFGs and QGs, with a more rapid increase in SFGs.

Key words: galaxies: active / galaxies: evolution / galaxies: general / galaxies: luminosity function, mass function / radio continuum: galaxies

© The Authors 2024

Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

Active galactic nuclei (AGNs) are believed to play an important role in the growth of galaxies through high kinetic or radiative power (AGN feedback; Fiore et al. 2017; Matzeu et al. 2023; Fabian 2012, for a general review). AGNs emit energy across the whole spectrum, which can be identified in multiple wavelengths, including optical (Baldwin et al. 1981; Kauffmann et al. 2003; Kewley et al. 2013), mid-infrared (MIR; Lacy et al. 2004; Donley et al. 2012), X-ray (e.g., Xue et al. 2016; Luo et al. 2017), and radio bands (e.g., Del Moro et al. 2013). Approaches through optical, MIR, and X-ray surveys are sensitive to selecting AGNs with relatively high Eddington ratios, while selection through radio bands can find more AGNs with weak nuclear activities and can help build a more complete AGN sample (Delvecchio et al. 2017; Radcliffe et al. 2021). Radio-AGNs have been found to preferentially reside in massive quiescent galaxies (QGs; Condon & Dressel 1978; Brown et al. 2011; Vaddi et al. 2016; Ho 2008, for a general review), and in dense environments (Best et al. 2005; Worpel et al. 2013; Pasini et al. 2022). They may have a significant effect on the evolution of their host galaxies and environments through the mechanical or jet or radio mode feedback (Fabian 2012; Hardcastle & Croston 2020; Kondapally et al. 2023; Magliocchetti 2022, for reviews). Constructing a large and complete radio-AGN sample over a broad luminosity range and a wide redshift range can greatly increase our understanding of their abundance and occupation fraction, and improve our knowledge about the impact AGN feedback has on host galaxies and their environments.

One effective approach to selecting radio-AGNs is through selecting objects with radio emission exceeding that expected from star formation, especially for less luminous radio-AGNs. In the star formation process, IR emission is expected to be correlated with radio emission due to their mutual origins from activities of massive stars (Condon 1992; Dubner & Giacani 2015), called the IR–radio correlation (IRRC), which has been established by many observational studies (e.g., Helou et al. 1985; Condon 1992; Yun et al. 2001; Appleton et al. 2004; Ivison et al. 2010; Sargent et al. 2010a,b; Delhaize et al. 2017; Molnár et al. 2021). This IRRC can be used to identify radio-AGNs (also called radio-excess AGNs; Donley et al. 2005; Park et al. 2008; Del Moro et al. 2013; Calistro Rivera et al. 2017) because AGNs may have extra radio emission related to nuclear activities, such as radio jets, AGN-driven outflows, and the innermost accretion disk coronal activities (Panessa et al. 2019, for a review). This method can select many AGNs that are usually missed by other established techniques, such as optical, X-ray, or MIR surveys (Park et al. 2008; Del Moro et al. 2013).

Separating radio-AGNs from star-forming galaxies (SFGs) requires both sensitive radio and far-infrared (FIR) observations, which is particularly important for low-luminosity radio-AGNs. In the local universe, it is known that SFGs dominate at the low radio luminosity end, while AGNs constitute a significant fraction of radio populations at the high-luminosity end (e.g., Seymour et al. 2008; Condon et al. 2019; Franzen et al. 2021; Matthews et al. 2021). At higher redshifts, however, large uncertainties remain due to the difficulty in simultaneously detecting faint and luminous objects with both IR and radio facilities (e.g., Novak et al. 2018). The situation has been significantly improved with the accomplishment of deep and wide radio surveys, for example the 3 GHz VLA COSMOS survey (Smolčić et al. 2017b,a; Novak et al. 2017; Delvecchio et al. 2022), the Low Frequency Array Two-metre Sky Survey (LoTSS) Deep Fields (Best et al. 2023; Sabater et al. 2021; Tasse et al. 2021), and the MeerKAT International GHz Tiered Extragalactic Exploration (MIGHTEE) Survey (Jarvis et al. 2016). Moreover, even deeper VLA (Owen 2018; Alberts et al. 2020) and far-infrared surveys (Liu et al. 2018, Wang et al. in prep.) are available in the GOODS fields; these improvements greatly enlarge the radio sample with more faint objects. Furthermore, thanks to the detailed de-blended photometry for the FIR/submillimeter imaging data in the GOODS and COSMOS fields (Liu et al. 2018; Jin et al. 2018, Wang et al. in prep.), more precise estimates for the IR luminosity are available to further improve the ability to distinguish between different radio populations. Combining the power of these wide and deep surveys has great potential to constrain the abundance and physical properties of radio-AGNs down to lower luminosities and up to higher redshifts.

In addition to the abundances of radio-AGNs, the occupation fraction of radio-AGNs in different galaxy populations is also essential in order to constrain models of AGN accretion and feedback. While the most luminous radio-AGNs are primarily located in massive quiescent galaxies in the local universe (e.g., Peacock & Nicholson 1991; Magliocchetti et al. 2004; Donoso et al. 2010; Worpel et al. 2013; Kolwa et al. 2019; Dullo et al. 2023), it remains unclear whether this is the case at high redshifts (e.g., z ∼ 2; Malavasi et al. 2015), and for the less luminous galaxies (Uchiyama et al. 2022). Radio luminosities of AGNs may fundamentally depend on the efficiency of gas accretion onto the black hole (BH), which can be reflected in the BH fundamental plane (L_R − L_X − M_BH relation; e.g., Merloni et al. 2003; Bonchi et al. 2013; Xie & Yuan 2017; Bariuan et al. 2022). Many works have systematically investigated the dependence of the radio-loud AGN (RL AGN) fraction on AGN radio luminosity, stellar mass (BH mass), and galaxy types of host galaxies, especially in the local universe. RL AGN fractions are found to increase with stellar mass in the local universe (0.03 < z < 0.3; Best et al. 2005; Sabater et al. 2019) and at higher redshifts (∼2; Williams & Röttgering 2015). This stellar-mass dependence may decrease with redshift (from z ∼ 0.3 to z ∼ 2; Williams & Röttgering 2015; Zhu et al. 2023). Moreover, the RL AGN fraction may decrease with radio luminosity (Best et al. 2005; Sabater et al. 2019) and increase with redshift (Donoso et al. 2009; Zheng et al. 2022). Galaxy colors and galaxy populations have a significant effect on the RL AGN fraction (Janssen et al. 2012; Kondapally et al. 2022), which may indicate that supermassive black holes (SMBHs) are fuelled by different mechanisms in different galaxy populations (Kondapally et al. 2022). However, at higher redshift and for fainter radio-AGNs the dependence of radio-AGN fraction on diverse physical properties of galaxies is still unclear. Thanks to the deeper and wider radio surveys (e.g., Smolčić et al. 2017b; Owen 2018; Alberts et al. 2020; Best et al. 2023), studying the cosmic evolution of both weak and powerful radio-AGNs has become possible, which will greatly enrich our understanding of the physical properties of radio-AGNs and their feedback. Moreover, given that different galaxy populations (e.g., SFGs and QGs) may present different associations with AGN feedback (Fiore et al. 2017; Delvecchio et al. 2022; Matzeu et al. 2023; Fabian 2012, for a general review), systematic studies of the radio-AGN fractions in different galaxy populations are required to help us further understand the detailed effects of AGN feedback. In addition, most of the previous studies usually estimated the AGN fraction as a function of solely stellar mass, luminosity, or redshift. To systematically investigate the possible physical properties affecting the radio activities of AGNs, a unified quantitative relation describing radio-AGN fractions as a function of stellar mass, radio luminosity, and redshift is required. We may expect this unified relation to serve as an important complement to AGN feedback mode in simulations, such as IllustrisTNG (Weinberger et al. 2017; Pillepich et al. 2018) and SIMBA (Davé et al. 2019), in the future.

In this work we first use the UV-optical-MIR surveys in the GOODS-N and COSMOS/UltraVISTA fields (totaling 1.55 deg²) to construct a large galaxy sample at 0.1 < z < 4 (Sect. 2). Then we cross-match this UV-optical-MIR catalog with the deep and large radio surveys and the de-blended IR luminosity catalogs in these fields (Sect. 2). Next we calculate IR and radio luminosities for individual objects (Sect. 3), and use the IR-radio-luminosity-ratio to select radio-excess AGNs at 0.1 < z < 4 (Sect. 4). Further, we investigate the cosmic evolution and physical properties of radio-excess AGNs through the following two aspects: (1) constructing radio luminosity functions for SFGs and radio-excess AGNs, repectively, and studying their evolution with redshift (Sect. 5); (2) calculating the radio-excess AGN fraction as a function of stellar mass, radio luminosity, and redshift in different galaxy populations such as SFGs and QGs (Sect. 6). The interpretation of our results is discussed in Sect. 7. We summarize our conclusions in Sect. 8. Throughout this paper we assume a Chabrier (2003) initial mass function (IMF) and a flat cosmology with the following parameters: Ω_m = 0.3, Ω_Λ = 0.7, and H₀ = 70 km s⁻¹ Mpc⁻¹.

2. Data and samples

Our sample was selected from two fields: the Great Observatories Origins Deep Survey North (GOODS-N; Barro et al. 2019) and the UltraVISTA survey in the Cosmic Evolution Survey (COSMOS; Scoville et al. 2007; Weaver et al. 2022). The GOODS-N survey is part of the deep fields in the Cosmic Assembly Near-infrared Deep Extragalactic Legacy Survey (CANDELS; Grogin et al. 2011; Koekemoer et al. 2011) and is suitable to study the faint objects in this work. The GOODS-N field has a region size of 171 arcmin², which is not large enough to detect many luminous objects. Therefore, we also used observational data from the COSMOS/UltraVISTA field with a region size of 1.5 deg². The COSMOS/UltraVISTA field is within the full COSMOS field, which is the largest HST survey imaging a 2 deg² field (Scoville et al. 2007; Weaver et al. 2022). Here we only used the data from the COSMOS/UltraVISTA field in order to utilize the de-blended FIR-submillimeter photometry catalog from Jin et al. (2018). In addition, the GOODS-N field has deeper radio (Owen 2018) and IR (Liu et al. 2018) surveys than the COSMOS field (Smolčić et al. 2017b; Jin et al. 2018). Therefore, summing up the data from these two fields helped us construct a large sample that includes both faint and luminous objects. The available multi-wavelength data and total samples from the GOODS-N and COSMOS fields in this work are summarized in Fig. 1. Next we show the detailed matching processes and sample selections in each field (the overall flowchart is summarized in Fig. 1).

Fig. 1.

Flowchart of the matching processes and sample selections (see details in Sect. 2). ⁽¹⁾Details about broadband SED fitting can be found in Sect. 3.1 and Appendix A. ⁽²⁾Details about IR-radio correlation (IRRC) analysis can be found in Sect. 4. ⁽³⁾Details about the radio luminosity function can be found in Sect. 5. ⁽⁴⁾Details about the maximum-likelihood (ML) fitting method used for calculating the probability of hosting a radio-excess AGN can be found in Sect. 6. ^(*) L_TIR: rest-frame total IR luminosity between 8 and 1000 μm.

2.1. GOODS-N field

First we collected UV-optical-MIR data in the GOODS-N field from Barro et al. (2019; hereafter B19). This UV-optical-MIR catalog covers the wavelength range between 0.4 and 8 μm, and contains 35 445 sources over 171 arcmin². We selected 29 267 sources as our All Galaxies Sample according to two criteria: 0.1 < z ≤ 4.0; signal-to-noise ratio (S/N) of F160W band ≥5.

Next, the FIR-submillimeter-radio data in the GOODS-N field were derived from Liu et al. (2018; hereafter L18). This FIR-submillimeter-radio catalog was established using Spitzer 24 μm or VLA 20 cm detected sources as priors for FIR-submillimeter photometry that have potentially significant source confusion. This “super de-blended” photometry method provides more accurate estimates for FIR-submillimeter flux of each individual source. This catalog contains 3306 sources with photometry at MIR (Spitzer 16 and 24 μm), FIR (Herschel 100, 160, 250, 350, and 500 μm), submillimeter (SCUBA2 850 μm and AzTEC+MAMBO 1160 μm), and radio (VLA 1.4 GHz) bands. The 1.4 GHz data are primarily from Owen (2018) with rms noise (σ) in the radio image center of 2.2 μJy, supplemented by Morrison et al. (2010) with σ in the radio image center of 3.9 μJy beam⁻¹. For sources weaker than the 5σ detection limit of Owen (2018), Liu et al. (2018) performed prior-extraction photometry and Monte Carlo simulations on the radio images of both Owen (2018) and Morrison et al. (2010).

Then we cross-matched the FIR-submillimeter-radio catalog (L18) with the UV-optical-MIR catalog (B19) by a match radius of 1.5 arcsec. The match radius is defined according to the average angular resolution of 1.4 GHz VLA survey in the GOODS-N field, which is 1.6 arcsec (Owen 2018). A total of 2584 (out of 3306) objects in L18 have UV-optical-MIR counterparts in B19. The remaining 722 sources in L18 are out of the CANDELS F160W mosaic region over 171 arcmin² in the CANDELS/GOODS-N field, and are not used in this work.

Finally, among the above-mentioned 2584 sources, we selected 509 sources as our Radio Sources Sample according to two criteria: 0.1 < z ≤ 4.0; S/N of 1.4 GHz radio flux ≥5. The Radio Sources Sample has multi-wavelength data (UV-optical-IR-submillimeter-radio bands) that can be used to estimate various galaxy properties (such as 8–1000 μm IR luminosity L_TIR and stellar mass M_⋆) through broadband spectral energy distribution (SED) fitting (see details in Sect. 3.1).

2.2. COSMOS/UltraVISTA field

The UV-optical-MIR catalog in the COSMOS/UltraVISTA field was derived from the full COSMOS field (Weaver et al. 2022, hereafter W22). The full COSMOS survey covers the wavelength range between 0.2 and 8 μm, and contains 1 720 700 sources over 2 deg². The COSMOS/UltraVISTA field is within the central 1.5 deg² of the full COSMOS survey, which contains about 888 705 objects. From these 888 705 objects, 405 408 objects were selected as our All Galaxies Sample in the COSMOS/UltraVISTA field according to two criteria: 0.1 < z ≤ 4.0; S/N of Ks band ≥5.

Next, we collected FIR-submillimeter-radio data from Jin et al. (2018; hereafter J18). Similar to the analysis in the GOODS-N, J18 provided a de-blended FIR (Herschel 100, 160, 250, 350, and 500 μm) and (sub)millimeter (SCUBA2 850 μm, AzTEC 1.1 mm, and MAMBO 1.2 mm) photometric catalog for 191 624 objects in the COSMOS/UltraVISTA field. This catalog also contains observational data at Spitzer 24 μm band and radio detections at 3 GHz from the VLA-COSMOS 3 GHz project with an rms noise in the field center of 2.3 μJy beam⁻¹ (Smolčić et al. 2017b).

Then we cross-matched the FIR-submillimeter-radio catalog (J18) with the UV-optical-MIR catalog (W22) by a radius of 0.5 arcsec, and 178 494 of 191 624 objects in J18 have UV-optical-MIR counterparts in W22. The match radius is defined according to the average angular resolution of the 3 GHz VLA-COSMOS survey, which is 0.75 arcsec (Smolčić et al. 2017b).

Finally, from the above-mentioned 178 494 objects, we selected 7006 objects as our Radio Sources Sample in the COSMOS/UltraVISTA field. We followed two criteria: 0.1 < z ≤ 4.0; S/N of 3 GHz radio flux ≥5.

3. Galaxy property estimation

3.1. Rest-frame 8–1000 μm infrared luminosity

GOODS-N field. For the 2584 objects in the GOODS-N field that have multi-wavelength data from UV to FIR-submillimeter to radio bands, we estimated the total IR luminosity from 8 to 1000 μm through SED fitting with Code Investigating GALaxy Emission (CIGALE 2022.0; Burgarella et al. 2005; Noll et al. 2009; Boquien et al. 2019; Yang et al. 2020, 2022). We refer to Appendix A for more details about the SED fitting.

COSMOS/UltraVISTA field. For the sources in the COSMOS/UltraVISTA field, we used the 8–1000 μm luminosity estimated by the SED fitting in J18. Their SED components consist of a stellar component from Bruzual & Charlot (2003) with a Small Magellanic Cloud attenuation law, dust continuum emission based on Draine & Li (2007), a MIR AGN torus component from Mullaney et al. (2011), and a power-law radio component. In addition, our work used the same Chabrier (2003) IMF as was used in J18.

All the objects of the Radio Sources Sample in the GOODS-N field have available IR luminosity estimates. A total of 6985 (out of 7006) objects of the Radio Sources Sample in the COSMOS/UltraVISTA field have L_TIR estimates from J18 (∼99.7%). Moreover, we verified that different methods give similar IR luminosity measurements. For example, at 0.1 < z < 1, averaged IR luminosities of galaxies with 10^10.5 < M_⋆ < 10^11.5 M_⊙ in the GOODS-N and COSMOS/UltraVISTA fields are 3.2 × 10³⁷ and 2.6 × 10³⁷ W, repectively.

3.2. Rest-frame 1.4 GHz radio luminosity

Throughout this work the radio spectrum of each radio source is assumed to follow a simple power-law shape of S_ν ∝ ν^α, where S_ν is the flux density at frequency ν and α is the spectral index which is assumed to be −0.8 in this work. This simple power law assumed for radio spectrum is a widely used approximation due to insufficient radio data, and it is worth noting that adopting a single value of α may bring a scatter or bias to the calculation of radio luminosity (e.g., Gim et al. 2019). The rest-frame 1.4 GHz luminosity (in the unit of W Hz⁻¹) of the Radio Sources Sample can be calculated by

$\begin{matrix} L_{1.4 GHz, rest} = \frac{4 π D_{L}^{2}}{{(1 + z)}^{1 + α}} {(\frac{1.4 GHz}{ν_{obs}})}^{α} S_{ν, obs}, \end{matrix}$ $\begin{aligned} L_{1.4 \,\mathrm {GHz}, \mathrm {rest}} = \frac{4\pi D_{\rm L}^2}{(1+z)^{1+\alpha }} \left(\frac{1.4\, \mathrm {GHz}}{\nu _{\rm obs}}\right)^{\alpha } S_{\nu , \mathrm {obs}}, \end{aligned}$ (1)

where D_L is the luminosity distance (in the unit of meter), z is the redshift, ν_obs is the observed frequency (in the unit of GHz), and S_ν,obs is the observed integrated flux densities at the observed frequency ν_obs (Ceraj et al. 2018). For the GOODS-N field, we used the observed 1.4 GHz flux to calculate the rest-frame 1.4 GHz luminosity, while for the COSMOS/UltraVISTA field, we used the observed flux at 3 GHz. We show the distribution of rest-frame 1.4 GHz luminosity of our Radio Sources Sample with redshift in Fig. 2.

Fig. 2.

Distribution of rest-frame 1.4 GHz luminosity of the Radio Sources Sample (blue circles) and the Radio-excess AGN Sample (orange diamonds) with redshift for the GOODS-N (top) and COSMOS/UltraVISTA (bottom) fields. The 5σ detection limit for the COSMOS/UltraVISTA field with σ = 2.3 μJy beam⁻¹ at 3 GHz (Smolčić et al. 2017b) is shown as a gray dashed line in the top panel and as a black solid line in the bottom panel. The 5σ detection limit for the GOODS-N field with σ = 2.2 μJy at 1.4 GHz (Owen 2018) is shown as a black solid line in the top panel and as a gray dashed line in the bottom panel. All the detection limits are calculated by Eq. (1) assuming a fixed spectral index of −0.8 (see details in Sect. 3.2).

3.3. Stellar mass and UVJ magnitude

For the All Galaxies Sample in the GOODS-N field, their stellar masses and rest-frame UVJ magnitudes were derived from Barro et al. (2019), where stellar masses are estimated by SED fitting with the codes FAST (Kriek et al. 2009) and Synthesizer (Pérez-González et al. 2005, 2008), and rest-frame UVJ luminosities are estimated by EAZY (Brammer et al. 2008). For the All Galaxies Sample in the COSMOS/UltraVISTA field, their stellar masses and rest-frame UVJ magnitudes were derived from Weaver et al. (2022), where stellar masses are estimated by LePhare (Arnouts et al. 2002; Ilbert et al. 2006) and rest-frame UVJ luminosities are estimated by EAZY (Brammer et al. 2008). Our work and the works mentioned here used the same Chabrier (2003) IMF. In addition, to examine whether these different methods will bring measurement bias for stellar mass, we also compared stellar masses given by the above methods ( $M_{⋆}^{FAST / LePhare}$ $M_{\star}^\mathrm{\mathtt{FAST}/\mathtt{LePhare}}$ ) with those obtained by the broadband SED fitting ( $M_{⋆}^{SED}$ $M_{\star}^{\mathrm{SED}}$ ; see details in Appendix A for the GOODS-N field and in Jin et al. 2018 for the COSMOS field). The logarithm of the ratio of $M_{⋆}^{FAST / LePhare}$ $M_{\star}^\mathrm{\mathtt{FAST}/\mathtt{LePhare}}$ to $M_{⋆}^{SED}$ $M_{\star}^{\mathrm{SED}}$ ( $log [M_{⋆}^{FAST / LePhare} / M_{⋆}^{SED}]$ $\log [M_{\star}^\mathrm{\mathtt{FAST}/\mathtt{LePhare}}/M_{\star}^{\mathrm{SED}}]$ ) is around 0 ± 0.3 dex. This result indicates that different methods give consistent measurements for the stellar mass. Given that only the Radio Sources Sample and Radio-excess AGN Sample have available broadband SED fitting results, while the All Galaxies Sample do not have these measures, we used the stellar mass measured by FAST or LePhare in this work.

3.4. Galaxy type: Star-forming and quiescent galaxies

Color-magnitude and color-color criteria are effective methods of selecting SFGs and QGs, and do not require accurate measurements for the star formation rate and stellar mass (e.g., Williams et al. 2009). For this work we used the UVJ selection criteria in Schreiber et al. (2015) to separate all our samples into SFGs and QGs at all redshift and stellar masses:

$\begin{matrix} quiescent = {\begin{matrix} U - V > 1.3 and \\ V - J < 1.6 and \\ U - V > 0.88 \times (V - J) + 0.49 \end{matrix}, \end{matrix}$ $\begin{aligned} \mathrm {quiescent} = {\left\{ \begin{array}{ll} U-V > 1.3\ \mathrm {and} \\ V-J < 1.6\ \mathrm {and} \\ U-V > 0.88\times (V-J)+0.49 \\ \end{array}\right.} ,\end{aligned}$ (2)

$\begin{matrix} star-forming = {\begin{matrix} U - V \leq 1.3 or \\ V - J \geq 1.6 or \\ U - V \leq 0.88 \times (V - J) + 0.49 . \end{matrix} \end{matrix}$ $\begin{aligned} {\text{ star-forming}}\ = {\left\{ \begin{array}{ll} U-V \le 1.3\ \mathrm {or} \\ V-J \ge 1.6\ \mathrm {or} \\ U-V \le 0.88\times (V-J)+0.49. \\ \end{array}\right.} \end{aligned}$ (3)

4. Radio-excess active galactic galaxies

The excess radio emission from AGNs makes these types of sources deviate from the IRRC of SFGs. Therefore, AGN activities can be identified by these excess radio emission. Following Helou et al. (1985), we define the logarithmic ratio between total IR emission (8–1000 μm; TIR) and radio emission as

$\begin{matrix} q_{TIR} = log (\frac{L_{TIR}}{L_{1.4 GHz, rest} \times 3.75 \times 10^{12} Hz}), \end{matrix}$ $\begin{aligned} q_{\rm {TIR}} = \log \left( \frac{L_{\rm {TIR}}}{L_{1.4 \mathrm {GHz}, \mathrm {rest}} \times 3.75 \times 10^{12}\ \mathrm {Hz}} \right), \end{aligned}$ (4)

where L_TIR is the total IR luminosity and 3.75 × 10¹² Hz is the frequency at the center of FIR band (λ = 80μm). In this work we consider a broad redshift range from 0.1 to 4, which we divided into six redshift bins: 0.1 < z ≤ 0.5, 0.5 < z ≤ 1.0, 1.0 < z ≤ 1.5, 1.5 < z ≤ 2.0, 2.0 < z ≤ 3.0, and 3.0 < z ≤ 4.0. In this section we focus on the Radio Sources Sample and then use the q_TIR distribution to select AGNs with excess radio emission, hereafter referred to as radio-excess AGNs.

4.1. q_TIR distribution and radio-excess AGN selection

In this analysis we merge the radio sources in the GOODS-N and COSMOS/UltraVISTA fields together for the following reasons. First, source numbers in each redshift bin for the GOODS-N field are small, which may result in a large uncertainty in determining the AGN selection threshold. Second, in each redshift bin, the GOODS-N and COSMOS fields have similar q_TIR distributions and ${\bar{q}}_{TIR}$ $\bar{q}_\mathrm{{TIR}}$ values (see Fig. B.1). The q_TIR distribution consists of a Gaussian component peaked at ${\bar{q}}_{TIR}$ $\bar{q}_\mathrm{{TIR}}$ and an extended tail toward lower q_TIR, which is consistent with previous works (e.g., Del Moro et al. 2013). The Gaussian component is attributed to the IRRC of SFGs, while the extended tail is thought to be associated with the extra radio emission from AGN. Therefore, to ensure the same selection criteria and reduce the uncertainty of the AGN selection threshold, we combine the data from the two fields together. The consistent radio luminosity functions for these two fields also indicate that our selection method is plausible (see details in Sect. 5). Next we analyze the q_TIR distribution in SFGs and QGs.

The entire q_TIR distribution can be described by multiple Gaussian models (see the yellow dotted curves in Fig. 3). We used f-test probability at ≥95% confidence level to decide how many models were required. The final best-fit models are shown in Fig. 3 with black solid curves. As we mentioned above, the highest Gaussian peak and the extended low-q_TIR tail are thought to be contributed by SFGs and AGNs, respectively. Therefore, here we define the cross point between the highest Gaussian and the second-highest Gaussian components as the threshold to separate the SFGs and radio-excess AGNs (see the q_AGN values and red dash-dotted lines in Fig. 3). It means that a radio source will be selected as a radio-excess AGN if its q_TIR value is lower than this cross point (q_AGN). It is worth noting that the q_TIR distributions of QGs also show a highest peak similar to SFGs. However, QGs are not expected to have this peak due to little or no star formation. One possible reason may be that UVJ color selection might classify a part of the massive SFGs with relatively low star formation rates (SFRs) as QGs (Popesso et al. 2023). We also found that nearly half of UVJ-selected QGs (without AGNs) in the Radio Sources Sample are located very close to the dividing line in the UVJ diagram, which indicates that UVJ selection might not accurately classify these sources. In addition, it is worth noting that the QGs in the All Galaxies Sample do not show such a distribution in the UVJ diagram. Therefore, another possible reason may be that radio detection favors sources with higher star formation, even among QGs. Radio-excess AGNs hosted by QGs do not show a prominent population close to the dividing line in the UVJ diagram, which may be due to the radio-excess AGNs being selected by a threshold (q_AGN) that is significantly lower than the mean q_TIR of SFGs. This may explain why QGs show a highest q_TIR peak in Fig. 3. In addition, the difference between the second-highest and the third-highest Gaussian component in the q_TIR distribution might be associated with different types of radio-AGNs, which is beyond the scope of this work and will be studied in our future works.

Fig. 3.

Distribution of IR-radio-ratio (q_TIR) in SFGs (left column) and QGs (right column) at different redshift bins. The q_TIR distributions were obtained by merging the Radio Sources Sample from the GOODS-N and COSMOS/UltraVISTA fields. Each yellow dotted curve represents a single-Gaussian model, while the black solid curve represents the best-fit model to the entire q_TIR distribution. The vertical black dashed line (corresponding to the ${\bar{q}}_{TIR}$ $\bar{q}_{\mathrm{TIR}}$ value in the figure) represents the peak position of the highest Gaussian component. The vertical red dash-dotted line (corresponding to the q_AGN value in the figure) represents the cross point between the highest Gaussian component and the second-highest Gaussian component, which is used as the threshold to separate SFGs and radio-excess AGNs in this work. For QGs at 2.0 < z ≤ 4.0, it is difficult to find a best-fit model due to the low source number, so here the threshold at 0.1 < z ≤ 0.5 was used, which does not have a significant effect on the results.

Finally, we selected 983 radio-excess AGNs at 0.1 < z ≤ 4 as our Radio-excess AGN Sample. This sample contains 102 sources in the GOODS-N field and 881 sources in the COSMOS/UltraVISTA field (see Fig. 2).

4.2. q_TIR evolution

As Fig. 4 shows, ${\bar{q}}_{TIR}$ $\bar{q}_\mathrm{{TIR}}$ presents a weak evolution with redshift, which is in the form of

Fig. 4.

Distribution of q_TIR with redshift. The gray points show the Radio Sources Sample in the GOODS-N and COSMOS/UltraVISTA fields (see details in Sect. 2). The blue points represent the ${\bar{q}}_{TIR}$ $\bar{q}_{\mathrm{TIR}}$ values of SFGs, which are shown in Fig. 3 and in Sect. 4.1. The blue dashed line represents the best-fit model to the data of SFGs, while the shaded blue region is the associated 1σ uncertainty of the best-fit model. The green dash-dotted line represents the best-fit model from Magnelli et al. (2015) using data in the GOODS-N, GOODS-S, ECDFS, and COSMOS fields. Both the red dotted line (Delhaize et al. 2017) and the yellow dash-dotted line (Novak et al. 2017) represent the best-fit model obtained with the VLA-COSMOS 3 GHz project but using different source selection criteria. The dark red dash-dotted line represents the best-fit model from Enia et al. (2022) based on data from the GOODS-N field. The black star represents the q_TIR value for SFGs in the local universe from Bell (2003).

$\begin{matrix} {\bar{q}}_{TIR} = (2.62 \pm 0.08) \times {(1 + z)}^{- 0.08 \pm 0.03} . \end{matrix}$ $\begin{aligned} \bar{q}_{\rm {TIR}} = (2.62\pm 0.08)\times (1+z)^{-0.08\pm 0.03}. \end{aligned}$ (5)

Our result is almost consistent with previous works (see Fig. 4; Bell 2003; Magnelli et al. 2015; Delhaize et al. 2017; Novak et al. 2017; Enia et al. 2022). Bell (2003) obtained q_TIR value for SFGs in the local universe. Magnelli et al. (2015) combined data from multiple fields of the GOODS-N, GOODS-S, ECDFS, and COSMOS up to z ∼ 2 with rms noise greater than 3.9 μJy at 1.4 GHz. Delhaize et al. (2017) and Novak et al. (2017) both used the VLA-COSMOS 3 GHz project, but used different source selection criteria. Enia et al. (2022) focused on SFGs in the GOODS-N field. Compared to their works, we use simultaneously deeper and larger radio surveys, and more precise IR luminosity based on de-blended IR photometry. The slight difference between our results and previous works might be due to the different selection method for sources. In addition, the averaged standard deviation (1σ error of q_TIR in Fig. 4) of q_TIR distribution for SFGs is around 0.18 which is consistent with the previous works (e.g., Del Moro et al. 2013).

More recently, Delvecchio et al. (2021) found that the IRRC of SFGs strongly depends on stellar mass (M_⋆), while An et al. (2021) found a weak dependence. To examine whether considering M_⋆-dependence or not affects our results, we performed some simple tests. For the Radio Sources Sample in the COSMOS/UltraVISTA field, at 0.1 < z < 1, ${\bar{q}}_{TIR}$ $\bar{q}_{\mathrm{TIR}}$ is around 2.55 and 2.53 for low-mass (10¹⁰ < M_⋆ < 10^10.5 M_⊙) and high-mass (10^10.5 < M_⋆ < 10¹¹ M_⊙) SFGs, respectively. For the Radio Sources Sample in the GOODS-N field, at 0.1 < z < 1, ${\bar{q}}_{TIR}$ $\bar{q}_{\mathrm{TIR}}$ is around 2.56 and 2.51 for low-mass (10¹⁰ < M_⋆ < 10^10.5 M_⊙) and high-mass (10^10.5 < M_⋆ < 10¹¹ M_⊙) SFGs, respectively. These results indicate that for our samples, the ${\bar{q}}_{TIR}$ $\bar{q}_{\mathrm{TIR}}$ of SFGs shows no significant dependence on stellar mass within the uncertainties of q_TIR. Therefore, disregarding the M_*-dependent IRRC may have no significant effect on our results.

5. Radio luminosity functions and their cosmic evolution

In this section we construct the radio luminosity function (RLF) for SFGs and radio-excess AGNs, which are separated in Sect. 4. First, we estimate the combined 1.4 GHz RLF from the GOODS-N and COSMOS/UltraVISTA fields (Sect. 5.1). We then show the fitting procedures and results for the evolution of SFG RLFs (Sect. 5.2) and AGN RLFs (Sect. 5.3). In this work we do not consider the star-formation contamination for the 1.4 GHz radio luminosity of AGN. We verified that for most of our Radio-excess AGN Sample, radio luminosity from the star formation is not significant compared to luminosity from AGN.

5.1. Constructing radio luminosity functions out to z ∼ 4

Using the 1/V_max method (Schmidt 1968), the RLF in each 1.4 GHz luminosity bin and each redshift bin is calculated by

$\begin{matrix} Φ (L, z) = \frac{1}{Δ log L} \sum_{i} \frac{1}{\frac{Ω}{4 π} \times V_{max, i} \times w_{i}}, \end{matrix}$ $\begin{aligned} \Phi (L,z) = \frac{1}{\Delta \log L} \sum _i \frac{1}{\frac{\Omega }{4 \pi } \times V_{\max ,i} \times w_i}, \end{aligned}$ (6)

where ΔlogL is the size of the 1.4 GHz luminosity bin, Ω is the observed area (171 arcmin² for the GOODS-N field and 1.5 deg² for the COSMOS/UltraVISTA field), V_max, i is the co-moving volume of the ith source that is defined as V_max, i = V_{z_max, i} − V_{z_min}, and w_i is the completeness and bias correction factor of the ith source. Further, V_{z_max, i} is the co-moving volume at the maximum redshift where the ith source can be observed given the 1.4 GHz flux detection limit (the maximum value of z_max is equal to the upper limit of each redshift bin) and V_{z_min} is the co-moving volume at the lower boundary of each redshift bin. The parameter w_i is the flux density completeness of our catalogs, which takes into account the effects of the sensitivity limit. The completeness and bias corrections of the COSMOS and GOODS-N fields are derived from Smolčić et al. (2017b) and Enia et al. (2022), respectively, which are estimated by Monte Carlo simulations where mock sources are inserted in and retrieved from the mosaic. The uncertainty of RLF in each luminosity bin and each redshift bin can be defined as

$\begin{matrix} σ_{Φ} (L, z) = \frac{1}{Δ log L} \sqrt{\sum_{i} {(\frac{1}{\frac{Ω}{4 π} \times V_{max, i} \times w_{i}})}^{2}} . \end{matrix}$ $\begin{aligned} \sigma _{\Phi } (L,z) = \frac{1}{\Delta \log L} \sqrt{\sum _i \left(\frac{1}{\frac{\Omega }{4 \pi } \times V_{\max ,i} \times w_i}\right)^2}. \end{aligned}$ (7)

We used Eqs. (6) and (7) to calculate RLFs and their uncertainties for the sources in the GOODS-N and COSMOS/UltraVISTA, respectively (see Fig. C.1). For both SFGs and radio-excess AGNs, the RLFs in these fields present generally consistent results. At the faintest end of the lowest redshift bin (0.1 < z ≤ 0.5), SFG RLF of the COSMOS/UltraVISTA field presents a slight decline compared to that of the GOODS-N field (see Fig. C.1). This may be due to the higher spatial resolution of the radio survey in the COSMOS/UltraVISTA field (see details in Sect. 2) resulting in more flux losses for the extended sources, especially for the faint populations at low redshift. In order to well constrain the RLFs, in the following sections we use the averaged SFG RLFs in the GOODS-N and COSMOS/UltraVISTA fields, and use the averaged AGN RLFs in the two fields to make further analysis (see Fig. 5). In addition, we also correct the classification purity for the radio-excess AGN in these three fields (see details in Appendix B).

Fig. 5.

1.4 GHz radio luminosity function of SFGs (blue circles) and radio-excess AGNs (red squares) in different redshift bins. For SFGs, the blue solid line represents the best-fit PLE model, which is the best model to describe the data in this work. For radio-excess AGNs, the red solid line represents the best-fit PDE model, which is the best model to describe the data. The magenta vertical dash-dotted line represents the turnover luminosity between the radio luminosity function of SFGs and radio-excess AGNs (see details in Sect. 5.4). The light blue circles show the results for SFGs from the VLA-COSMOS 3 GHz project (Novak et al. 2017), while the gold squares show the results for radio-AGNs from the VLA-COSMOS 3 GHz project (Smolčić et al. 2017a). The redshift bins in this work have different definitions than in Novak et al. (2017) and Smolčić et al. (2017a).

5.2. Fitting the evolution of a star-forming galaxy’s radio luminosity function

We assume a modified-Schechter function (Saunders et al. 1990; Smolčić et al. 2009; Gruppioni et al. 2013; Novak et al. 2017) to describe the SFG RLF. The local RLF of SFGs (Novak et al. 2017) is defined as

$\begin{matrix} Φ_{0}^{SFG} (L) = Φ_{⋆, 0} {(\frac{L}{L_{⋆, 0}})}^{1 - α} exp [- \frac{1}{2 σ^{2}} {log}^{2} (1 + \frac{L}{L_{⋆, 0}})], \end{matrix}$ $\begin{aligned} \Phi _{0}^\mathrm{SFG} (L) = \Phi _{\star ,0}\left(\frac{L}{L_{\star ,0}}\right)^{1-\alpha }\exp \left[-\frac{1}{2\sigma ^2}\log ^2\left(1+\frac{L}{L_{\star ,0}}\right)\right], \end{aligned}$ (8)

where L_⋆,0 = 1.85 × 10²¹ W Hz⁻¹ is the turnover position of the local SFG RLF, Φ_⋆,0 = 3.55 × 10⁻³ Mpc⁻³ dex⁻¹ is the local turnover normalization, and α = 1.22 and σ = 0.63 are used to fit the distribution in the faint and bright ends, respectively.

The RLFs of SFGs have been shown to evolve with redshift (e.g., Novak et al. 2017; Ocran et al. 2020; Cochrane et al. 2023). Given that it is difficult to simultaneously constrain the evolutions of all the parameters, we fix the values of α and σ at all cosmic times to those of the local RLF (Novak et al. 2017) in the fit. This means that we assume an unchanged RLF shape at all cosmic times, and only allow the turnover position (L_⋆) and turnover normalization (Φ_⋆) to change with redshift. In reality, α and σ might change with redshift.

The L_⋆ evolution and Φ_⋆ evolution can be described by simple power laws:

$\begin{matrix} L_{⋆} = {(1 + z)}^{α_{L}} \times L_{⋆, 0} \end{matrix}$ $\begin{aligned} L_{\star } = (1+z)^{\alpha _{\rm L}} \times L_{\star ,0}\ \end{aligned}$ (9)

and

$\begin{matrix} Φ_{⋆} = {(1 + z)}^{α_{D}} \times Φ_{⋆, 0} . \end{matrix}$ $\begin{aligned} \Phi _{\star } = (1+z)^{\alpha _{\rm D}} \times \Phi _{\star ,0}. \end{aligned}$ (10)

Therefore, the redshift-evolved SFG RLF (luminosity and density evolution, LDE) is defined as

$\begin{matrix} Φ_{LDE}^{SFG} (L, z, α_{L}, α_{D}) = {(1 + z)}^{α_{D}} \times Φ_{0}^{SFG} (\frac{L}{{(1 + z)}^{α_{L}}}), \end{matrix}$ $\begin{aligned} \Phi ^\mathrm{SFG}_{\rm LDE}(L, z, \alpha _{\rm L}, \alpha _{\rm D}) = (1+z)^{\alpha _{\rm D}}\times \Phi _{0}^\mathrm{SFG}\left(\frac{L}{(1+z)^{\alpha _{\rm L}}}\right), \\ \end{aligned}$ (11)

where α_D is the parameter of pure density evolution (PDE; vertical shift of RLF), α_L is the parameter of pure luminosity evolution (PLE; horizontal shift of RLF), and $Φ_{0}^{SFG} (L)$ $\Phi_{0}^{\mathrm{SFG}}(L)$ is the local SFG RLF given by Eq. (8). If we assume that there is no evolution with redshift for Φ_⋆ and only L_⋆ evolves with redshift (α_D = 0), we can define a SFG PLE model as

$\begin{matrix} Φ_{PLE}^{SFG} (L, z, α_{L}) = Φ_{0}^{SFG} (\frac{L}{{(1 + z)}^{α_{L}}}) . \end{matrix}$ $\begin{aligned} \Phi ^\mathrm{SFG}_{\rm PLE}(L, z, \alpha _{\rm L}) = \Phi _{0}^\mathrm{SFG}\left(\frac{L}{(1+z)^{\alpha _{\rm L}}}\right). \end{aligned}$ (12)

Similarly, the SFG PDE model only allows Φ_⋆ to evolve with redshift (α_L = 0), which is defined as

$\begin{matrix} Φ_{PDE}^{SFG} (L, z, α_{D}) = {(1 + z)}^{α_{D}} \times Φ_{0}^{SFG} (L) . \end{matrix}$ $\begin{aligned} \Phi ^\mathrm{SFG}_{\rm PDE}(L, z, \alpha _{\rm D}) = (1+z)^{\alpha _{\rm D}}\times \Phi _{0}^\mathrm{SFG}(L). \end{aligned}$ (13)

Then we used the Markov chain Monte Carlo (MCMC) algorithm in the Python package EMCEE (Foreman-Mackey et al. 2013) to fit the data (see Table D.1) with these three models. The averaged reduced χ² over all the redshift bins of the best-fit LDE, PLE, and PDE models are 3.9, 4.8, and 526.4, respectively. First, we ignored the PDE model because it has the worst fit to the observational data. The best-fit LDE and PLE models show significant differences at higher redshift (z ≳ 1.5) where our data points only sample the bright end of RLF. This results in a degeneracy in estimation of α_D and α_L, preventing a precise calculation to the turnover position and turnover normalization. Therefore, for simplicity, we only consider the PLE model (see Eq. (12)) here.

The best-fit α_L of the PLE model in each redshift bin is summarized in Table 1. The parameter α_L shows an evolution with redshift, which can be described by α_L = A_L × z + B_L. We obtain the best-fit results with A_L = −0.34 ± 0.01 and B_L = 3.57 ± 0.01 (see Panel A of Fig. 6), which are also summarized in Table 1. The turnover position L_⋆ is estimated by Eq. (9), which is also listed in Table 1. Our results are generally consistent with those in Novak et al. (2017) with the VLA-COSMOS 3 GHz project (see Figs. 5 and 6), although our results show a slightly larger B_L (see Panel A of Fig. 6) and a slightly larger L_⋆ (see Panel B of Fig. 6), which may be due to the different selection criteria adopted for radio SFGs. Novak et al. (2017) separated the entire radio sources sample into SFGs and radio-excess AGNs subsets, while we first separate our entire radio sources sample into SFGs and QGs according to the UVJ selection criteria, and further distinguish the radio-excess AGNs from them.

Table 1.

Best-fit evolution parameters of the radio luminosity function for SFGs and radio-excess AGNs.

Fig. 6.

Evolution of the best-fit parameters for the 1.4 GHz radio radio luminosity function with redshift. Shown are the evolution parameter (α_L) of the pure luminosity evolution (PLE) for SFGs (Panel A); turnover position (L_⋆) of the PLE for SFGs, calculated by Eq. (9) (Panel B); evolution parameter (α_D) of the pure density evolution (PDE) for radio-excess AGN (Panel C); turnover normalization (Φ_⋆) of the PDE for radio-excess AGN, calculated by Eq. (10) (Panel D); and turnover luminosity (L_D) above which the number density of radio-AGNs surpasses that of SFGs (L_D) (Panel E). The gray dash-dotted lines in Panels A and B represent the results for SFGs from Novak et al. (2017), while the gray regions represent the corresponding 1σ uncertainty. The gray dash-dotted lines in Panels C and D represent the results for radio-excess AGNs from Smolčić et al. (2017a), while the gray regions represent the corresponding 1σ uncertainty. The stars in Panels B and D represent the local values obtained by Novak et al. (2017) and Mauch & Sadler (2007), respectively. The magenta star in Panel E represents the local L_D obtained by the cross point between the local SFGs RLF from Novak et al. (2017) and the local AGN RLF from Mauch & Sadler (2007). The best-fit values for α_L, α_D, L_⋆, Φ_⋆, and L_D are summarized in Table 1.

5.3. Fitting the evolution of an AGN’s radio luminosity function

The AGN RLF can be assumed as a double power-law shape (Mauch & Sadler 2007). The local RLF of AGN (Mauch & Sadler 2007) is defined as

$\begin{matrix} Φ_{0}^{AGN} (L) = \frac{Φ_{⋆, 0}}{{(L_{⋆, 0} / L)}^{α} + {(L_{⋆, 0} / L)}^{β}}, \end{matrix}$ $\begin{aligned} \Phi _{0}^\mathrm{AGN} (L) = \frac{\Phi _{\star ,0}}{(L_{\star ,0}/L)^{\alpha }+(L_{\star ,0}/L)^{\beta }}, \end{aligned}$ (14)

where $Φ_{⋆, 0} = \frac{1}{0.5} 10^{- 5.5} {Mpc}^{- 3} {dex}^{- 1}$ $\Phi_{\star,0}=\frac{1}{0.5}10^{-5.5}\ \rm{Mpc}^{-3}\ \rm{dex}^{-1}$ is the turnover normalization of the local AGN RLF, L_⋆,0 = 10^24.59 W Hz⁻¹ is the local turnover position, and α = −1.27 and β = −0.49 are the indices at the bright and faint end, respectively.

Similar to the SFG RLF, we assume that AGN RLF has an unchanged shape at all cosmic times. Thus, similar to Eqs. (11)–(13), the LDE, PLE, and PDE models to describe AGN RLFs can be respectively defined as

$\begin{matrix} Φ_{LDE}^{AGN} (L, z, α_{L}, α_{D}) = {(1 + z)}^{α_{D}} \times Φ_{0}^{AGN} (\frac{L}{{(1 + z)}^{α_{L}}}), \end{matrix}$ $\begin{aligned} \Phi ^\mathrm{AGN}_{\rm {LDE}}(L, z, \alpha _{\rm L}, \alpha _{\rm D}) = (1+z)^{\alpha _{\rm D}}\times \Phi _{0}^\mathrm{AGN}\left(\frac{L}{(1+z)^{\alpha _{\rm L}}}\right), \end{aligned}$ (15)

$\begin{matrix} Φ_{PLE}^{AGN} (L, z, α_{L}) = Φ_{0}^{AGN} (\frac{L}{{(1 + z)}^{α_{L}}}), \end{matrix}$ $\begin{aligned} \Phi ^\mathrm{AGN}_{\rm {PLE}}(L, z, \alpha _{\rm L}) = \Phi _{0}^\mathrm{AGN}\left(\frac{L}{(1+z)^{\alpha _{\rm L}}}\right), \end{aligned}$ (16)

$\begin{matrix} Φ_{PDE}^{AGN} (L, z, α_{D}) = {(1 + z)}^{α_{D}} \times Φ_{0}^{AGN} (L) . \end{matrix}$ $\begin{aligned} \Phi ^\mathrm{AGN}_{\rm {PDE}}(L, z, \alpha _{\rm D}) = (1+z)^{\alpha _{\rm D}}\times \Phi _{0}^\mathrm{AGN}(L). \end{aligned}$ (17)

We also use the EEMCEE package to fit the observational data for radio-excess AGNs (see Table D.2) with these three models. The averaged reduced χ² over all the redshift bins of the best-fit LDE, PDE, and PLE models are 2.3, 2.2, and 2.5, respectively. Here we only consider the PDE model (see Eq. (16)), which has the best fit to the data.

The best-fit α_D of the PDE model in each redshift bin is summarized in Table 1, which also shows a redshift evolution in a form of α_D = A_D × z + B_D. The best-fit A_D and B_D are −0.77 ± 0.06 and 2.69 ± 0.09, respectively (see Table 1). As we assumed an unchanged AGN RLF shape and a pure density evolution, the evolution of Φ_⋆ (calculated by Eq. (10)) with redshift can simply represent the number density distribution of AGN, which peaked at z ∼ 1.5 (see Panel D of Fig. 6). Our result is generally consistent with that in Smolčić et al. (2017a) based on the VLA-COSMOS 3 GHz project (see Figs. 5 and 6), and is consistent with the AGN accretion rate density history obtained with X-ray surveys (Aird et al. 2010, and references therein).

We also study the RLF of the Radio-excess AGN Sample hosted by different galaxy populations (SFGs and QGs; see Fig. E.1). At z < 1 the radio-excess AGN population is mainly dominated by the sources hosted by QGs, while radio-excess AGNs hosted by SFGs only have a comparable space densities at the faint end. At z > 1 the space densities of radio-excess AGNs hosted by QGs decrease toward higher redshift, while radio-excess AGNs hosted by SFGs show the opposite evolution trend. In addition, radio-excess AGNs hosted by SFGs completely dominate those hosted by QGs at z > 1.5. These trends are consistent with those for low-excitation radio galaxies (LERGs) in Kondapally et al. (2022). LERGs in Kondapally et al. (2022) are defined as objects showing powerful radio emission from AGN but not identified in other bands (e.g., IR and X-ray), which represent 89% of their entire radio-excess AGNs (11% for high-excitation radio galaxies). For our radio-excess AGNs, nearly 80% of objects do not have X-ray detections or have intrinsic 2–10 keV X-ray luminosity lower than 10⁴² erg s⁻¹. Our radio-excess AGNs may be otherwise identified by IR or other methods, so LERG fraction in our Radio-excess AGN Sample is expected to be lower than 80%, while the fraction of high-excitation radio galaxies (HERGs) will be higher than 20%. In this work our Radio-excess AGN Sample includes both LERGs and HERGs subsets, while Kondapally et al. (2022) mainly focused on LERGs. Even so, both their work and our work found that radio-AGNs in SFGs and QGs have significantly different cosmic evolutions. Given that SFGs and QGs may have different fuelling mechanisms toward central SMBHs (e.g., Kauffmann & Heckman 2009; Kondapally et al. 2022; Ni et al. 2023), these results indicate that the radio activities of central engines may depend on the fuelling mechanisms of their host galaxies. In Sect. 6, we further investigate this topic by studying the probability of radio-excess AGNs hosted by different galaxy populations.

5.4. Crossover luminosity between SFG RLF and AGN RLF

The comparison between the SFG RLF and AGN RLF shows that AGNs dominate the radio populations beyond a certain luminosity (hereafter turnover luminosity L_D; see the magenta dash-dotted lines in Fig. 5). In the local universe this crossover luminosity is around 10²³ W Hz⁻¹, which is usually used to select RL AGNs (e.g., Best et al. 2005; Kukreti et al. 2023). Toward higher redshifts, given that both SFG RLFs and AGN RLFs show evolutions with redshift, this crossover luminosity is not expected to be constant at all cosmic time. Our results indicate that the redshift evolution of this crossover luminosity can be described by

$\begin{matrix} L_{D} = L_{D, 0} \times {(1 + z)}^{α_{LD}}, \end{matrix}$ $\begin{aligned} L_{\rm D} = L_{\rm D,0} \times (1+z)^{\alpha _{\rm LD}}, \end{aligned}$ (18)

where L_D, 0 is the crossover luminosity in the local universe and α_LD is the evolution index of the crossover luminosity. Here, L_D, 0 = 10^22.9 W Hz⁻¹, is obtained by the cross point between the local SFG RLF (Novak et al. 2017) and local AGN RLF (Mauch & Sadler 2007). Further, α_LD also shows a weak evolution with redshift:

$\begin{matrix} α_{LD} = (0.06 \pm 0.01) \times z + (3.05 \pm 0.02) . \end{matrix}$ $\begin{aligned} \alpha _{\rm LD} = (0.06\pm 0.01)\times z + (3.05\pm 0.02). \end{aligned}$ (19)

These results provide us with another way to select powerful radio-AGNs at different redshifts solely through radio survey. This means that if IR data are not available, a radio source at redshift z can be selected as a radio-AGN when its 1.4 GHz radio luminosity is higher than L_D = 10^22.9 × (1 + z)^{0.06 × z + 3.05} W Hz⁻¹. This crossover luminosity describes a luminosity threshold where radio-AGNs begin to dominate the entire radio population. This result is consistent with those from previous works (e.g., McAlpine et al. 2013), and also consistent with the radio luminosity threshold proposed for radio-AGN selection in the past (e.g., Magliocchetti et al. 2017).

6. The probability of hosting radio-excess AGNs in SFGs and QGs

In this section we study the probability of different galaxy populations hosting a radio-excess AGN. Here we use the combined All Galaxies Sample (see Sect. 2 and Fig. 1) and combined Radio-excess AGN Sample (see Sect. 4 and Fig. 1) in the GOODS-N and COSMOS/UltraVISTA fields. The All Galaxies Sample and the host galaxies of the Radio-excess AGN Sample had been divided into SFGs and QGs according to the UVJ selection method (see Sect. 3.4). In Sects. 6.1 and 6.2, we use the observational data to calculate the probability of a SFG or a QG with the stellar mass M_⋆ and at redshift z hosting a radio-excess AGN with the 1.4 GHz luminosity L_R in each redshift bin. In Sect. 6.3 we use the method of Aird et al. (2012) to quantitatively calculate the probability of a SFG or a QG hosting a radio-excess AGN as a function of stellar mass, radio luminosity, and redshift. This method is based on the maximum-likelihood approach, which does not require data binning for stellar mass, luminosity, and redshift (Aird et al. 2012). Here we do not consider the star-formation contamination for the AGN radio luminosity, while for most of our Radio-excess AGN Sample the contributions from star formation to the radio luminosity have been verified as not significant compared to that from AGN.

6.1. Calculating the probability for observational data

In order to conveniently compare the data and the best-fit model, we divided our All Galaxies Sample and Radio-excess AGN Sample into six redshift bins (the same as the redshift bins in Sect. 4). Within each redshift bin we subdivided our All Galaxies Sample and Radio-excess AGN Sample into different stellar mass bins. The conditional probability density function p(L_R | M_⋆, z) describes the probability of a galaxy with stellar mass M_⋆ and at redshift z hosting a radio-excess AGN with 1.4 GHz luminosity L_R. Following Aird et al. (2012), here p(L_R | M_⋆, z) is defined as the probability density per logarithmic luminosity interval (units are dex⁻¹). Thus, p(L_R | M_⋆, z) can be converted to AGN fraction (f_AGN) according to

$\begin{matrix} f_{AGN} (M_{⋆}, z) = \int_{log L_{\lim}}^{\infty} p (L_{R} | M_{⋆}, z) d log L_{R} . \end{matrix}$ $\begin{aligned} f_{\rm AGN}(M_{\star },\ z) = \int _{{\log L_{\rm lim}}}^{\infty }p(L_{\rm R}\ |\ M_{\star },\ z)\ d \log L_{\rm R}. \end{aligned}$ (20)

Here f_AGN is defined as the fraction of galaxies with stellar mass M_⋆ and at redshift z that host a radio-excess AGN, and L_lim is the lower limit of the radio luminosity. For the combined observational data from the GOODS-N and COSMOS/UltraVISTA fields, in each redshift bin, p(L_R | M_⋆, z) in the mth stellar mass bin and the nth luminosity bin is defined by

$\begin{matrix} p (L_{n} | M_{m}, z) = \frac{N_{AGN, m n}}{N_{gal, m} Δ log L_{n}} = \frac{f_{AGN} (M_{⋆}, z)}{Δ log L_{n}}, \end{matrix}$ $\begin{aligned} p(L_{n}\ |\ M_{m},\ z) = \frac{N_{\mathrm{AGN},mn}}{N_{\mathrm{gal},m}\Delta \log L_{n}} = \frac{f_{\rm AGN}(M_{\star },\ z)}{\Delta \log L_{n}}, \end{aligned}$ (21)

where N_AGN, mn is the number of the Radio-excess AGN Sample in the mth stellar mass bin and the nth luminosity bin, N_gal, m is the number of the All Galaxies Sample in the mth stellar mass bin, and ΔlogL_n is the width of the nth luminosity bin. It is worth noting that when calculating the probability of SFGs (or QGs, or all galaxies) hosting a radio-excess AGN, here N_AGN, mn refers to all the radio-excess AGNs hosted by SFGs (or QGs, or all galaxies) and N_gal refers to all the SFGs (or QGs, or all galaxies) in the All Galaxies Sample. The estimates of p(L_R | M_⋆, z) for SFGs in the six redshift bins are shown as colored symbols in Figs. F.1 and F.2. The radio-excess AGNs in both QGs and all galaxies exhibit similar trends to SFGs, so we do not show their details here. The probability of a galaxy hosting a radio-excess AGN with a given L_R increases with M_⋆ at all redshift bins (see Fig. F.1), which is consistent with the previous works (e.g., Best et al. 2005; Sabater et al. 2019). The probability of a galaxy with a given M_⋆ hosting a radio-excess AGN decreases with L_R, which is also consistent with previous works (e.g., Best et al. 2005). These trends of radio-excess AGNs in this work are consistent with those of X-ray AGNs (e.g., Haggard et al. 2010; Aird et al. 2012; Wang et al. 2017).

6.2. Simple χ² fits in each redshift bin

To quantitatively study the trends of p(L_R | M_⋆, z), we apply χ² fits to the results calculated by Eq. (21) (see the data points in Figs. F.1 and F.2). Next we take the analysis for SFGs as an example to show the detailed analysis procedures. At each fixed L_R, we assume a simple power-law relation for p(L_R | M_⋆, z) as a function of M_⋆,

$\begin{matrix} log [p (L_{R} | M_{⋆}, z)] = a + b log [\frac{M_{⋆}}{10^{11} M_{⊙}}], \end{matrix}$ $\begin{aligned} \log [p(L_{\rm R}\ |\ M_{\star },\ z)] = a + b \log \left[\frac{M_{\star }}{10^{11}\ M_{\odot }}\right], \end{aligned}$ (22)

where a and b are the intercept and slope of the relation, respectively. The best-fit values for a and b in each L_R bin and each redshift bin are presented in the right column of Fig. F.1. The averaged best-fit reduced-χ² over all the redshift bins and all the L_R bins is 1.22. In all redshift bins, the intercept a decreases with higher L_R, while the slope b does not show significant changes. Similarly, for each fixed M_⋆, p(L_R | M_⋆, z) as a function of L_R can be defined as

$\begin{matrix} log [p (L_{R} | M_{⋆}, z)] = c + d log [\frac{L_{R}}{10^{23} W {Hz}^{- 1}}], \end{matrix}$ $\begin{aligned} \log [p(L_{\rm R}\ |\ M_{\star },\ z)] = c + d \log \left[\frac{L_{\rm R}}{10^{23}\ \mathrm {W}\ \mathrm {Hz}^{-1}}\right], \end{aligned}$ (23)

where c and d are the intercept and slope of the relation, respectively. The best-fit parameters in each M_⋆ bin and each redshift bin are shown in the right column of Fig. F.2. The averaged best-fit reduced-χ² over all the redshift bins and all the M_⋆ bins is 1.76. In all redshift bins, the intercept c decreases with higher M_⋆, while the slope d nearly keeps constant. These results indicate that the slope b (or the slope d) of p(L_R | M_⋆, z) as a function of M_⋆ (or L_R) is independent of L_R (or M_⋆). The QGs and all galaxies show similar trends but different best-fit parameters compared to SFGs. Given that we mainly focus on the trends in this section, for brevity we do not discuss their results in details here. The main differences among different galaxy populations are discussed in Sect. 6.3.4 and shown in Fig. 7.

Fig. 7.

Best-fit parameters of the maximum-likelihood fitting for p(L_R | M_⋆, z) of SFGs (blue circles and blue dashed lines), QGs (red squares and red dash-dotted lines), and all galaxies (black stars and black solid lines). The colored symbols represent the best-fit parameters obtained by Eqs. (24)–(27) in each redshift bin (see details in Sects. 6.3.1 and 6.3.2). The colored lines show the best-fit results over the entire redshift range obtained by Eqs. (28)–(29) and (25)–(27) (see details in Sect. 6.3.3). The dotted gray lines represent the results for the X-ray AGNs in Aird et al. (2012).

6.3. Redshift evolution

The above analysis starts by dividing sources into different redshift bins, then binning the sources according to their M_⋆ and L_R within each redshift bin, and finally fitting the binned data. Therefore, the best-fit result may be significantly affected by the bin size of redshift, M_⋆, and L_R, and the information of each individual galaxy cannot be fully utilized due to the binning. Therefore, next we utilize the maximum-likelihood fitting approach in Aird et al. (2012) to measure the dependence of p(L_R | M_⋆, z) on both M_⋆ and L_R. First, in Sect. 6.3.1, we introduce the fitting approach applied to our samples. Then we test this approach in each redshift bin in Sect. 6.3.2. Finally, in Sect. 6.3.3 we apply this approach incorporating the redshift evolution in order to remove the effects brought by the redshift binning.

6.3.1. Maximum-likelihood fitting

According to the results in Sect. 6.2, the slope b (or d) in Eq. (22) (or (23)) of p(L_R | M_⋆, z) as a function of M_⋆ (or L_R) is independent of L_R (or M_⋆). It means that p(L_R | M_⋆, z) at a fixed redshift can be expressed as a separable function of M_⋆ and L_R in the form of

$\begin{matrix} p (L_{R} | M_{⋆}, z) d log L_{R} = K {(\frac{M_{⋆}}{M_{0}})}^{γ_{M}} {(\frac{L_{R}}{L_{0}})}^{γ_{L}} d log L_{R}, \end{matrix}$ $\begin{aligned} p(L_{\rm R}\ |\ M_{\star },\ z)\ d\log L_{\rm R} = K\left(\frac{M_{\star }}{M_0}\right)^{\gamma _{\rm M}}\left(\frac{L_{\rm R}}{L_0}\right)^{\gamma _{\rm L}} d \log L_{\rm R},\\ \end{aligned}$ (24)

where the scaling factor M₀ is set to 10¹¹ M_⊙, the scaling factor L₀ is set to 10²³ W Hz⁻¹, K is the normalization, γ_M and γ_L are indices. The best-fit parameters are found through maximizing the log-likelihood function,

$\begin{matrix} ln L = - N + \sum_{k = 1}^{N_{i}^{AGN}} ln p_{k}, \end{matrix}$ $\begin{aligned} \ln \mathcal{L} = - \mathcal{N} + \sum _{k = 1}^{N_i^\mathrm{AGN}} \ln p_k, \end{aligned}$ (25)

where $N_{i}^{AGN}$ $N_i^{\rm AGN}$ is the number of radio-excess AGNs in the GOODS-N and COSMOS/UltraVISTA fields in the ith redshift bin, p_k is the probability of a galaxy with stellar mass M_k hosting a radio-excess AGN with 1.4 GHz luminosity L_k:

$\begin{matrix} p_{k} = p (L_{k} | M_{k}, z_{k}) = K {(\frac{M_{k}}{M_{0}})}^{γ_{M}} {(\frac{L_{k}}{L_{0}})}^{γ_{L}} . \end{matrix}$ $\begin{aligned} p_k = p(L_k\ |\ M_k,\ z_k) = K\left(\frac{M_k}{M_0}\right)^{\gamma _{\rm M}}\left(\frac{L_k}{L_0}\right)^{\gamma _{\rm L}}. \end{aligned}$ (26)

The expected number of radio-excess AGNs, 𝒩, is defined as

$\begin{matrix} N = \sum_{j = 1}^{N_{i}^{gal, gn}} \int_{log L_{\lim, gn, j}}^{\infty} p (L_{R} | M_{j}, z_{j}) d log L_{R} \\ + \sum_{j = 1}^{N_{i}^{gal, cs}} \int_{log L_{\lim, cs, j}}^{\infty} p (L_{R} | M_{j}, z_{j}) d log L_{R} \end{matrix}$ $\begin{aligned}&\mathcal{N} = \sum _{j = 1}^{N_i^\mathrm{gal,gn}} \int _{\log L_{\mathrm{lim,gn},j}}^{\infty } p(L_{\rm R}\ |\ M_j,\ z_j)\ d \log L_{\rm R} \nonumber \\&\qquad + \sum _{j = 1}^{N_i^\mathrm{gal,cs}} \int _{\log L_{\mathrm{lim,cs},j}}^{\infty } p(L_{\rm R}\ |\ M_j,\ z_j)\ d \log L_{\rm R} \end{aligned}$ (27)

$\begin{matrix} = \sum_{j = 1}^{N_{i}^{gal, gn}} K {(\frac{M_{j}}{M_{0}})}^{γ_{M}} \int_{log L_{\lim, gn, j}}^{\infty} {(\frac{L_{R}}{L_{0}})}^{γ_{L}} d log L_{R}, \\ + \sum_{j = 1}^{N_{i}^{gal, cs}} K {(\frac{M_{j}}{M_{0}})}^{γ_{M}} \int_{log L_{\lim, cs, j}}^{\infty} {(\frac{L_{R}}{L_{0}})}^{γ_{L}} d log L_{R}, \end{matrix}$ $\begin{aligned}&\ \ = \sum _{j = 1}^{N_i^\mathrm{gal,gn}} K \left(\frac{M_j}{M_0}\right)^{\gamma _{\rm M}} \int _{\log L_{\mathrm{lim,gn},j}}^{\infty } \left(\frac{L_{\rm R}}{L_0}\right)^{\gamma _{\rm L}}\ d\log L_{\rm R}, \nonumber \\&\qquad + \sum _{j = 1}^{N_i^\mathrm{gal,cs}} K \left(\frac{M_j}{M_0}\right)^{\gamma _{\rm M}} \int _{\log L_{\mathrm{lim,cs},j}}^{\infty } \left(\frac{L_{\rm R}}{L_0}\right)^{\gamma _{\rm L}}\ d\log L_{\rm R},\nonumber \end{aligned}$ (28)

where $N_{i}^{gal, gn}$ $N_i^{\rm gal,gn}$ and $N_{i}^{gal, cs}$ $N_i^{\rm gal,cs}$ are the numbers of galaxies in the ith redshift bin for the GOODS-N and COSMOS/UltraVISTA fields, respectively. The L_lim of the GOODS-N field (L_{lim, gn, j}) corresponds to the 5σ detection limit of the radio survey at redshift z_j (see Fig. 2), while for the COSMOS (L_{lim, cs, j}) field, it corresponds to the L_1.4 GHz threshold for correcting the AGN classification purity (see details in Appendix B). The 1σ error of each model parameter is estimated by the maximum projection of Δ𝒮 = 1.0 where 𝒮 = −2lnℒ according to Aird et al. (2012).

6.3.2. Results in each redshift bin

In each redshift bin we use the above maximum-likelihood fitting approach to constrain p(L_R | M_⋆, z) for SFGs, QGs, and all galaxies. The best-fit values for γ_M, γ_L, and K in each redshift bin are shown in the form of colored symbols in Fig. 7. For both SFGs and QGs, we find that γ_M and γ_L almost keep constant at all redshift bins except at z ≲ 0.5 (see detailed discussions in Sect. 6.3.3). For all galaxies, γ_M and γ_L do not show any significant changes with redshift (0.1 < z ≤ 4.0), which is consistent with the trends of X-ray AGNs in Aird et al. (2012) at 0.2 < z ≤ 1.0. For SFGs, QGs, and all galaxies, the normalization K significantly increases with redshift. In addition, SFGs and QGs show significantly different values of γ_M, γ_L, and K, which will be discussed in detail in Sect. 6.3.4. Overall, the probability of a galaxy with a given stellar mass hosting a radio-excess AGN with a given radio luminosity increases toward higher redshift. This indicates higher AGN activities and potentially more significant AGN feedbacks toward higher redshift.

6.3.3. Incorporating redshift evolution into the maximum-likelihood fitting

According to the results of Sect. 6.3.2, both γ_M and γ_L are nearly constant over the entire redshift range, which indicates that the redshift evolution of p(L_R | M_⋆, z) can be assumed to be independent of L_R and M_⋆. Only the normalization K shows an evident evolution with redshift, so p(L_R | M_⋆, z) incorporating the redshift evolution for K with a simple power law can be rewritten as

$\begin{matrix} p (L_{R} | M_{⋆}, z) d log L_{R} = A {(\frac{1 + z}{1 + z_{0}})}^{γ_{z}} {(\frac{M_{⋆}}{M_{0}})}^{γ_{M}} {(\frac{L_{R}}{L_{0}})}^{γ_{L}} d log L_{R}, \end{matrix}$ $\begin{aligned} p(L_{\rm R}\ |\ M_{\star },\ z)\ d\log L_{\rm R} = A\left(\frac{1+z}{1+z_0}\right)^{\gamma _z} \left(\frac{M_{\star }}{M_0}\right)^{\gamma _{\rm M}}\left(\frac{L_{\rm R}}{L_0}\right)^{\gamma _{\rm L}} \ d \log L_{\rm R}, \end{aligned}$ (29)

where the scaling factor z₀ is set to be 2.0 (median of the redshift range for our sample), A is the overall normalization. Thus,

$\begin{matrix} K = A {(\frac{1 + z}{1 + z_{0}})}^{γ_{z}} . \end{matrix}$ $\begin{aligned} K = A\left(\frac{1+z}{1+z_0}\right)^{\gamma _z}. \end{aligned}$ (30)

Then we modify p(L_R | M_⋆, z) in Eqs. (25)–(27) with Eq. (28) to perform maximum-likelihood fitting over the entire redshift range (0.1 < z ≤ 4.0). The best-fit parameters are listed in Table 2 and plotted in Fig. 7 as the colored lines. The redshift evolution model over the entire redshift range is almost consistent with the result obtained in each redshift bin (colored symbols in Fig. 7).

Table 2.

Best-fit parameters from maximum-likelihood fitting with redshift evolution.

To compare the model with the observational data, we use the N_obs/N_mdl method of Aird et al. (2012; see Eq. (14) in their paper) to scale the observed radio-excess AGN number (N_obs) by the expected number estimated by the model (N_mdl) in each redshift bin (see Eq. (28)). Overall, the binned estimates for the observational data in each redshift bin are well described by the best-fit model obtained with maximum-likelihood fitting over the entire redshift range (see results for SFGs in Fig. 8, for QGs in Fig. G.1, for all galaxies in Fig. G.2). However, the best-fit model cannot fully explain the data in the following cases: low-luminosity radio-excess AGNs in SFGs at 0.1 < z ≤ 1.0 (L_1.4GHz < 10²³ W Hz⁻¹; see Panel A1 in Fig. 8) and at 2 < z ≤ 4 (L_1.4GHz < 10^23.5 W Hz⁻¹; see Panel D1 in Fig. 8); radio-excess AGNs in massive SFGs (M_⋆ > 10^11.5 M_⊙) at 0.1 < z ≤ 1.0 (see Panel A2 in Fig. 8); radio-excess AGNs in low-mass QGs (M_⋆ < 10¹⁰ M_⊙) at 0.1 < z ≤ 1.0 (see Panel A2 in Fig. G.1). The following reasons might explain the slight differences between the observational data and model predictions. In the analysis we assume that γ_M (or γ_L) is independent of L_R (or M_⋆) and redshift (see Sects. 6.2 and 6.3.2). On the one hand, a mild dependence of γ_M (or γ_L) on L_R (or M_⋆) cannot be fully ruled out. In addition, γ_M may be different between low-mass and high-mass galaxies (Williams & Röttgering 2015; Zhu et al. 2023), while γ_L may also change from low-luminosity to high-luminosity populations (Best et al. 2005). On the other hand, a weak evolution of both γ_M and γ_L with redshift might exist (see colored symbols in Fig. 7 and discussion in Sect. 7). Some works had found that γ_M decreases from the local universe to z ∼ 2 (Williams & Röttgering 2015; Zhu et al. 2023), while Kondapally et al. (2022) found a weak positive evolution with redshift (for SFGs only). Even so, on the whole our model can still well describe the redshift evolution of radio-excess AGN fraction across the entire redshift range.

Fig. 8.

Estimates of p(L_R | M_⋆, z) as a function of M_⋆ and L_R based on our maximum-likelihood fitting results for SFGs. The colored dashed lines in the top and bottom rows represent our best-fit model from the unbinned maximum-likelihood fitting (see Eq. (28) and Table 2) through combining data in the GOODS-N and COSMOS/UltraVISTA fields over the entire redshift range (0.1 < z < 4) evaluated at the center of each redshift bin. The binned data points are scaled with the probability estimated by the model (see Eq. (28)) using the N_obs/N_mdl method of Aird et al. (2012) (see details in Sect. 6.3.3). In the top panel, the different colors represent different M_⋆ bins. In the bottom panel, the different colors show different L_R bins.

6.3.4. Quantitative relation of radio-excess AGN fraction in SFGs and QGs

To sum up, we obtain quantitative relations for the probability of SFGs, QGs, or all galaxies hosting a radio-excess AGN as a function of M_⋆, L_R, and redshift, respectively:

$\begin{matrix} p {(L_{R} | M_{⋆}, z)}_{SFG} = 10^{- 0.79} {(\frac{M_{⋆}}{M_{0}})}^{1.06} {(\frac{L_{R}}{L_{0}})}^{- 0.77} {(\frac{1 + z}{1 + z_{0}})}^{3.08}, \\ p {(L_{R} | M_{⋆}, z)}_{QG} = 10^{- 0.70} {(\frac{M_{⋆}}{M_{0}})}^{1.41} {(\frac{L_{R}}{L_{0}})}^{- 0.60} {(\frac{1 + z}{1 + z_{0}})}^{2.47}, \\ p {(L_{R} | M_{⋆}, z)}_{All} = 10^{- 0.88} {(\frac{M_{⋆}}{M_{0}})}^{1.10} {(\frac{L_{R}}{L_{0}})}^{- 0.67} {(\frac{1 + z}{1 + z_{0}})}^{2.24} . \end{matrix}$ $\begin{aligned}&p(L_{\rm R}\ |\ M_{\star },\ z)_{\rm SFG} = 10^{-0.79} \left(\frac{M_{\star }}{M_0}\right)^{1.06} \left(\frac{L_{\rm R}}{L_0}\right)^{-0.77} \left(\frac{1+z}{1+z_0}\right)^{3.08},\nonumber \\&p(L_{\rm R}\ |\ M_{\star },\ z)_{\rm QG} = 10^{-0.70} \left(\frac{M_{\star }}{M_0}\right)^{1.41} \left(\frac{L_{\rm R}}{L_0}\right)^{-0.60} \left(\frac{1+z}{1+z_0}\right)^{2.47},\nonumber \\&p(L_{\rm R}\ |\ M_{\star },\ z)_{\rm All} = 10^{-0.88} \left(\frac{M_{\star }}{M_0}\right)^{1.10} \left(\frac{L_{\rm R}}{L_0}\right)^{-0.67}\ \left(\frac{1+z}{1+z_0}\right)^{2.24}. \end{aligned}$ (31)

Here M₀ is set to 10¹¹ M_⊙, L₀ is set to 10²³ W Hz⁻¹, and z₀ is set to 2.0. These best-fit parameters and their uncertainties are summarized in Table 2. The SFGs and QGs have significantly different evolution trends for radio-excess AGN fraction (see Fig. 7). The QGs have larger γ_M and γ_L than the SFGs, which indicates that radio-excess AGNs in QGs are most often found to reside in more massive galaxies with higher radio luminosity compared to those in SFGs. Both SFGs and QGs have an increasing radio-excess AGN fraction with redshift (see Fig. 7). SFGs show a larger γ_z than QGs (see Table 2), which implies that radio-excess AGN fractions in SFGs increase more rapidly toward higher redshifts than those in QGs.

7. Discussion: Comparison with previous works

As Eq. (30) shows, we obtain a smaller γ_M (∼1) than that in Best et al. (2005) (γ_M ∼ 2.5; 0.03 < z < 0.3). Our result is consistent with that in Williams & Röttgering (2015) at 1.5 < z < 2 (γ_M ∼ 1.0), while their results showed a decreasing γ_M trend from the local universe (γ_M ∼ 2.7) to z ∼ 2. Zhu et al. (2023) also found a decreasing stellar-mass dependence from z ∼ 0.3 to z ∼ 2.3. In addition, the dependence of radio-AGN fraction on the stellar mass may change when M_* > 10¹¹ M_⊙ (Williams & Röttgering 2015; Zhu et al. 2023). The differences between our results and theirs may be due to the different radio-AGN selections and different calculation methods. Previous works usually focused on powerful radio-AGNs (L_1.4 GHz > 10²³ or 10²⁴ W Hz⁻¹), while our sample additionally includes some faint radio-AGNs with L_1.4 GHz < 10²³ W Hz⁻¹. Moreover, the different evolution trends of radio-AGN in SFGs and QGs indicate that disregarding the galaxy types may have a nonnegligible effect on statistical results. Kondapally et al. (2022) proposed for the first time a radio-AGNs study for QGs and SFGs separately. They found that low-excitation radio galaxies in QGs have a larger γ_M (∼2.5) than those in SFGs (γ_M ∼ 1.37). We found a similar trend but lower γ_M values than their results (γ_M ∼ 1.4 for QGs and γ_M ∼ 1.0 for SFGs in our work). One possible reason for this difference may be that we focus on different radio-AGNs samples compared to Kondapally et al. (2022) (see details in Sect. 5.3). Our results show that γ_M has no significant evolution with redshift at 0.1 < z ≤ 4.0, while for LERGs in Kondapally et al. (2022) at 0.3 < z ≤ 1.5, γ_M shows no significant evolution for QGs and a positive evolution for SFGs.

Given that we use the same calculation method as Aird et al. (2012) (see details in Sect. 6.3), we also compare the results for radio-excess AGNs with those for X-ray AGNs from Aird et al. (2012), which are plotted as gray dotted lines in Fig. 7. Compared to X-ray AGNs (Aird et al. 2012), the evolutions of radio-excess AGN fraction (see black solid lines in Fig. 7) have a larger γ_M, which indicates that radio-excess AGNs tend to reside in more massive galaxies (see Fig. 7). X-ray AGN fractions in red sequence galaxies exhibit a more rapid redshift evolution than those in blue cloud galaxies (Aird et al. 2012; Wang et al. 2017), while the radio-excess AGN fraction in SFGs increases more rapidly toward higher redshift than those in QGs. These different evolution trends between radio-excess AGNs and X-ray AGNs suggest that black holes with different accretion states may influence their host galaxies in different modes.

In addition, we also examine whether M_⋆-dependent IRRC (Delvecchio et al. 2021) affects the results about the probability of hosting a radio-excess AGN. Following van der Vlugt et al. (2022), we select radio-excess AGNs with q_TIR deviating more than 3σ from the M_⋆-dependent IRRC of Delvecchio et al. (2021). Then we recalculate the probability of a radio-excess AGN hosted by SFGs using Eq. (28) (see the dark blue region in Fig. H.1), which is consistent with the results based on the radio-excess-AGN selection method in Sect. 4 (see the green region in Fig. H.1). This result indicates that using the M_*-dependent IRRC (Delvecchio et al. 2021) does not alter the results in this work.

8. Summary

We used the optical to MIR surveys to select an All Galaxies Sample with ∼400 000 sources at 0.1 < z < 4 in GOODS-N and COSMOS/UltraVISTA fields (totaling 1.6 deg²; Sect. 2). After cross-matching with the deep–large radio surveys (1.4 GHz or 3 GHz) in these fields, we selected 7494 radio sources with S/N of radio flux ≥5 at 0.1 < z < 4 as our Radio Sources Sample (Sect. 2). Combining this sample with the de-blended IR photometry in these fields, we calculated the infrared-to-radio ratio (q_TIR) distributions of these radio sources to select a Radio-excess AGN Sample with 983 sources (Sect. 4.1). Next we subdivided all our samples into SFGs and QGs using the UVJ method. Based on the Radio Sources Sample, we constructed 1.4 GHz radio luminosity functions (RLFs) for SFGs and radio-excess AGNs at 0.1 < z < 4 (Sect. 5.1), respectively, and studied their evolutions with redshift. Based on the All Galaxies Sample and Radio-excess AGN Sample, we further investigated the probability of different galaxies populations (SFGs and QGs) hosting a radio-excess AGN as a function of stellar mass, radio luminosity, and redshift. The main conclusions we obtain are shown as follows:

The q_TIR value of SFGs shows a weak evolution with redshift: ${\bar{q}}_{TIR} = 2.62 \times {(1 + z)}^{- 0.08}$ $\bar{q}_\mathrm{{TIR}} = 2.62\times (1+z)^{-0.08}$ (Sect. 4.2), which is generally consistent with previous works.
The evolution of RLFs with redshift for SFGs can be well described by the pure luminosity evolution model of L_⋆ ∝ (1 + z)^{−0.34 × z + 3.57} (Sect. 5.2). The evolution of AGN RLFs follows the pure density evolution model of Φ_⋆ ∝ (1 + z)^{−0.77 × z + 2.69} (Sect. 5.3).
The evolution of crossover luminosity between SFG RLFs and AGN RLFs is shown as L_D = 10^22.9 × (1 + z)^{0.06 × z + 3.05} W Hz⁻¹, which can be used to select powerful radio-AGNs at different redshifts solely through radio surveys (Sect. 5.4). This result also indicates a decreasing contribution of AGNs to entire radio populations toward higher redshift.
The probability of a galaxy hosting a radio-excess AGN is shown as a function of stellar mass, radio luminosity, and redshift: $p (L_{R} | M_{⋆}, z) = 10^{- 0.79} {(\frac{M_{⋆}}{M_{0}})}^{1.06} {(\frac{L_{R}}{L_{0}})}^{- 0.77} {(\frac{1 + z}{1 + z_{0}})}^{3.08}$ $p(L_{\mathrm{R}}\ |\ M_{\star},\ z) = 10^{-0.79} \left(\frac{M_{\star}}{M_0}\right)^{1.06} \left(\frac{L_{\mathrm{R}}}{L_0}\right)^{-0.77} \left(\frac{1+z}{1+z_0}\right)^{3.08}$ for SFGs, $p (L_{R} | M_{⋆}, z) = 10^{- 0.70} {(\frac{M_{⋆}}{M_{0}})}^{1.41} {(\frac{L_{R}}{L_{0}})}^{- 0.60} {(\frac{1 + z}{1 + z_{0}})}^{2.47}$ $p(L_{\mathrm{R}}\ |\ M_{\star},\ z) = 10^{-0.70} \left(\frac{M_{\star}}{M_0}\right)^{1.41} \left(\frac{L_{\mathrm{R}}}{L_0}\right)^{-0.60} \left(\frac{1+z}{1+z_0}\right)^{2.47}$ for QGs (Sect. 6.3). It indicates that radio-excess AGNs in QGs reside in more massive galaxies with higher radio luminosity than those in SFGs. The fractions of radio-excess AGNs in both SFGs and QGs are increasing from the local universe to the higher redshift. In addition, this increasing trend in SFGs is more significant than in QGs.

The above studies can lay the foundation for further investigations with the upcoming revolutionary radio facilities, such as the Square Kilometer Array (SKA; Dewdney et al. 2009; Norris et al. 2013; McAlpine et al. 2015) and the Next Generation Very Large Array (ngVLA; Hughes et al. 2015).

Acknowledgments

We thank the referee for constructive comments that greatly improved this paper. This work is supported by the National Natural Science Foundation of China (Project No. 12173017 and Key Project No. 12141301).

References

Aird, J., Nandra, K., Laird, E. S., et al. 2010, MNRAS, 401, 2531 [Google Scholar]
Aird, J., Coil, A. L., Moustakas, J., et al. 2012, ApJ, 746, 90 [CrossRef] [Google Scholar]
Alberts, S., Rujopakarn, W., Rieke, G. H., Jagannathan, P., & Nyland, K. 2020, ApJ, 901, 168 [NASA ADS] [CrossRef] [Google Scholar]
An, F., Vaccari, M., Smail, I., et al. 2021, MNRAS, 507, 2643 [NASA ADS] [CrossRef] [Google Scholar]
Appleton, P. N., Fadda, D. T., Marleau, F. R., et al. 2004, ApJS, 154, 147 [Google Scholar]
Arnouts, S., Moscardini, L., Vanzella, E., et al. 2002, MNRAS, 329, 355 [Google Scholar]
Baldwin, J. A., Phillips, M. M., & Terlevich, R. 1981, PASP, 93, 5 [Google Scholar]
Bariuan, L. G. C., Snios, B., Sobolewska, M., Siemiginowska, A., & Schwartz, D. A. 2022, MNRAS, 513, 4673 [NASA ADS] [CrossRef] [Google Scholar]
Barro, G., Pérez-González, P. G., Cava, A., et al. 2019, ApJS, 243, 22 [NASA ADS] [CrossRef] [Google Scholar]
Bell, E. F. 2003, ApJ, 586, 794 [Google Scholar]
Best, P. N., Kauffmann, G., Heckman, T. M., et al. 2005, MNRAS, 362, 25 [Google Scholar]
Best, P. N., Kondapally, R., Williams, W. L., et al. 2023, MNRAS, 523, 1729 [NASA ADS] [CrossRef] [Google Scholar]
Bonchi, A., La Franca, F., Melini, G., Bongiorno, A., & Fiore, F. 2013, MNRAS, 429, 1970 [CrossRef] [Google Scholar]
Boquien, M., Burgarella, D., Roehlly, Y., et al. 2019, A&A, 622, A103 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Brammer, G. B., van Dokkum, P. G., & Coppi, P. 2008, ApJ, 686, 1503 [Google Scholar]
Brown, M. J. I., Jannuzi, B. T., Floyd, D. J. E., & Mould, J. R. 2011, ApJ, 731, L41 [Google Scholar]
Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000 [NASA ADS] [CrossRef] [Google Scholar]
Burgarella, D., Buat, V., & Iglesias-Páramo, J. 2005, MNRAS, 360, 1413 [NASA ADS] [CrossRef] [Google Scholar]
Calistro Rivera, G., Williams, W. L., Hardcastle, M. J., et al. 2017, MNRAS, 469, 3468 [Google Scholar]
Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682 [NASA ADS] [CrossRef] [Google Scholar]
Ceraj, L., Smolčić, V., Delvecchio, I., et al. 2018, A&A, 620, A192 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Chabrier, G. 2003, ApJ, 586, L133 [NASA ADS] [CrossRef] [Google Scholar]
Cochrane, R. K., Kondapally, R., Best, P. N., et al. 2023, MNRAS, 523, 6082 [NASA ADS] [CrossRef] [Google Scholar]
Condon, J. J. 1992, ARA&A, 30, 575 [Google Scholar]
Condon, J. J., & Dressel, L. L. 1978, ApJ, 221, 456 [NASA ADS] [CrossRef] [Google Scholar]
Condon, J. J., Matthews, A. M., & Broderick, J. J. 2019, ApJ, 872, 148 [NASA ADS] [CrossRef] [Google Scholar]
Dale, D. A., Helou, G., Magdis, G. E., et al. 2014, ApJ, 784, 83 [Google Scholar]
Davé, R., Anglés-Alcázar, D., Narayanan, D., et al. 2019, MNRAS, 486, 2827 [Google Scholar]
Delhaize, J., Smolčić, V., Delvecchio, I., et al. 2017, A&A, 602, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Del Moro, A., Alexander, D. M., Mullaney, J. R., et al. 2013, A&A, 549, A59 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Delvecchio, I., Smolčić, V., Zamorani, G., et al. 2017, A&A, 602, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Delvecchio, I., Daddi, E., Sargent, M. T., et al. 2021, A&A, 647, A123 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Delvecchio, I., Daddi, E., Sargent, M. T., et al. 2022, A&A, 668, A81 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Dewdney, P. E., Hall, P. J., Schilizzi, R. T., & Lazio, T. J. L. W. 2009, IEEE Proc., 97, 1482 [Google Scholar]
Donley, J. L., Rieke, G. H., Rigby, J. R., & Pérez-González, P. G. 2005, ApJ, 634, 169 [Google Scholar]
Donley, J. L., Koekemoer, A. M., Brusa, M., et al. 2012, ApJ, 748, 142 [Google Scholar]
Donoso, E., Best, P. N., & Kauffmann, G. 2009, MNRAS, 392, 617 [Google Scholar]
Donoso, E., Li, C., Kauffmann, G., Best, P. N., & Heckman, T. M. 2010, MNRAS, 407, 1078 [NASA ADS] [CrossRef] [Google Scholar]
Draine, B. T., & Li, A. 2007, ApJ, 657, 810 [CrossRef] [Google Scholar]
Dubner, G., & Giacani, E. 2015, A&A Rev., 23, 3 [NASA ADS] [CrossRef] [Google Scholar]
Dullo, B. T., Knapen, J. H., Beswick, R. J., et al. 2023, MNRAS, 522, 3412 [NASA ADS] [CrossRef] [Google Scholar]
Enia, A., Talia, M., Pozzi, F., et al. 2022, ApJ, 927, 204 [NASA ADS] [CrossRef] [Google Scholar]
Fabian, A. C. 2012, ARA&A, 50, 455 [Google Scholar]
Fiore, F., Feruglio, C., Shankar, F., et al. 2017, A&A, 601, A143 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
Franzen, T. M. O., Seymour, N., Sadler, E. M., et al. 2021, PASA, 38, e041a [NASA ADS] [CrossRef] [Google Scholar]
Gim, H. B., Yun, M. S., Owen, F. N., et al. 2019, ApJ, 875, 80 [NASA ADS] [CrossRef] [Google Scholar]
Grogin, N. A., Kocevski, D. D., Faber, S. M., et al. 2011, ApJS, 197, 35 [NASA ADS] [CrossRef] [Google Scholar]
Gruppioni, C., Pozzi, F., Rodighiero, G., et al. 2013, MNRAS, 432, 23 [Google Scholar]
Haggard, D., Green, P. J., Anderson, S. F., et al. 2010, ApJ, 723, 1447 [NASA ADS] [CrossRef] [Google Scholar]
Hardcastle, M. J., & Croston, J. H. 2020, New A Rev., 88, 101539a [NASA ADS] [CrossRef] [Google Scholar]
Helou, G., Soifer, B. T., & Rowan-Robinson, M. 1985, ApJ, 298, L7 [Google Scholar]
Ho, L. C. 2008, ARA&A, 46, 475 [Google Scholar]
Hughes, A. M., Beasley, A., & Carilli, C. 2015, IAU Gener. Assembly, 29, 2255106 [Google Scholar]
Ilbert, O., Arnouts, S., McCracken, H. J., et al. 2006, A&A, 457, 841 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Ivison, R. J., Magnelli, B., Ibar, E., et al. 2010, A&A, 518, L31 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Janssen, R. M. J., Röttgering, H. J. A., Best, P. N., & Brinchmann, J. 2012, A&A, 541, A62 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Jarvis, M., Taylor, R., Agudo, I., et al. 2016, in MeerKAT Science: On the Pathway to the SKA, 6 [Google Scholar]
Jin, S., Daddi, E., Liu, D., et al. 2018, ApJ, 864, 56 [Google Scholar]
Kauffmann, G., & Heckman, T. M. 2009, MNRAS, 397, 135 [NASA ADS] [CrossRef] [Google Scholar]
Kauffmann, G., Heckman, T. M., Tremonti, C., et al. 2003, MNRAS, 346, 1055 [Google Scholar]
Kewley, L. J., Maier, C., Yabe, K., et al. 2013, ApJ, 774, L10 [Google Scholar]
Koekemoer, A. M., Faber, S. M., Ferguson, H. C., et al. 2011, ApJS, 197, 36 [NASA ADS] [CrossRef] [Google Scholar]
Kolwa, S., Jarvis, M. J., McAlpine, K., & Heywood, I. 2019, MNRAS, 482, 5156 [CrossRef] [Google Scholar]
Kondapally, R., Best, P. N., Cochrane, R. K., et al. 2022, MNRAS, 513, 3742 [NASA ADS] [CrossRef] [Google Scholar]
Kondapally, R., Best, P. N., Raouf, M., et al. 2023, MNRAS, 523, 5292 [NASA ADS] [CrossRef] [Google Scholar]
Kriek, M., van Dokkum, P. G., Labbé, I., et al. 2009, ApJ, 700, 221 [Google Scholar]
Kukreti, P., Morganti, R., Tadhunter, C., & Santoro, F. 2023, A&A, 674, A198 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Lacy, M., Storrie-Lombardi, L. J., Sajina, A., et al. 2004, ApJS, 154, 166 [Google Scholar]
Liu, D., Daddi, E., Dickinson, M., et al. 2018, ApJ, 853, 172 [Google Scholar]
Luo, B., Brandt, W. N., Xue, Y. Q., et al. 2017, ApJS, 228, 2 [Google Scholar]
Magliocchetti, M. 2022, A&A Rev., 30, 6 [NASA ADS] [CrossRef] [Google Scholar]
Magliocchetti, M., Maddox, S. J., Hawkins, E., et al. 2004, MNRAS, 350, 1485 [NASA ADS] [CrossRef] [Google Scholar]
Magliocchetti, M., Popesso, P., Brusa, M., et al. 2017, MNRAS, 464, 3271 [NASA ADS] [CrossRef] [Google Scholar]
Magnelli, B., Ivison, R. J., Lutz, D., et al. 2015, A&A, 573, A45 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Malavasi, N., Bardelli, S., Ciliegi, P., et al. 2015, A&A, 576, A101 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Matthews, A. M., Condon, J. J., Cotton, W. D., & Mauch, T. 2021, ApJ, 914, 126 [CrossRef] [Google Scholar]
Matzeu, G. A., Brusa, M., Lanzuisi, G., et al. 2023, A&A, 670, A182 [Google Scholar]
Mauch, T., & Sadler, E. M. 2007, MNRAS, 375, 931 [Google Scholar]
McAlpine, K., Jarvis, M. J., & Bonfield, D. G. 2013, MNRAS, 436, 1084 [CrossRef] [Google Scholar]
McAlpine, K., Prandoni, I., Jarvis, M., et al. 2015, in Advancing Astrophysics with the Square Kilometre Array (AASKA14), 83 [Google Scholar]
Merloni, A., Heinz, S., & di Matteo, T. 2003, MNRAS, 345, 1057 [Google Scholar]
Molnár, D. C., Sargent, M. T., Leslie, S., et al. 2021, MNRAS, 504, 118 [CrossRef] [Google Scholar]
Morrison, G. E., Owen, F. N., Dickinson, M., Ivison, R. J., & Ibar, E. 2010, ApJS, 188, 178 [Google Scholar]
Mullaney, J. R., Alexander, D. M., Goulding, A. D., & Hickox, R. C. 2011, MNRAS, 414, 1082 [Google Scholar]
Ni, Q., Aird, J., Merloni, A., et al. 2023, MNRAS, 524, 4778 [CrossRef] [Google Scholar]
Noll, S., Burgarella, D., Giovannoli, E., et al. 2009, A&A, 507, 1793 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Norris, R. P., Afonso, J., Bacon, D., et al. 2013, PASA, 30, e020a [NASA ADS] [CrossRef] [Google Scholar]
Novak, M., Smolčić, V., Delhaize, J., et al. 2017, A&A, 602, A5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Novak, M., Smolčić, V., Schinnerer, E., et al. 2018, A&A, 614, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Ocran, E. F., Taylor, A. R., Vaccari, M., et al. 2020, MNRAS, 491, 5911 [NASA ADS] [CrossRef] [Google Scholar]
Owen, F. N. 2018, ApJS, 235, 34 [Google Scholar]
Panessa, F., Baldi, R. D., Laor, A., et al. 2019, Nat. Astron., 3, 387 [Google Scholar]
Park, S. Q., Barmby, P., Fazio, G. G., et al. 2008, ApJ, 678, 744 [NASA ADS] [CrossRef] [Google Scholar]
Pasini, T., Brüggen, M., Hoang, D. N., et al. 2022, A&A, 661, A13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Peacock, J. A., & Nicholson, D. 1991, MNRAS, 253, 307 [NASA ADS] [CrossRef] [Google Scholar]
Pérez-González, P. G., Rieke, G. H., Egami, E., et al. 2005, ApJ, 630, 82 [Google Scholar]
Pérez-González, P. G., Rieke, G. H., Villar, V., et al. 2008, ApJ, 675, 234 [Google Scholar]
Pillepich, A., Springel, V., Nelson, D., et al. 2018, MNRAS, 473, 4077 [Google Scholar]
Popesso, P., Concas, A., Cresci, G., et al. 2023, MNRAS, 519, 1526 [Google Scholar]
Radcliffe, J. F., Barthel, P. D., Thomson, A. P., et al. 2021, A&A, 649, A27 [EDP Sciences] [Google Scholar]
Sabater, J., Best, P. N., Hardcastle, M. J., et al. 2019, A&A, 622, A17 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Sabater, J., Best, P. N., Tasse, C., et al. 2021, A&A, 648, A2 [EDP Sciences] [Google Scholar]
Sargent, M. T., Schinnerer, E., Murphy, E., et al. 2010a, ApJS, 186, 341 [Google Scholar]
Sargent, M. T., Schinnerer, E., Murphy, E., et al. 2010b, ApJ, 714, L190 [Google Scholar]
Saunders, W., Rowan-Robinson, M., Lawrence, A., et al. 1990, MNRAS, 242, 318 [Google Scholar]
Schmidt, M. 1968, ApJ, 151, 393 [Google Scholar]
Schreiber, C., Pannella, M., Elbaz, D., et al. 2015, A&A, 575, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Scoville, N., Aussel, H., Brusa, M., et al. 2007, ApJS, 172, 1 [Google Scholar]
Seymour, N., Dwelly, T., Moss, D., et al. 2008, MNRAS, 386, 1695 [NASA ADS] [CrossRef] [Google Scholar]
Smolčić, V., Schinnerer, E., Zamorani, G., et al. 2009, ApJ, 690, 610 [CrossRef] [Google Scholar]
Smolčić, V., Delvecchio, I., Zamorani, G., et al. 2017a, A&A, 602, A2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Smolčić, V., Novak, M., Bondi, M., et al. 2017b, A&A, 602, A1 [Google Scholar]
Stalevski, M., Fritz, J., Baes, M., Nakos, T., & Popović, L. Č. 2012, MNRAS, 420, 2756 [Google Scholar]
Stalevski, M., Ricci, C., Ueda, Y., et al. 2016, MNRAS, 458, 2288 [Google Scholar]
Tasse, C., Shimwell, T., Hardcastle, M. J., et al. 2021, A&A, 648, A1 [EDP Sciences] [Google Scholar]
Uchiyama, H., Yamashita, T., Nagao, T., et al. 2022, ApJ, 934, 68 [NASA ADS] [CrossRef] [Google Scholar]
Vaddi, S., O’Dea, C. P., Baum, S. A., et al. 2016, ApJ, 818, 182 [NASA ADS] [CrossRef] [Google Scholar]
van der Vlugt, D., Hodge, J. A., Algera, H. S. B., et al. 2022, ApJ, 941, 10 [NASA ADS] [CrossRef] [Google Scholar]
Wang, T., Elbaz, D., Alexander, D. M., et al. 2017, A&A, 601, A63 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Weaver, J. R., Kauffmann, O. B., Ilbert, O., et al. 2022, ApJS, 258, 11 [NASA ADS] [CrossRef] [Google Scholar]
Weinberger, R., Springel, V., Hernquist, L., et al. 2017, MNRAS, 465, 3291 [Google Scholar]
Williams, R. J., Quadri, R. F., Franx, M., van Dokkum, P., & Labbé, I. 2009, ApJ, 691, 1879 [NASA ADS] [CrossRef] [Google Scholar]
Williams, W. L., & Röttgering, H. J. A. 2015, MNRAS, 450, 1538 [NASA ADS] [CrossRef] [Google Scholar]
Worpel, H., Brown, M. J. I., Jones, D. H., Floyd, D. J. E., & Beutler, F. 2013, ApJ, 772, 64 [NASA ADS] [CrossRef] [Google Scholar]
Xie, F.-G., & Yuan, F. 2017, ApJ, 836, 104 [NASA ADS] [CrossRef] [Google Scholar]
Xue, Y. Q., Luo, B., Brandt, W. N., et al. 2016, ApJS, 224, 15 [Google Scholar]
Yang, G., Boquien, M., Buat, V., et al. 2020, MNRAS, 491, 740 [Google Scholar]
Yang, G., Boquien, M., Brandt, W. N., et al. 2022, ApJ, 927, 192 [NASA ADS] [CrossRef] [Google Scholar]
Yun, M. S., Reddy, N. A., & Condon, J. J. 2001, ApJ, 554, 803 [Google Scholar]
Zheng, X. C., Röttgering, H. J. A., van der Wel, A., & Duncan, K. 2022, A&A, 665, A114 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Zhu, S., Brandt, W. N., Zou, F., et al. 2023, MNRAS, 522, 3506 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Spectral energy distribution fitting

CIGALE uses a series of AGN and galaxies templates to efficiently model the observed multi-wavelength data from the X-rays and far-ultraviolet (FUV) to the far-infrared (FIR) and radio bands. Through a Bayesian-like approach, this code estimates many fundamental physical properties, such as star formation rate, stellar mass, dust luminosity, AGN contribution, and other quantities. In the fitting process the galaxy templates used in this work include the following six modules: star formation history, single stellar population (Bruzual & Charlot 2003), dust attenuation (Calzetti et al. 2000), dust emission (Dale et al. 2014), and AGN module (Stalevski et al. 2012, 2016). These modules and their parameters are summarized in Table A.1. In Fig. A.1 we present the best-fit SEDs of four radio sources in the GOODS-N field as examples.

Table A.1.

Model parameters for SED fitting with CIGALE

Fig. A.1.

Best-fit SED model from CIGALE fitting for four radio sources as examples. Panels (A), (B), (C), and (D) show radio sources from the GOODS-N field with redshift z of 0.504, 1.76, 2.33, and 3.48, respectively.

Appendix B: Correcting the classification purity for the radio-excess AGNs

We use the cross point between the highest and the second-highest Gaussian components of q_TIR distribution as the threshold to select radio-excess AGNs (see details in Section 4.1). Radio-excess AGNs with q_TIR around the selection threshold (which corresponds to the q_AGN in Fig. 3 and the vertical black dash-dotted lines in Fig. B.1) are partially contaminated by the SFGs. To correct the classification purity for radio-excess AGN, we follow Delvecchio et al. (2022) to define f_AGN as the 1 − N_SFG/N_TOT, where N_SFG represents the best-fit Gaussian model for the highest q_TIR peak (see the highest dotted curves in Fig. B.1), and N_TOT represents the best-fit model for the entire q_TIR distribution (see the purple solid curves in Fig. B.1). At the threshold q_AGN (see the black dash-dotted line in Fig. B.1), f_AGN is usually below 50% for samples in the COSMOS field, while f_AGN is around 100% for the GOODS-N field. It indicates that the Radio-excess AGN Sample in the GOODS-N have a nearly 100% AGN purity under the selection threshold used in this work. Therefore, we only make AGN purity corrections for the objects in the COSMOS fields. Following Delvecchio et al. (2022), we use the L_1.4GHz threshold converted from the median IR luminosity in each redshift bin to make the purity corrections.

Fig. B.1.

Distribution of IR-to-radio ratio (q_TIR) for SFGs in the Radio Sources Sample (upper row) and fractions of radio-excess AGNs (hosted by SFGs) in all the SFGs of the Radio Sources Sample (bottom row). Upper row: At each redshift bin (from left to right) the q_TIR distributions of the GOODS-N field and the COSMOS/UltraVISTA field are shown as the yellow histogram and the violet histogram, respectively. Each dotted curve represents each single-Gaussian model, while the solid curve represents the best-fit model to the entire q_TIR distribution. The vertical dashed line corresponds to the peak position of the highest Gaussian component. The vertical black dash-dotted line represents the threshold used in this work to separate SFGs and radio-excess AGNs (which corresponds to q_AGN in Fig. 3; see details in Section 4.1). Bottom row: At each redshift bin (from left to right) fractions of the radio-excess AGNs for the GOODS-N field and the COSMOS/UltraVISTA field are shown as the yellow dashed curve and the violet dotted curve, respectively. The parameter f_AGN is defined as 1 − N_SFG/N_TOT, where N_SFG represents the height of the highest Gaussian model (the highest dotted curve) at a certain q_TIR value, and N_TOT represents the height of the entire best-fit model (entire solid curve) at a certain q_TIR value. The vertical black dash-dotted line represents the threshold used in this work to separate SFGs and radio-excess AGNs (as in the upper panels).

Appendix C: Radio luminosity function in the GOODS-N and COSMOS/UltraVISTA fields

We obtain the radio luminosity function in each field following the same procedure described in Section 5. For both SFGs and radio-excess AGNs, RLFs in these three fields present generally consistent results (see detailed discussions in Section 5).

Fig. C.1.

1.4 GHz radio luminosity function of SFGs (top row) and radio-excess AGNs (bottom row) in each field (green circles: GOODS-N; blue triangles: COSMOS/UltraVISTA). The black solid curves in the top row represent the best-fit model (PLE) for SFGs. The black solid curves in the bottom row represent the best-fit model (PDE) for radio-excess AGNs.

Appendix D: Tables of radio luminosity functions for SFGs and radio-excess AGNs

Table D.1.

Radio luminosity function of SFGs obtained with the 1/V_max method.

Table D.2.

Radio luminosity function of radio-excess AGNs obtained with the 1/V_max method.

Appendix E: Radio luminosity function for radio-excess AGNs split into the subsets hosted by SFGs and QGs

We obtain the radio luminosity function for radio-excess AGNs hosted by different populations (SFGs and QGs) following the same procedure described in Section 5. We found a similar evolution trend for radio-excess AGNs hosted by different populations to that for low-excitation radio galaxies in Kondapally et al. (2022) (see detailed discussions in Section 5.3).

Fig. E.1.

1.4 GHz radio luminosity function of all the radio-excess AGNs (black squares) split into the subsets hosted by SFGs (blue triangles) and QGs (red circles). The grey solid curves represent the best-fit PDE model for all the radio-excess AGNs and the gray regions represent the 1σ uncertainty of the best-fit model.

Appendix F: Simple χ² fits to the probability of hosting a radio-excess AGN in each redshift bin

As Equations 22 and 23 show, at each fixed L_R (or M_⋆), we assume a simple power-law relation for p(L_R | M_⋆, z) as a function of M_⋆ (or L_R). Then we apply χ² fits to the data (see colored symbols in Fig. F.1 and Fig. F.2). The best-fit model in each L_R (or M_⋆) bin is shown as the colored line in the left column of Fig. F.1 (or Fig. F.2). The best-fit parameters in each L_R (or M_⋆) bin are shown as colored points in the right column of Fig. F.1 (or Fig. F.2). For more details about analysis and discussion, we refer to Section 6.3.3.

Fig. F.1.

Probability density function p(L_R | M_⋆, z) of a SFG hosting a radio-excess AGN with a given 1.4 GHz luminosity L_R as a function of stellar mass M_⋆ (left column) and the corresponding best-fit parameters (right column). Left column: At each redshift bin (from top to bottom) the different colors represent different L_R bins. Each dashed line represents the best-fit power-law relation for p(L_R | M_⋆, z) as a function of M_⋆ in each L_R bin (defined by Equation 22). Right column: Best-fit slope b (upper panel at each redshift bin) and intercept a (bottom panel at each redshift bin) of the above power-law relation. The different colors represent different L_R bins (same as those in the left column). In the upper right panel of each redshift bin, the gray dashed line shows the mean of slopes over all the L_R bins ( $\bar{b}$ $\bar{b}$ ), while the gray region represents $\bar{b} \pm 0.1$ $\bar{b} \pm 0.1$ . In the bottom right panel of each redshift bin, the pink dashed line shows the mean of intercepts over all the L_R bins ( $\bar{a}$ $\bar{a}$ ), while the pink region represents $\bar{a} \pm 0.1$ $\bar{a} \pm 0.1$ . The region scale (±0.1) is just a reference value to show the difference between each data point with the mean value.

Fig. F.2.

Probability density function p(L_R | M_⋆, z) of a SFG with a given stellar mass M_⋆ hosting a radio-excess AGN as a function of 1.4 GHz luminosity L_R (left column) and the corresponding best-fit parameters (right column). Left column: At each redshift bin (from top to bottom) the different colors represent different M_⋆ bins. Each dashed line represents the best-fit power-law relation for p(L_R | M_⋆, z) as a function of L_R in each M_⋆ bin (defined by Equation 23). Right column: Best-fit slope d (upper panel at each redshift bin) and intercept c (bottom panel at each redshift bin) of the above power-law relation. The different colors represent different M_⋆ bins (same as the left column). In the upper right panel of each redshift bin, the gray dashed line shows the mean of slopes over all the M_⋆ bins ( $\bar{d}$ $\bar{d}$ ), while the gray region represents $\bar{d} \pm 0.1$ $\bar{d} \pm 0.1$ . In the bottom right panel of each redshift bin, the pink dashed line shows the mean of intercepts over all the M_⋆ bins ( $\bar{c}$ $\bar{c}$ ), while the pink region represents $\bar{c} \pm 0.1$ $\bar{c} \pm 0.1$ . The region scale (±0.1) is just a reference value to show the difference between each data point with the mean value.

Appendix G: Probability of hosting radio-excess AGNs in QGs and all the galaxies

Similar to Section 6.3 (mainly for SFGs), here we show the probability of hosting a radio-excess AGN in QGs and all the galaxies based on maximum-likelihood fitting over the entire redshift range.

Fig. G.1.

Estimates of p(L_R | M_⋆, z) as a function of M_⋆ and L_R based on our maximum-likelihood fitting results for QGs. The colored dashed lines in the top and bottom rows represent our best-fit model from the unbinned maximum-likelihood fitting (see Equation 28 and Table 2) through combining data in the three fields (GOODS-N and COSMOS/UltraVISTA) over the entire redshift range (0.1 < z < 4) evaluated at the center of each redshift bin. The binned data points are scaled with the probability estimated by the model (see Equation 28) using the N_obs/N_mdl method of Aird et al. (2012) (see details in Section 6.3.3). In the top panel, the different colors represent different M_⋆ bins. In the bottom panel, the different colors indicate different L_R bins.

Fig. G.2.

Same as Fig. 8 and Fig. G.1, but for all galaxies.

Appendix H: Evolution of p(L_R | M_⋆, z) in SFGs considering the M_⋆-dependent IRRC

In order to examine whether the M_⋆-dependent IRRC will affect our results, we follow van der Vlugt et al. (2022) to select radio-excess AGNs with q_TIR deviating more than 3σ from the M_⋆-dependent IRRC of Delvecchio et al. (2021). To compare with the results shown in Section 6.3 (green region in Fig. H.1), we also use Equation 28 to recalculate the probability of SFGs hosting a radio-excess AGN (see the dark blue region in Fig. H.1). Figure H.1 shows that using the M_*-dependent IRRC (Delvecchio et al. 2021) does not alter the results in this work.

Fig. H.1.

Best-fit parameters of the maximum-likelihood fitting for p(L_R | M_⋆, z) in SFGs over the entire redshift range. The green regions correspond to the best-fit parameters based on the sample selected by the method in Section 4 (same as the blue region in Fig. 7). The dark blue regions correspond to the best-fit parameters based on the sample selected by the M_*-dependent IRRC (Delvecchio et al. 2021, see details in Section 7).

All Tables

Table 1.

Best-fit evolution parameters of the radio luminosity function for SFGs and radio-excess AGNs.

In the text

Table 2.

Best-fit parameters from maximum-likelihood fitting with redshift evolution.

In the text

Table A.1.

Model parameters for SED fitting with CIGALE

In the text

Table D.1.

Radio luminosity function of SFGs obtained with the 1/V_max method.

In the text

Table D.2.

Radio luminosity function of radio-excess AGNs obtained with the 1/V_max method.

In the text

All Figures

Fig. 1.

Flowchart of the matching processes and sample selections (see details in Sect. 2). ⁽¹⁾Details about broadband SED fitting can be found in Sect. 3.1 and Appendix A. ⁽²⁾Details about IR-radio correlation (IRRC) analysis can be found in Sect. 4. ⁽³⁾Details about the radio luminosity function can be found in Sect. 5. ⁽⁴⁾Details about the maximum-likelihood (ML) fitting method used for calculating the probability of hosting a radio-excess AGN can be found in Sect. 6. ^(*) L_TIR: rest-frame total IR luminosity between 8 and 1000 μm.

In the text

Fig. 2.

Distribution of rest-frame 1.4 GHz luminosity of the Radio Sources Sample (blue circles) and the Radio-excess AGN Sample (orange diamonds) with redshift for the GOODS-N (top) and COSMOS/UltraVISTA (bottom) fields. The 5σ detection limit for the COSMOS/UltraVISTA field with σ = 2.3 μJy beam⁻¹ at 3 GHz (Smolčić et al. 2017b) is shown as a gray dashed line in the top panel and as a black solid line in the bottom panel. The 5σ detection limit for the GOODS-N field with σ = 2.2 μJy at 1.4 GHz (Owen 2018) is shown as a black solid line in the top panel and as a gray dashed line in the bottom panel. All the detection limits are calculated by Eq. (1) assuming a fixed spectral index of −0.8 (see details in Sect. 3.2).

In the text

Fig. 3.

Distribution of IR-radio-ratio (q_TIR) in SFGs (left column) and QGs (right column) at different redshift bins. The q_TIR distributions were obtained by merging the Radio Sources Sample from the GOODS-N and COSMOS/UltraVISTA fields. Each yellow dotted curve represents a single-Gaussian model, while the black solid curve represents the best-fit model to the entire q_TIR distribution. The vertical black dashed line (corresponding to the ${\bar{q}}_{TIR}$ $\bar{q}_{\mathrm{TIR}}$ value in the figure) represents the peak position of the highest Gaussian component. The vertical red dash-dotted line (corresponding to the q_AGN value in the figure) represents the cross point between the highest Gaussian component and the second-highest Gaussian component, which is used as the threshold to separate SFGs and radio-excess AGNs in this work. For QGs at 2.0 < z ≤ 4.0, it is difficult to find a best-fit model due to the low source number, so here the threshold at 0.1 < z ≤ 0.5 was used, which does not have a significant effect on the results.

In the text

Fig. 4.

Distribution of q_TIR with redshift. The gray points show the Radio Sources Sample in the GOODS-N and COSMOS/UltraVISTA fields (see details in Sect. 2). The blue points represent the ${\bar{q}}_{TIR}$ $\bar{q}_{\mathrm{TIR}}$ values of SFGs, which are shown in Fig. 3 and in Sect. 4.1. The blue dashed line represents the best-fit model to the data of SFGs, while the shaded blue region is the associated 1σ uncertainty of the best-fit model. The green dash-dotted line represents the best-fit model from Magnelli et al. (2015) using data in the GOODS-N, GOODS-S, ECDFS, and COSMOS fields. Both the red dotted line (Delhaize et al. 2017) and the yellow dash-dotted line (Novak et al. 2017) represent the best-fit model obtained with the VLA-COSMOS 3 GHz project but using different source selection criteria. The dark red dash-dotted line represents the best-fit model from Enia et al. (2022) based on data from the GOODS-N field. The black star represents the q_TIR value for SFGs in the local universe from Bell (2003).

In the text

Fig. 5.

1.4 GHz radio luminosity function of SFGs (blue circles) and radio-excess AGNs (red squares) in different redshift bins. For SFGs, the blue solid line represents the best-fit PLE model, which is the best model to describe the data in this work. For radio-excess AGNs, the red solid line represents the best-fit PDE model, which is the best model to describe the data. The magenta vertical dash-dotted line represents the turnover luminosity between the radio luminosity function of SFGs and radio-excess AGNs (see details in Sect. 5.4). The light blue circles show the results for SFGs from the VLA-COSMOS 3 GHz project (Novak et al. 2017), while the gold squares show the results for radio-AGNs from the VLA-COSMOS 3 GHz project (Smolčić et al. 2017a). The redshift bins in this work have different definitions than in Novak et al. (2017) and Smolčić et al. (2017a).

In the text

Fig. 6.

Evolution of the best-fit parameters for the 1.4 GHz radio radio luminosity function with redshift. Shown are the evolution parameter (α_L) of the pure luminosity evolution (PLE) for SFGs (Panel A); turnover position (L_⋆) of the PLE for SFGs, calculated by Eq. (9) (Panel B); evolution parameter (α_D) of the pure density evolution (PDE) for radio-excess AGN (Panel C); turnover normalization (Φ_⋆) of the PDE for radio-excess AGN, calculated by Eq. (10) (Panel D); and turnover luminosity (L_D) above which the number density of radio-AGNs surpasses that of SFGs (L_D) (Panel E). The gray dash-dotted lines in Panels A and B represent the results for SFGs from Novak et al. (2017), while the gray regions represent the corresponding 1σ uncertainty. The gray dash-dotted lines in Panels C and D represent the results for radio-excess AGNs from Smolčić et al. (2017a), while the gray regions represent the corresponding 1σ uncertainty. The stars in Panels B and D represent the local values obtained by Novak et al. (2017) and Mauch & Sadler (2007), respectively. The magenta star in Panel E represents the local L_D obtained by the cross point between the local SFGs RLF from Novak et al. (2017) and the local AGN RLF from Mauch & Sadler (2007). The best-fit values for α_L, α_D, L_⋆, Φ_⋆, and L_D are summarized in Table 1.

In the text

Fig. 7.

Best-fit parameters of the maximum-likelihood fitting for p(L_R | M_⋆, z) of SFGs (blue circles and blue dashed lines), QGs (red squares and red dash-dotted lines), and all galaxies (black stars and black solid lines). The colored symbols represent the best-fit parameters obtained by Eqs. (24)–(27) in each redshift bin (see details in Sects. 6.3.1 and 6.3.2). The colored lines show the best-fit results over the entire redshift range obtained by Eqs. (28)–(29) and (25)–(27) (see details in Sect. 6.3.3). The dotted gray lines represent the results for the X-ray AGNs in Aird et al. (2012).

In the text

Fig. 8.

Estimates of p(L_R | M_⋆, z) as a function of M_⋆ and L_R based on our maximum-likelihood fitting results for SFGs. The colored dashed lines in the top and bottom rows represent our best-fit model from the unbinned maximum-likelihood fitting (see Eq. (28) and Table 2) through combining data in the GOODS-N and COSMOS/UltraVISTA fields over the entire redshift range (0.1 < z < 4) evaluated at the center of each redshift bin. The binned data points are scaled with the probability estimated by the model (see Eq. (28)) using the N_obs/N_mdl method of Aird et al. (2012) (see details in Sect. 6.3.3). In the top panel, the different colors represent different M_⋆ bins. In the bottom panel, the different colors show different L_R bins.

In the text

	Fig. A.1. Best-fit SED model from CIGALE fitting for four radio sources as examples. Panels (A), (B), (C), and (D) show radio sources from the GOODS-N field with redshift z of 0.504, 1.76, 2.33, and 3.48, respectively.
In the text

Fig. B.1.

Distribution of IR-to-radio ratio (q_TIR) for SFGs in the Radio Sources Sample (upper row) and fractions of radio-excess AGNs (hosted by SFGs) in all the SFGs of the Radio Sources Sample (bottom row). Upper row: At each redshift bin (from left to right) the q_TIR distributions of the GOODS-N field and the COSMOS/UltraVISTA field are shown as the yellow histogram and the violet histogram, respectively. Each dotted curve represents each single-Gaussian model, while the solid curve represents the best-fit model to the entire q_TIR distribution. The vertical dashed line corresponds to the peak position of the highest Gaussian component. The vertical black dash-dotted line represents the threshold used in this work to separate SFGs and radio-excess AGNs (which corresponds to q_AGN in Fig. 3; see details in Section 4.1). Bottom row: At each redshift bin (from left to right) fractions of the radio-excess AGNs for the GOODS-N field and the COSMOS/UltraVISTA field are shown as the yellow dashed curve and the violet dotted curve, respectively. The parameter f_AGN is defined as 1 − N_SFG/N_TOT, where N_SFG represents the height of the highest Gaussian model (the highest dotted curve) at a certain q_TIR value, and N_TOT represents the height of the entire best-fit model (entire solid curve) at a certain q_TIR value. The vertical black dash-dotted line represents the threshold used in this work to separate SFGs and radio-excess AGNs (as in the upper panels).

In the text

	Fig. C.1. 1.4 GHz radio luminosity function of SFGs (top row) and radio-excess AGNs (bottom row) in each field (green circles: GOODS-N; blue triangles: COSMOS/UltraVISTA). The black solid curves in the top row represent the best-fit model (PLE) for SFGs. The black solid curves in the bottom row represent the best-fit model (PDE) for radio-excess AGNs.
In the text

	Fig. E.1. 1.4 GHz radio luminosity function of all the radio-excess AGNs (black squares) split into the subsets hosted by SFGs (blue triangles) and QGs (red circles). The grey solid curves represent the best-fit PDE model for all the radio-excess AGNs and the gray regions represent the 1σ uncertainty of the best-fit model.
In the text

Fig. F.1.

Probability density function p(L_R | M_⋆, z) of a SFG hosting a radio-excess AGN with a given 1.4 GHz luminosity L_R as a function of stellar mass M_⋆ (left column) and the corresponding best-fit parameters (right column). Left column: At each redshift bin (from top to bottom) the different colors represent different L_R bins. Each dashed line represents the best-fit power-law relation for p(L_R | M_⋆, z) as a function of M_⋆ in each L_R bin (defined by Equation 22). Right column: Best-fit slope b (upper panel at each redshift bin) and intercept a (bottom panel at each redshift bin) of the above power-law relation. The different colors represent different L_R bins (same as those in the left column). In the upper right panel of each redshift bin, the gray dashed line shows the mean of slopes over all the L_R bins ( $\bar{b}$ $\bar{b}$ ), while the gray region represents $\bar{b} \pm 0.1$ $\bar{b} \pm 0.1$ . In the bottom right panel of each redshift bin, the pink dashed line shows the mean of intercepts over all the L_R bins ( $\bar{a}$ $\bar{a}$ ), while the pink region represents $\bar{a} \pm 0.1$ $\bar{a} \pm 0.1$ . The region scale (±0.1) is just a reference value to show the difference between each data point with the mean value.

In the text

Fig. F.2.

Probability density function p(L_R | M_⋆, z) of a SFG with a given stellar mass M_⋆ hosting a radio-excess AGN as a function of 1.4 GHz luminosity L_R (left column) and the corresponding best-fit parameters (right column). Left column: At each redshift bin (from top to bottom) the different colors represent different M_⋆ bins. Each dashed line represents the best-fit power-law relation for p(L_R | M_⋆, z) as a function of L_R in each M_⋆ bin (defined by Equation 23). Right column: Best-fit slope d (upper panel at each redshift bin) and intercept c (bottom panel at each redshift bin) of the above power-law relation. The different colors represent different M_⋆ bins (same as the left column). In the upper right panel of each redshift bin, the gray dashed line shows the mean of slopes over all the M_⋆ bins ( $\bar{d}$ $\bar{d}$ ), while the gray region represents $\bar{d} \pm 0.1$ $\bar{d} \pm 0.1$ . In the bottom right panel of each redshift bin, the pink dashed line shows the mean of intercepts over all the M_⋆ bins ( $\bar{c}$ $\bar{c}$ ), while the pink region represents $\bar{c} \pm 0.1$ $\bar{c} \pm 0.1$ . The region scale (±0.1) is just a reference value to show the difference between each data point with the mean value.

In the text

Fig. G.1.

Estimates of p(L_R | M_⋆, z) as a function of M_⋆ and L_R based on our maximum-likelihood fitting results for QGs. The colored dashed lines in the top and bottom rows represent our best-fit model from the unbinned maximum-likelihood fitting (see Equation 28 and Table 2) through combining data in the three fields (GOODS-N and COSMOS/UltraVISTA) over the entire redshift range (0.1 < z < 4) evaluated at the center of each redshift bin. The binned data points are scaled with the probability estimated by the model (see Equation 28) using the N_obs/N_mdl method of Aird et al. (2012) (see details in Section 6.3.3). In the top panel, the different colors represent different M_⋆ bins. In the bottom panel, the different colors indicate different L_R bins.

In the text

	Fig. G.2. Same as Fig. 8 and Fig. G.1, but for all galaxies.
In the text

Fig. H.1.

Best-fit parameters of the maximum-likelihood fitting for p(L_R | M_⋆, z) in SFGs over the entire redshift range. The green regions correspond to the best-fit parameters based on the sample selected by the method in Section 4 (same as the blue region in Fig. 7). The dark blue regions correspond to the best-fit parameters based on the sample selected by the M_*-dependent IRRC (Delvecchio et al. 2021, see details in Section 7).

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.

[1] Aird, J., Nandra, K., Laird, E. S., et al. 2010, MNRAS, 401, 2531 [Google Scholar]

[2] Aird, J., Coil, A. L., Moustakas, J., et al. 2012, ApJ, 746, 90 [CrossRef] [Google Scholar]

[3] Alberts, S., Rujopakarn, W., Rieke, G. H., Jagannathan, P., & Nyland, K. 2020, ApJ, 901, 168 [NASA ADS] [CrossRef] [Google Scholar]

[4] An, F., Vaccari, M., Smail, I., et al. 2021, MNRAS, 507, 2643 [NASA ADS] [CrossRef] [Google Scholar]

[5] Appleton, P. N., Fadda, D. T., Marleau, F. R., et al. 2004, ApJS, 154, 147 [Google Scholar]

[6] Arnouts, S., Moscardini, L., Vanzella, E., et al. 2002, MNRAS, 329, 355 [Google Scholar]

[7] Baldwin, J. A., Phillips, M. M., & Terlevich, R. 1981, PASP, 93, 5 [Google Scholar]

[8] Bariuan, L. G. C., Snios, B., Sobolewska, M., Siemiginowska, A., & Schwartz, D. A. 2022, MNRAS, 513, 4673 [NASA ADS] [CrossRef] [Google Scholar]

[9] Barro, G., Pérez-González, P. G., Cava, A., et al. 2019, ApJS, 243, 22 [NASA ADS] [CrossRef] [Google Scholar]

[10] Bell, E. F. 2003, ApJ, 586, 794 [Google Scholar]

[11] Best, P. N., Kauffmann, G., Heckman, T. M., et al. 2005, MNRAS, 362, 25 [Google Scholar]

[12] Best, P. N., Kondapally, R., Williams, W. L., et al. 2023, MNRAS, 523, 1729 [NASA ADS] [CrossRef] [Google Scholar]

[13] Bonchi, A., La Franca, F., Melini, G., Bongiorno, A., & Fiore, F. 2013, MNRAS, 429, 1970 [CrossRef] [Google Scholar]

[14] Boquien, M., Burgarella, D., Roehlly, Y., et al. 2019, A&A, 622, A103 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[15] Brammer, G. B., van Dokkum, P. G., & Coppi, P. 2008, ApJ, 686, 1503 [Google Scholar]

[16] Brown, M. J. I., Jannuzi, B. T., Floyd, D. J. E., & Mould, J. R. 2011, ApJ, 731, L41 [Google Scholar]

[17] Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000 [NASA ADS] [CrossRef] [Google Scholar]

[18] Burgarella, D., Buat, V., & Iglesias-Páramo, J. 2005, MNRAS, 360, 1413 [NASA ADS] [CrossRef] [Google Scholar]

[19] Calistro Rivera, G., Williams, W. L., Hardcastle, M. J., et al. 2017, MNRAS, 469, 3468 [Google Scholar]

[20] Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682 [NASA ADS] [CrossRef] [Google Scholar]

[21] Ceraj, L., Smolčić, V., Delvecchio, I., et al. 2018, A&A, 620, A192 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[22] Chabrier, G. 2003, ApJ, 586, L133 [NASA ADS] [CrossRef] [Google Scholar]

[23] Cochrane, R. K., Kondapally, R., Best, P. N., et al. 2023, MNRAS, 523, 6082 [NASA ADS] [CrossRef] [Google Scholar]

[24] Condon, J. J. 1992, ARA&A, 30, 575 [Google Scholar]

[25] Condon, J. J., & Dressel, L. L. 1978, ApJ, 221, 456 [NASA ADS] [CrossRef] [Google Scholar]

[26] Condon, J. J., Matthews, A. M., & Broderick, J. J. 2019, ApJ, 872, 148 [NASA ADS] [CrossRef] [Google Scholar]

[27] Dale, D. A., Helou, G., Magdis, G. E., et al. 2014, ApJ, 784, 83 [Google Scholar]

[28] Davé, R., Anglés-Alcázar, D., Narayanan, D., et al. 2019, MNRAS, 486, 2827 [Google Scholar]

[29] Delhaize, J., Smolčić, V., Delvecchio, I., et al. 2017, A&A, 602, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[30] Del Moro, A., Alexander, D. M., Mullaney, J. R., et al. 2013, A&A, 549, A59 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[31] Delvecchio, I., Smolčić, V., Zamorani, G., et al. 2017, A&A, 602, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[32] Delvecchio, I., Daddi, E., Sargent, M. T., et al. 2021, A&A, 647, A123 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[33] Delvecchio, I., Daddi, E., Sargent, M. T., et al. 2022, A&A, 668, A81 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[34] Dewdney, P. E., Hall, P. J., Schilizzi, R. T., & Lazio, T. J. L. W. 2009, IEEE Proc., 97, 1482 [Google Scholar]

[35] Donley, J. L., Rieke, G. H., Rigby, J. R., & Pérez-González, P. G. 2005, ApJ, 634, 169 [Google Scholar]

[36] Donley, J. L., Koekemoer, A. M., Brusa, M., et al. 2012, ApJ, 748, 142 [Google Scholar]

[37] Donoso, E., Best, P. N., & Kauffmann, G. 2009, MNRAS, 392, 617 [Google Scholar]

[38] Donoso, E., Li, C., Kauffmann, G., Best, P. N., & Heckman, T. M. 2010, MNRAS, 407, 1078 [NASA ADS] [CrossRef] [Google Scholar]

[39] Draine, B. T., & Li, A. 2007, ApJ, 657, 810 [CrossRef] [Google Scholar]

[40] Dubner, G., & Giacani, E. 2015, A&A Rev., 23, 3 [NASA ADS] [CrossRef] [Google Scholar]

[41] Dullo, B. T., Knapen, J. H., Beswick, R. J., et al. 2023, MNRAS, 522, 3412 [NASA ADS] [CrossRef] [Google Scholar]

[42] Enia, A., Talia, M., Pozzi, F., et al. 2022, ApJ, 927, 204 [NASA ADS] [CrossRef] [Google Scholar]

[43] Fabian, A. C. 2012, ARA&A, 50, 455 [Google Scholar]

[44] Fiore, F., Feruglio, C., Shankar, F., et al. 2017, A&A, 601, A143 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[45] Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]

[46] Franzen, T. M. O., Seymour, N., Sadler, E. M., et al. 2021, PASA, 38, e041a [NASA ADS] [CrossRef] [Google Scholar]

[47] Gim, H. B., Yun, M. S., Owen, F. N., et al. 2019, ApJ, 875, 80 [NASA ADS] [CrossRef] [Google Scholar]

[48] Grogin, N. A., Kocevski, D. D., Faber, S. M., et al. 2011, ApJS, 197, 35 [NASA ADS] [CrossRef] [Google Scholar]

[49] Gruppioni, C., Pozzi, F., Rodighiero, G., et al. 2013, MNRAS, 432, 23 [Google Scholar]

[50] Haggard, D., Green, P. J., Anderson, S. F., et al. 2010, ApJ, 723, 1447 [NASA ADS] [CrossRef] [Google Scholar]

[51] Hardcastle, M. J., & Croston, J. H. 2020, New A Rev., 88, 101539a [NASA ADS] [CrossRef] [Google Scholar]

[52] Helou, G., Soifer, B. T., & Rowan-Robinson, M. 1985, ApJ, 298, L7 [Google Scholar]

[53] Ho, L. C. 2008, ARA&A, 46, 475 [Google Scholar]

[54] Hughes, A. M., Beasley, A., & Carilli, C. 2015, IAU Gener. Assembly, 29, 2255106 [Google Scholar]

[55] Ilbert, O., Arnouts, S., McCracken, H. J., et al. 2006, A&A, 457, 841 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[56] Ivison, R. J., Magnelli, B., Ibar, E., et al. 2010, A&A, 518, L31 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[57] Janssen, R. M. J., Röttgering, H. J. A., Best, P. N., & Brinchmann, J. 2012, A&A, 541, A62 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[58] Jarvis, M., Taylor, R., Agudo, I., et al. 2016, in MeerKAT Science: On the Pathway to the SKA, 6 [Google Scholar]

[59] Jin, S., Daddi, E., Liu, D., et al. 2018, ApJ, 864, 56 [Google Scholar]

[60] Kauffmann, G., & Heckman, T. M. 2009, MNRAS, 397, 135 [NASA ADS] [CrossRef] [Google Scholar]

[61] Kauffmann, G., Heckman, T. M., Tremonti, C., et al. 2003, MNRAS, 346, 1055 [Google Scholar]

[62] Kewley, L. J., Maier, C., Yabe, K., et al. 2013, ApJ, 774, L10 [Google Scholar]

[63] Koekemoer, A. M., Faber, S. M., Ferguson, H. C., et al. 2011, ApJS, 197, 36 [NASA ADS] [CrossRef] [Google Scholar]

[64] Kolwa, S., Jarvis, M. J., McAlpine, K., & Heywood, I. 2019, MNRAS, 482, 5156 [CrossRef] [Google Scholar]

[65] Kondapally, R., Best, P. N., Cochrane, R. K., et al. 2022, MNRAS, 513, 3742 [NASA ADS] [CrossRef] [Google Scholar]

[66] Kondapally, R., Best, P. N., Raouf, M., et al. 2023, MNRAS, 523, 5292 [NASA ADS] [CrossRef] [Google Scholar]

[67] Kriek, M., van Dokkum, P. G., Labbé, I., et al. 2009, ApJ, 700, 221 [Google Scholar]

[68] Kukreti, P., Morganti, R., Tadhunter, C., & Santoro, F. 2023, A&A, 674, A198 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[69] Lacy, M., Storrie-Lombardi, L. J., Sajina, A., et al. 2004, ApJS, 154, 166 [Google Scholar]

[70] Liu, D., Daddi, E., Dickinson, M., et al. 2018, ApJ, 853, 172 [Google Scholar]

[71] Luo, B., Brandt, W. N., Xue, Y. Q., et al. 2017, ApJS, 228, 2 [Google Scholar]

[72] Magliocchetti, M. 2022, A&A Rev., 30, 6 [NASA ADS] [CrossRef] [Google Scholar]

[73] Magliocchetti, M., Maddox, S. J., Hawkins, E., et al. 2004, MNRAS, 350, 1485 [NASA ADS] [CrossRef] [Google Scholar]

[74] Magliocchetti, M., Popesso, P., Brusa, M., et al. 2017, MNRAS, 464, 3271 [NASA ADS] [CrossRef] [Google Scholar]

[75] Magnelli, B., Ivison, R. J., Lutz, D., et al. 2015, A&A, 573, A45 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[76] Malavasi, N., Bardelli, S., Ciliegi, P., et al. 2015, A&A, 576, A101 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[77] Matthews, A. M., Condon, J. J., Cotton, W. D., & Mauch, T. 2021, ApJ, 914, 126 [CrossRef] [Google Scholar]

[78] Matzeu, G. A., Brusa, M., Lanzuisi, G., et al. 2023, A&A, 670, A182 [Google Scholar]

[79] Mauch, T., & Sadler, E. M. 2007, MNRAS, 375, 931 [Google Scholar]

[80] McAlpine, K., Jarvis, M. J., & Bonfield, D. G. 2013, MNRAS, 436, 1084 [CrossRef] [Google Scholar]

[81] McAlpine, K., Prandoni, I., Jarvis, M., et al. 2015, in Advancing Astrophysics with the Square Kilometre Array (AASKA14), 83 [Google Scholar]

[82] Merloni, A., Heinz, S., & di Matteo, T. 2003, MNRAS, 345, 1057 [Google Scholar]

[83] Molnár, D. C., Sargent, M. T., Leslie, S., et al. 2021, MNRAS, 504, 118 [CrossRef] [Google Scholar]

[84] Morrison, G. E., Owen, F. N., Dickinson, M., Ivison, R. J., & Ibar, E. 2010, ApJS, 188, 178 [Google Scholar]

[85] Mullaney, J. R., Alexander, D. M., Goulding, A. D., & Hickox, R. C. 2011, MNRAS, 414, 1082 [Google Scholar]

[86] Ni, Q., Aird, J., Merloni, A., et al. 2023, MNRAS, 524, 4778 [CrossRef] [Google Scholar]

[87] Noll, S., Burgarella, D., Giovannoli, E., et al. 2009, A&A, 507, 1793 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[88] Norris, R. P., Afonso, J., Bacon, D., et al. 2013, PASA, 30, e020a [NASA ADS] [CrossRef] [Google Scholar]

[89] Novak, M., Smolčić, V., Delhaize, J., et al. 2017, A&A, 602, A5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[90] Novak, M., Smolčić, V., Schinnerer, E., et al. 2018, A&A, 614, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[91] Ocran, E. F., Taylor, A. R., Vaccari, M., et al. 2020, MNRAS, 491, 5911 [NASA ADS] [CrossRef] [Google Scholar]

[92] Owen, F. N. 2018, ApJS, 235, 34 [Google Scholar]

[93] Panessa, F., Baldi, R. D., Laor, A., et al. 2019, Nat. Astron., 3, 387 [Google Scholar]

[94] Park, S. Q., Barmby, P., Fazio, G. G., et al. 2008, ApJ, 678, 744 [NASA ADS] [CrossRef] [Google Scholar]

[95] Pasini, T., Brüggen, M., Hoang, D. N., et al. 2022, A&A, 661, A13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[96] Peacock, J. A., & Nicholson, D. 1991, MNRAS, 253, 307 [NASA ADS] [CrossRef] [Google Scholar]

[97] Pérez-González, P. G., Rieke, G. H., Egami, E., et al. 2005, ApJ, 630, 82 [Google Scholar]

[98] Pérez-González, P. G., Rieke, G. H., Villar, V., et al. 2008, ApJ, 675, 234 [Google Scholar]

[99] Pillepich, A., Springel, V., Nelson, D., et al. 2018, MNRAS, 473, 4077 [Google Scholar]

[100] Popesso, P., Concas, A., Cresci, G., et al. 2023, MNRAS, 519, 1526 [Google Scholar]

[101] Radcliffe, J. F., Barthel, P. D., Thomson, A. P., et al. 2021, A&A, 649, A27 [EDP Sciences] [Google Scholar]

[102] Sabater, J., Best, P. N., Hardcastle, M. J., et al. 2019, A&A, 622, A17 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[103] Sabater, J., Best, P. N., Tasse, C., et al. 2021, A&A, 648, A2 [EDP Sciences] [Google Scholar]

[104] Sargent, M. T., Schinnerer, E., Murphy, E., et al. 2010a, ApJS, 186, 341 [Google Scholar]

[105] Sargent, M. T., Schinnerer, E., Murphy, E., et al. 2010b, ApJ, 714, L190 [Google Scholar]

[106] Saunders, W., Rowan-Robinson, M., Lawrence, A., et al. 1990, MNRAS, 242, 318 [Google Scholar]

[107] Schmidt, M. 1968, ApJ, 151, 393 [Google Scholar]

[108] Schreiber, C., Pannella, M., Elbaz, D., et al. 2015, A&A, 575, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[109] Scoville, N., Aussel, H., Brusa, M., et al. 2007, ApJS, 172, 1 [Google Scholar]

[110] Seymour, N., Dwelly, T., Moss, D., et al. 2008, MNRAS, 386, 1695 [NASA ADS] [CrossRef] [Google Scholar]

[111] Smolčić, V., Schinnerer, E., Zamorani, G., et al. 2009, ApJ, 690, 610 [CrossRef] [Google Scholar]

[112] Smolčić, V., Delvecchio, I., Zamorani, G., et al. 2017a, A&A, 602, A2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[113] Smolčić, V., Novak, M., Bondi, M., et al. 2017b, A&A, 602, A1 [Google Scholar]

[114] Stalevski, M., Fritz, J., Baes, M., Nakos, T., & Popović, L. Č. 2012, MNRAS, 420, 2756 [Google Scholar]

[115] Stalevski, M., Ricci, C., Ueda, Y., et al. 2016, MNRAS, 458, 2288 [Google Scholar]

[116] Tasse, C., Shimwell, T., Hardcastle, M. J., et al. 2021, A&A, 648, A1 [EDP Sciences] [Google Scholar]

[117] Uchiyama, H., Yamashita, T., Nagao, T., et al. 2022, ApJ, 934, 68 [NASA ADS] [CrossRef] [Google Scholar]

[118] Vaddi, S., O’Dea, C. P., Baum, S. A., et al. 2016, ApJ, 818, 182 [NASA ADS] [CrossRef] [Google Scholar]

[119] van der Vlugt, D., Hodge, J. A., Algera, H. S. B., et al. 2022, ApJ, 941, 10 [NASA ADS] [CrossRef] [Google Scholar]

[120] Wang, T., Elbaz, D., Alexander, D. M., et al. 2017, A&A, 601, A63 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[121] Weaver, J. R., Kauffmann, O. B., Ilbert, O., et al. 2022, ApJS, 258, 11 [NASA ADS] [CrossRef] [Google Scholar]

[122] Weinberger, R., Springel, V., Hernquist, L., et al. 2017, MNRAS, 465, 3291 [Google Scholar]

[123] Williams, R. J., Quadri, R. F., Franx, M., van Dokkum, P., & Labbé, I. 2009, ApJ, 691, 1879 [NASA ADS] [CrossRef] [Google Scholar]

[124] Williams, W. L., & Röttgering, H. J. A. 2015, MNRAS, 450, 1538 [NASA ADS] [CrossRef] [Google Scholar]

[125] Worpel, H., Brown, M. J. I., Jones, D. H., Floyd, D. J. E., & Beutler, F. 2013, ApJ, 772, 64 [NASA ADS] [CrossRef] [Google Scholar]

[126] Xie, F.-G., & Yuan, F. 2017, ApJ, 836, 104 [NASA ADS] [CrossRef] [Google Scholar]

[127] Xue, Y. Q., Luo, B., Brandt, W. N., et al. 2016, ApJS, 224, 15 [Google Scholar]

[128] Yang, G., Boquien, M., Buat, V., et al. 2020, MNRAS, 491, 740 [Google Scholar]

[129] Yang, G., Boquien, M., Brandt, W. N., et al. 2022, ApJ, 927, 192 [NASA ADS] [CrossRef] [Google Scholar]

[130] Yun, M. S., Reddy, N. A., & Condon, J. J. 2001, ApJ, 554, 803 [Google Scholar]

[131] Zheng, X. C., Röttgering, H. J. A., van der Wel, A., & Duncan, K. 2022, A&A, 665, A114 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

[132] Zhu, S., Brandt, W. N., Zou, F., et al. 2023, MNRAS, 522, 3506 [NASA ADS] [CrossRef] [Google Scholar]

Cosmic evolution of radio-excess active galactic nuclei in quiescent and star-forming galaxies across 0 < z < 4

1. Introduction

2. Data and samples

2.1. GOODS-N field

2.2. COSMOS/UltraVISTA field

3. Galaxy property estimation

3.1. Rest-frame 8–1000 μm infrared luminosity

3.2. Rest-frame 1.4 GHz radio luminosity

3.3. Stellar mass and UVJ magnitude

3.4. Galaxy type: Star-forming and quiescent galaxies

4. Radio-excess active galactic galaxies

4.1. qTIR distribution and radio-excess AGN selection

4.2. qTIR evolution

5. Radio luminosity functions and their cosmic evolution

5.1. Constructing radio luminosity functions out to z ∼ 4

5.2. Fitting the evolution of a star-forming galaxy’s radio luminosity function

5.3. Fitting the evolution of an AGN’s radio luminosity function

5.4. Crossover luminosity between SFG RLF and AGN RLF

6. The probability of hosting radio-excess AGNs in SFGs and QGs

6.1. Calculating the probability for observational data

6.2. Simple χ2 fits in each redshift bin

6.3. Redshift evolution

6.3.1. Maximum-likelihood fitting

6.3.2. Results in each redshift bin

6.3.3. Incorporating redshift evolution into the maximum-likelihood fitting

6.3.4. Quantitative relation of radio-excess AGN fraction in SFGs and QGs

7. Discussion: Comparison with previous works

8. Summary

Acknowledgments

References

Appendix A: Spectral energy distribution fitting

Appendix B: Correcting the classification purity for the radio-excess AGNs

Appendix C: Radio luminosity function in the GOODS-N and COSMOS/UltraVISTA fields

Appendix D: Tables of radio luminosity functions for SFGs and radio-excess AGNs

Appendix E: Radio luminosity function for radio-excess AGNs split into the subsets hosted by SFGs and QGs

Appendix F: Simple χ2 fits to the probability of hosting a radio-excess AGN in each redshift bin

Appendix G: Probability of hosting radio-excess AGNs in QGs and all the galaxies

Appendix H: Evolution of p(LR | M⋆, z) in SFGs considering the M⋆-dependent IRRC

All Tables

All Figures

4.1. q_TIR distribution and radio-excess AGN selection

4.2. q_TIR evolution

6.2. Simple χ² fits in each redshift bin

Appendix F: Simple χ² fits to the probability of hosting a radio-excess AGN in each redshift bin

Appendix H: Evolution of p(L_R | M_⋆, z) in SFGs considering the M_⋆-dependent IRRC