A striking relationship between dust extinction and radio detection in DESI QSOs: evidence for a dusty blow-out phase in red QSOs

Introduction

Quasi-stellar objects (QSOs), also known as quasars, are the most powerful class of Active Galactic Nuclei (AGN), with extremely high bolometric luminosities ($10^{45}$–$10^{48}$ erg s$^{-1}$). These high luminosities are now known to be due to mass accretion onto a supermassive black-hole (SMBH; masses of $10^7$–$10^9$ M$_{\odot}$). The observed correlation between the masses of SMBHs at the centers of spheroidal galaxies and the bulge mass of the hosting galaxy suggests an intrinsic link between the growth of the SMBH and the surrounding host galaxy [68].

A canonical model for this growth phase invokes gas-rich major mergers [25, 3, 23] in which rapid gas inflow to the nucleus fuels a powerful QSO, obscured in its early stages by a high column of gas and dust. Through powerful winds and/or outflows, an unobscured “typical” blue QSO is eventually revealed. However, it is still unclear what type of object represents the transitional “blow-out” stage between the obscured and unobscured QSO phase. Potential candidates include; 1) the well-studied dust obscured “red QSOs” (e.g., [86, 85, 6, 128, 11]), selected as broad line QSOs with red optical or near-infrared (NIR) colours; 2) “Broad Absorption Line Quasars” (BALQSOs; e.g., [72, 46, 73]), which display broad absorption typically blueward of the C iv emission line; and 3) the rarer “Extremely Red Quasars” (ERQs; e.g., [1, 84]), selected to have red optical–mid-infrared (MIR) colours. All of these objects are thought to play an important role in the regulation of star formation in the host galaxy (otherwise known as “AGN feedback”; see reviews by [48] and [108]).

Over recent years, a number of groups have explored whether red QSOs represent a transitional blow-out phase, determining their properties and how they relate to typical blue QSOs [35, 85, 6, 5, 53, 14, 18, 10, 4, 26, 11, 30, 19, 29, 22, 51, 31, 76, 90]. In our recent work, we have shown that the majority of red QSOs from the Sloan Digital Sky Survey (SDSS; [52]), selected based on their optical $g-i$ colour, are red due to line-of-sight dust extinction ([31, 51]). Somewhat surprisingly we found that red QSOs exhibited enhanced radio emission, compared to blue QSOs, which is typically compact ($<$ $2$ kpc; [11, 29, 30, 19]) and peaks around the radio-quiet/radio-loud threshold [29, 30]. The observed differences in the radio properties of red and blue QSOs do not appear to be associated with any differences in the star formation ([29, 31, 94]), accretion [11, 51], or large scale structure [76] properties between the two populations. On the basis of these results, a potential self-consistent scenario that links the enhanced radio emission to the dust in red QSOs is that the radio emission is due to shocks produced by either an outflow or a jet interacting with a higher opacity interstellar medium (ISM)/circumnuclear environment. These outflows/jets could also potentially drive out the surrounding gas and dust (i.e., AGN feedback), eventually revealing a typical blue QSO; therefore, red QSOs are prime candidates for an intermediate stage in the evolution of QSOs.

It is worth noting that the measured dust extinctions towards optically selected SDSS red QSOs are modest ($A_V$ $<$ $0.7$ mag; [31, 51]) and they likely represent the tip of the iceberg of a larger, more heavily reddened QSO population. Indeed, other studies have utilized NIR photometry in order to select more heavily reddened red QSO samples (e.g., [4, 85, 6, 53]), with measured dust extinctions up to $A_V$ $\sim$ $5$ mag. For example, [53] used a WISE–2MASS colour selection to define a sample of red QSOs with dust extinctions $E(B-V)>0.25$ mag ($A_V\gtrsim0.8$ mag)¹. Utilizing radio data from both the Faint Images of the Radio Sky at Twenty-centimeters (FIRST; [34]) and the Very Large Array Sky Survey (VLASS; [20]), they found a higher radio detection fraction for their red QSOs compared to a blue QSO sample ($E(B-V)<0.25$ mag), in addition to a higher fraction of compact radio morphologies for the red QSOs in both radio datasets. These results are remarkably similar to our less reddened optically selected red QSOs [11, 29, 30].

A limitation of our previous red QSO studies is the shallow magnitude limit of SDSS which restricts us to the most luminous and optically brightest objects, with only modest amounts of extinction. We therefore require larger, fainter, and redder samples of QSOs to robustly test whether red QSOs represent a transitional phase in the evolution of QSOs: the Dark Energy Spectroscopic Instrument (DESI; [44, 57]) has the capability to do just this. DESI has completed the first year of a five year survey which aims to obtain optical spectra for $\sim$ $3$ million QSOs down to $r$ $\sim$ $23.2$ mag [55, 56, 123]. These data will be revolutionary in QSO science, pushing 1–2 mag fainter than SDSS and, consequently, identifying both more obscured and more typical systems (see Section 5 of [61] for diversity in the QSO spectra).

The DESI QSO target selection consists of both an optical colour selection and a Random Forest selection [77, 62]². Like any standard optical QSO selection, both the Random Forest and colour selection are biased towards optically blue QSOs since the majority of QSOs, and therefore the Random Forest training sets, are blue. This bias will inevitably select against the more “exotic” and dust reddened QSO population. In October 2020, a call for DESI proposals for “secondary target” programs was announced, with the aim of utilizing spare DESI fibres and expanding the DESI science beyond the key cosmological aims of the survey. Our group successfully proposed a program (PI: V. Fawcett) which targets dust-reddened QSOs that have been missed by the nominal DESI QSO selections, with aims to: 1) push to higher levels of dust extinction than can be identified by either the SDSS or the nominal DESI QSO selection, 2) determine whether the enhanced radio emission found in SDSS red QSOs extends to more heavily reddened systems, 3) explore the origin of the enhanced radio emission found in red QSOs, and 4) explore exotic sub-populations of QSOs within DESI.

In this first paper we introduce and present the first eight months of data from our DESI secondary target program and explore the connection between dust extinction and radio emission in DESI QSOs, using data from the LOFAR Two-metre Sky Survey (LoTSS) Data Release 2 ([54]). In § 2 we describe our secondary target program and in § 3 we describe how our secondary target program extends the nominal DESI QSO selection, the methods used in this paper, and the construction of our combined QSO sample which is used for the scientific analyses in this paper. In § 4 we present our results and in § 5 we discuss how our sample extends previous red QSO work and explore the origin of the radio emission. Unless specified otherwise, in this paper we use AB magnitudes. Throughout our work we adopt a flat \uplambda-cosmology with $H_0$ $=$ 70 km s$^{-1}$Mpc$^{-1}$, $\Omega$M $=$ 0.3 and \Omega_{\uplambda} $=$ 0.7. The data underlying the figures in this paper are available online: https://doi.org/10.5281/zenodo.8147342.

Secondary target program

In this section we describe our secondary target selection and the DESI data (§ 2.1) and spectral classification method (§ 2.2).

Defining our secondary target program selection

Our secondary target program aims to expand the DESI QSO colour space, utilizing MIR data from the Wide-field Infrared Survey Explorer (WISE; [12]). The main goal of this program is to build up a statistically significant sample of dust-reddened QSOs that may have otherwise been missed by the nominal DESI QSO selections. However, due to our optical colour selection (which we describe in more detail below) our program will also select extremely blue QSOs. Therefore it should be noted that, although the bitmask³ assigned to the objects observed by our program is QSO_RED, the sample will include a small fraction of optically blue QSOs, in addition to dust-reddened QSOs. For the analyses in this paper, our secondary target program is solely used to boost the source statistics of the rarer dust-reddened QSOs and in combination with a larger nominal sample of DESI QSOs, the majority of which will be typical blue QSOs by selection.

Our current sample (which we refer to as the “DESI secondary program (SP) sample”) is taken from the first $\sim$ $8$ months of observations, including all of the DESI Survey Validation (SV; [124]) data (Dec 2020–June 2021) and part of the ‘Main’ survey data (May 2021–July 2021).⁴ The DESI SV observations are split into three different targeting campaigns: ‘SV1’, ‘SV2’, and ‘SV3’; our DESI SP sample only utilizes SV1 and SV3 observations since no secondary targets were assigned during SV2 [65, 124]. During the SV3 phase, primary targets were prioritized for multiple observations before secondary targets were scheduled in order to obtain highly complete samples of DESI primary targets. Thus, a portion of our SP sample is biased towards the colours of primary DESI classes, particularly Luminous Red Galaxies (LRGs; [79]; $z$ $<$ $21.7$ limit in SV3) and QSOs ([62]; $r$ $<$ $23$ limit in SV3). This bias will gradually diminish over the course of the DESI Main survey, for which our SP targets are assigned an observational priority at least as high as any other primary target (see Section 5 of [129]). The flowchart in Figure 1 displays our secondary target selection process, which we describe below.

A key distinguishing component of our SP compared to the DESI QSO sample is the MIR WISE selection. We first started with the DESI DR9 legacy photometric catalogues [67]⁵, which contain optical photometry in the $g$- (4730 Å), $r$- (6420 Å), and $z$-bands (9260 Å), as well as four MIR bands (at 3.4, 4.6, 12, and 22 \upmum) from stacks of WISE/NEOWISE data [12, 87, 88], referred to as unWISE coadds ([130]; for more details see Section 5 of [67]). We first used the MASKBIT column in order to remove objects with flagged photometric errors. In our selection we remove objects with a MASKBIT corresponding to either BRIGHT (object touches the pixel of a bright star), ALLMASK_G, ALLMASK_R, ALLMASK_Z (object touches a bad pixel in either the $g$-, $r$-, or $z$-band images), BAILOUT (object touches a pixel in a region where source fitting failed), GALAXY (object touches a pixel in a large galaxy), or CLUSTER (object touches a pixel in a globular cluster).⁶ These were chosen to be consistent with the nominal DESI QSO SV target selection [62]. We then applied a morphology cut, requiring either point sources (“PSF” morphology) or a small relative $\chi^2$ difference between the PSF and more extended morphological models (defined in [62]) from the DR9 legacy imaging [67, 63].

To ensure we minimize non-QSO contaminants in our sample, we required a signal-to-noise (SNR) cut of $>$ $3$ in the WISE $W1$, $W2$, and $W3$ bands and applied the [16] “AGN wedge” (a $W1-W2$ and $W2-W3$ colour-colour selection). We chose to use the [16] wedge rather than more simple, liberal MIR cuts (e.g., [13]; $W1-W2\geq0.8$, Vega) in order to optimize the number of reliable QSOs observed by our SP, in addition to obtaining accurate $W3$ fluxes used to calculate the $6$ \upmum luminosity (see § 3.4). However, this may result in fewer sources observed which have MIR colours that lie outside of the AGN wedge (in particular, exotic QSOs with more extreme MIR colours or higher redshift QSOs for which the AGN wedge is known to be less reliable). After applying the MIR cuts, we then imposed an $r$-band magnitude cut of $r$ $<$ $23.5$ mag (matching the more liberal cut used for the SV QSO selection). We further removed any object that falls within the region of $g-r$ versus $r-z$ colour-colour space that is the predominant colour selection for nominal DESI QSOs ($g-r$ $<$ $1.3$ & $r-z$ $>$ $-0.4$ & $r-z$ $<$ $1.1$; see Figure 2), in order to expand the DESI QSO colour space. This resulted in predominantly red objects, although, as noted earlier, the selection also includes extremely blue objects (e.g., $\sim$ 12.3 per cent of our secondary program QSOs have an $r-z<-0.4$; see Figure 2). The nominal QSO selection includes a Random Forest selection that will also target QSOs outside this colour region, although these objects will be rare; see Figure 5 and [77]. The number of objects that both satisfy these criteria, and were specifically targeted by our program, is 3038.⁷ It is worth noting that if this program continues for the full five years of the DESI survey, we expect over $\sim$ 8 times this number of SP targets to be observed.

Spectral classification

The DESI spectrograph covers an observed wavelength range of $\sim$ 3600–9800 Å, with a spectral resolution of $\sim$ 2000–5000 (blue to red). In order to measure the redshift and determine the spectral type (hereafter “spectype”) of each spectrum, DESI makes use of three tools. redrock⁸ ([?]) is the standard DESI spectral template-fitting code that uses a set of templates to classify objects into one of three broad categories: QSO, GALAXY, and STAR. redrock uses a Bayesian approach and selects the best-fitting solutions based on the lowest reduced $\chi^{2}$ values. In addition to redrock, DESI also utilizes two “afterburner” codes that are applied after redrock is run: QuasarNET⁹ [60], a neural network-based QSO classifier that updates the spectype and redshift of missed high confidence QSOs, and the Mg ii afterburner, which is an emission-based code that updates the spectype of a classified GALAXY to QSO if there is significant broad Mg ii emission present in the spectrum (see Section 6.2 in [62]). We refer to the combination of these three tools as the “modified pipeline”.

On the basis of visual inspection (VI), the modified pipeline achieves a good redshift ($\Delta v$ $<$ $3000$ km s$^{-1}$)¹⁰ and spectral type purity of 99 per cent and 94 per cent, respectively for the main QSO survey (see Tables 8 and 9 in [61]). However, from a VI categorization of QSOs with incorrect redshifts and/or incorrect spectypes from the modified pipeline, [61] showed that many were found to either have: 1) a very red continuum, 2) strong host-galaxy features, or 3) strong broad absorption features.¹¹ Therefore, our DESI SP sample, the majority of which have very red optical colors (Figure 2) in addition to pushing to fainter magnitudes than the main DESI QSO selection ($r$ $<$ $23$ mag), is likely to result in the sort of objects that the DESI pipeline struggles to classify. Due to this, we undertook complete VI of our SP sample of 3038 objects in order to retrieve accurate redshifts and spectypes. Our VI approach is based on that described in Section 2.2 of [61] and is summarized below.

Each of the 3038 objects in the SP sample are evaluated using the Prospect¹² visualization tool; for a detailed description of the Prospect tool, see [66]. Prospect displays the un-smoothed spectrum, along with the noise spectrum and an indication of where the three arms of the DESI spectrograph overlap. Available tools on the page allow the user to smooth the spectrum and display common absorption and emission lines. The nine best-fitting redrock template solutions are also displayed and the best four solutions are available to overlay on the spectrum. The visual assessor is then able to manually: 1) adjust and correct the redshift, 2) correct the spectype, 3) add a specific comment, 4) flag any bad data, and 5) assess the quality of the spectrum. Examples of specific comments that were added include flagging Type 2 AGN, BALQSOs, and spectra that included two sources (i.e., those with two redshift solutions). A quality class (VI_quality) of 0 to 4 is provided by the visual assessor for each spectrum, where 0 is the worst quality and 4 is the best quality (the specific meaning of each VI_quality is detailed in [61]). Every spectrum is visually assessed by two visual assessors and then the files are merged by a VI lead (D. M. Alexander for this dataset), who reassess spectra with significant disagreements between the two assessors and then decides on the final classification and/or redshift.¹³ The resulting spectrum is considered “high quality” with a robust redshift and spectype if the average VI quality class is $\geq2.5$ (i.e., the lowest quality spectra that are included in our high-quality sample are those where one assessor assigns a quality flag of 2 and the other assigns a 3).

Spectral type	Percentage of high-quality sources
QSO (Type 1)	69.7 per cent (1852/2658)
QSO (Type 2)	11.4 per cent (302/2658)
GALAXY	13.0 per cent (345/2658)
STAR	6.0 per cent (159/2658)

From visually inspecting our DESI SP sample, we found that $\sim$ $89$ per cent (2700/3038) of the DESI spectra have high quality (VI quality $\geq$ $2.5$) redshifts and spectypes. Out of these high-quality objects we found that $\sim$ 1.6 per cent (42/2700) had two redshift solutions; inspecting the images confirmed that there were indeed two objects within the same fibre, and we therefore removed these objects from our analyses. Excluding the low-quality sources and systems with two objects in the same spectrum, we found that $\sim$ $70$ per cent (1852/2658) of our SP sample are confirmed as Type 1 QSOs, displaying clear broad emission lines in their optical spectra; this is consistent with the expected fraction of high-quality QSOs from the main survey QSO targets ($\sim$ $76$ per cent; see Figure 4 in [61]). Additionally, $\sim$ 11 per cent (209/1852) of the QSOs were flagged as BALQSOs, displaying broad absorption blueward of C iv (e.g., [73]). We also found that $\sim$ $11$ per cent (302/2658) of the high-quality spectra are Type 2 AGN, displaying strong [Ne v]\uplambda\uplambda3346,3426 but no broad lines, $\sim$ $13$ per cent (345/2658) are galaxies (Figure 15), and $\sim$ $6$ per cent (159/2658) are stars (Figure 16), consistent with the main DESI QSO survey [61]; see Table 1 and Figure 3. This demonstrates that the [16] AGN wedge is efficient at selecting QSOs and that our program is indeed targeting additional QSOs beyond the main DESI survey, rather than introducing a large number of non-QSO contaminants.

Figure 2 displays the optical colour-colour distribution for the different spectral classes within our secondary sample.¹⁴ We find that $\sim$ 87, $\sim$ 84, and $\sim$ 81 per cent of the QSOs, galaxies, and stars have red optical colours ($g-r$ $>$ $1.3$ mag or $r-z$ $>$ $1.1$ mag), respectively. The number of sources for each spectype identified from our VI is displayed in Table 1. For a detailed comparison of the performance of the modified pipeline compared to the VI results, see Appendix § 9. In this paper we focus on the QSOs in our sample (Figure 4); for more discussion on the galaxies and stars within our DESI SP sample, see Appendix § 10.

Data and methods

In this section, we compare the properties of the QSOs selected by our secondary program to those selected by the nominal DESI QSO programs (§ 3.1), our spectral stacking method (§ 3.2), our dust fitting method (§ 3.3), the radio data utilized in this paper (§ 3.4), and the construction of the combined QSO sample used for the main analyses in this paper (§ 3.5).

How does our SP extend the nominal DESI QSO selection?

For the analyses in this paper we focus on the broad line QSOs in our SP sample (excluding the 306 Type 2 AGN; see Table 1). We also require a comparison QSO sample in order to assess how our SP sample extends the nominal DESI QSO population. We refer to the QSOs from the nominal DESI QSO population as “nominal QSOs” and the QSOs from our SP as “SP QSOs”. Our selection process for the nominal and SP QSO samples is displayed in Figure 4, which we describe below.

To construct the nominal QSO sample, we first selected objects that were targeted by the nominal DESI QSO program (bitmask $=$ QSO). We then applied the same MIR cuts that were used to form our SP sample (§ 2.1; SNR $>3$ in the $W123$ bands and satisfying the [16] AGN wedge; see Figure 5). Finally, we restricted the selection to objects that were classified as QSOs by the modified pipeline and apply ZWARN $==0$ to ensure the redshifts are robust. The resulting nominal QSO sample consists of 88 101 QSOs. To construct the SP QSO sample, we selected the visually confirmed QSOs from our SP sample (§ 2.2) which had a VI_quality $\geq2.5$ (i.e., had a robust spectype and redshift) and additionally removed the Type 2 AGN (no broad lines present) and systems with two objects in the same spectrum; this resulted in 1852 QSOs in our SP QSO sample

Due to re-observations, many DESI targets will obtain multiple spectra. In these cases, for our analyses we use spectra that coadd all observations of a given target. In addition to this, since the DESI SV observations each have distinct target selections (for more details see [65]) there may be multiple coadded spectra across the different campaigns. In these cases, we use the coadded spectrum with the highest squared template signal-to-noise ratio (TSNR2) value. The calculation of the TSNR2 value is target-class dependent and so in this work we use the TSNR2_QSO statistic. For more details on the TSNR2 metric, see [82].

Figure 5 displays the optical colour-colour distribution for both the nominal and SP QSO samples. We find that, as expected, the majority of the nominal QSO sample lie inside the main colour space used to select DESI QSOs, with only $\sim$ 3 per cent (2367/88 101) of QSOs that lie outside; these QSOs were likely selected by the Random Forest QSO selection. Our SP QSO sample expands the colour parameter space of DESI QSOs, with the majority of additional QSOs displaying very red optical colours. It is worth noting that $\sim$ 40 per cent (745/1852) of the SP QSOs were also targeted by the nominal QSO program; these 745 QSOs met both our DESI SP selection criteria and one of the nominal QSO selections. Our DESI SP therefore provides $\sim$ 32 per cent more QSOs that lie outside the main region of colour space, i.e., 1107 SP QSOs solely targeted by our program compared to the 2367 nominal QSOs that fall outside the main colour region (see Figure 5), despite being only $\sim$ 2 per cent the size of the nominal QSO sample.

To understand how the 1852 DESI SP QSOs compare to the SDSS QSOs utilized in our previous studies, we first restricted the SP QSO sample to the SDSS footprint, reducing the sample to 1312 QSOs. We then matched to the DR16 QSO catalogue [50] using a $1''$ search radius; this resulted in 107 matches. Therefore, only 8 per cent (107/1312) of the SP QSO sample within the SDSS footprint have also been observed by SDSS and met the criteria used by [50] to be included in the SDSS DR16 QSO catalog. Out of the 107 sources, we found seven with a redshift disagreement between our VI and the SDSS DR16 QSO catalogue; additional VI of these objects clarified that the SDSS redshifts were significantly incorrect. Matching the DR16 QSO catalogue with the other high-quality sources from our SP sample, we found three matches; one was classified as a galaxy and two were classified as Type 2 AGN by our VI. We further investigate these disagreements in Appendix § 11.

Figure 6 displays the L_{\rm 6\,\upmu m} and $r$-band magnitude versus redshift distributions for the nominal and SP QSO samples. To compute the L6 \upmum we used a log-linear interpolation or extrapolation of the fluxes in the W2 and W3 bands, following the approach from [11]. The larger density of SP QSOs at $r$ $\lesssim$ $22.5$ mag is a result of the different target prioritization during SV3 (see § 2.1); SP sources that were also primary targets i.e., QSOs [62], LRGs [79], or Emission Line Galaxies (ELGs; [80]), will have had a far higher probability of being observed than SP targets that were not also primary targets. The selection of these primary targets impose brighter magnitude limits than that of our SP, and therefore the targets observed during SV3 are more likely to reside at these brighter magnitudes. The SP QSOs that were also observed by SDSS are shown by the blue diamonds. The highest redshift QSO in the SP and nominal samples are at $z$ $=$ $4.4$ and 4.6, respectively. On average, the SP QSOs are $\sim$ $1.7$ magnitudes optically fainter than the QSOs in the nominal sample, demonstrating that our program is pushing to QSOs with fainter magnitudes. Applying the two-sided Kolmogorov-Smirnov (K-S) test to the redshift and L_{\rm 6\,\upmu m} distributions of the nominal and SP QSO samples, we cannot rule out that the two samples are drawn from the same parent distribution at a $<1$ per cent significance level. This also demonstrates that QSOs with fainter optical magnitudes are not necessarily intrinsically fainter; their optical emission may be relatively fainter due to dust extinction (see § 4.1). We utilize L_{\rm 6\,\upmu m} as a tracer for the intrinsic QSO luminosity in § 4.2.

Composite construction

To characterize the basic spectral properties of our QSO samples we constructed composite spectra following the approach outlined in [51]; hereafter, F22 [51]. First, we trim the ends of the spectrum from each spectral arm to remove noisy data. We then correct for Galactic extinction, using the [37] map and the [38] Milky Way extinction law, and shift to rest-frame wavelengths using the VI redshifts for the SP QSO sample and the redshifts obtained from the modified pipeline for the nominal QSO sample. Each spectrum is then re-binned to a common wavelength grid with 1 Å per bin and normalized to 4500 Å. We then take the median composite, applying a minimum threshold of 30 spectra per bin.

Dust extinction fitting method

To quantify the amount of extinction present in each QSO, we fitted a blue QSO spectral template with varying amounts of dust extinction to the QSO spectra following the same method as F22 [51]. We briefly outline the approach here.

For our blue QSO template, we used the VLT/X-shooter control QSO composite from F22 [51]. This composite was constructed by taking the geometric mean of 28 blue QSOs at $1.45$ $<$ $z$ $<$ $1.65$ (spanning rest-frame wavelengths $1250$ $<$ \uplambda $<$ $9000$ Å) observed with X-shooter (Program ID 0101.B-0739; PI: Klindt). We chose to use this composite rather than constructing a blue composite from DESI spectra due to the broad wavelength coverage of X-shooter; to achieve a similar wavelength range with DESI, we would need to use objects at different redshifts which could result in a non-physical continuum shape (see Appendix D in F22 [51] for more discussion). We then masked the emission-line regions and smoothed the spectrum of the blue QSO X-shooter composite with a Gaussian filter (\upsigma=3). We fitted this unreddened template to the QSO spectra using a least-squares minimization code, which varied the normalization (avoiding emission-line regions) and $E(B-V)$ parameter (ranging from $-1.0$ $<$ $E(B-V)$ $<$ $2.5$ mag¹⁵), using a simple power-law (PL; A_V\propto\uplambda^{-1}, $R_V=4.0$; [51]) dust extinction curve which was found to describe the dust in dust-reddened QSOs more effectively than the commonly used Small Magellanic Cloud (SMC; [24]) extinction law [51].¹⁶ The measured $E(B-V)$ was allowed to be negative to account for the fact that the X-shooter blue QSO composite used for the fitting will contain a small amount of dust and, therefore, some DESI QSOs will be bluer than the composite. However, QSOs with negative $E(B-V)$ values should be interpreted with caution; they represent the natural spread in the power-law slope for unobscured QSOs (e.g., [35]), but do not have any physical interpretation.

In order to determine the quality of the fits, we calculated the Mean Absolute Deviation (MAD $=$ \frac{1}{n}\sum|\rm X-\upmu|; where $n$ is the sample size, X is difference between the model and data, and \upmu is the mean of the difference). For this work, we choose a conservative cut of MAD $<$ $4\times10^{-5}$ to define whether a fit is good or not. Although this will cut down our sample substantially, it should ensure the results in this paper are robust (see Appendix § 12 for a comparison of how different MAD cuts affect our results).

It has previously been noted that QSOs with low luminosities (log L_{\rm 6\upmu m}\lesssim45 erg s$^{-1}$) and redshifts ($z$ $\lesssim$ $1$) can have a substantial contribution of light from the host galaxy [36, 31, 109]. Therefore, for the QSOs at $z\leq1$, we first fitted and subtracted a host-galaxy model from the spectrum to avoid potential stellar emission biasing dust extinction estimates towards the QSO. To fit the host-galaxy component, we used the publicly available code PyQSOFit¹⁷, following the approach taken in F22 [51]. For each source, we first smoothed the spectra with a 10 pixel box-car and then globally fitted both the continuum and emission lines of the entire spectrum. To fit the continuum, we used a power-law, Balmer continuum (BC), third-order polynomial, continuum, and host-galaxy components. Following the same approach as F22 [51], we chose to apply the [40] template for \uplambda $>$ $2000$ Å (obtained via private communication), combined with the [39] template for \uplambda $<$ $2000$ Å (available on https://github.com/legolason/PyQSOFit/blob/master/pyqsofit/fe_uv.txt), to give a continuous template from 1000–12000 Å. The host galaxy and QSO contributions are decomposed using principal component analysis (PCA; [69, 70]). The PCA method is based on the assumption that the observed QSO spectrum is a combination of two independent sets of eigenspectra taken from pure galaxy and pure QSO samples. We used 5 galaxy and 20 QSO PCA components (following a similar approach as [42]), and adopted the stellar host model from [71]. We used this model instead of the default [69] model in PyQSOFit due to the wider wavelength coverage. A host-galaxy component is only included if $>$ $100$ pixels in the host-galaxy template are non-negative; we note the lower resolution of the [71] compared to the [69] template at higher redshifts may lead to fewer spectra successfully fitted by a host-galaxy template. Overall, we found that $\sim$ $10$ per cent (97/948) of the SP QSO sample at $z\leq1$ could successfully be fit with a host-galaxy component. Each of the other continuum components is described in detail in F22 [51]. For examples of our dust fitting method, see Figure 7.

Due to the limited wavelength range of our chosen blue composite ($1250$ $<$ \uplambda $<$ $9000$ Å), our dust extinction fitting method is less effective for sources at $z$ $>$ $2.5$ since this corresponds to a large portion of the DESI spectrum with rest-frame wavelengths $<$ $1250$ Å. The fitting procedure is most effective at $z$ $<$ $2.0$, where the whole DESI wavelength range is covered by the composite. However, to maximize the source statistics for our dust extinction analyses, we include QSOs at $2.0$ $<$ $z$ $<$ $2.5$ that have good fits.

Radio data

In this work we utilize radio data from the LOFAR Two Metre Sky Survey (LoTSS) Data Release 2 (LDR2; [54])¹⁸. This survey covers 5740 square degrees at 144 MHz, with a median sensitivity of 83 \upmuJy/beam and $6''$ resolution. LDR2 contains 4 496 228 sources with a $>$ $5$ $\sigma$ detection. At equivalent frequencies, for sources with typical synchrotron spectra, LDR2 is $\sim$ $9$ times deeper than the Faint Images of the Radio Sky at Twenty-centimeters (FIRST; [34]) 1.4 GHz survey; this makes LDR2 ideal to analyze the radio properties of DESI QSOs, which can have modest AGN luminosities.

To determine whether a QSO is ‘radio-quiet’ or ‘radio-loud’, we adopted the same “radio-loudness” parameter ($\mathcal{R}$) as that first used in [11], and adapted for 144 MHz data in [30], defined as the dimensionless quantity:

The $L_{\rm 144\,MHz}$ is calculated using the methodology described in [7], assuming a uniform radio spectral index of $\alpha$ $=$ $0.7$ for the $K$-correction. We used the rest-frame 6 \upmum luminosity to define radio loudness, rather than optical luminosity which is used in most classical definitions (e.g., [49]), since L6 \upmum is less sensitive to extinction. For example, even at the most extreme $E(B-V)$ for an optically selected QSO ($E(B-V)=1$ mag) at $z=1.5$, the difference in L6 \upmum due to loss of flux by dust extinction at rest-frame 6 \upmum is $\sim$ 0.15 dex. We scaled the radio-loud/radio-quiet threshold from 1.4 GHz ($\mathcal{R}$ $=$ $-4.2$) to our 144 MHz data by assuming a canonical spectral index of $\alpha$ $=$ 0.7 ($S_{\nu}\propto \nu^{-\alpha}$).¹⁹ This corresponded to a boundary at $\mathcal{R}$ $=$ $-4.5$, which is broadly consistent with the classical radio-quiet/radio-loud threshold (often defined using a 5 GHz-to-2500 Å flux ratio; [49]), but is less susceptible to obscuration from dust (see [11] for full details).

Constructing the combined QSO sample

For the scientific analyses in this paper, we constructed a combined QSO sample, which included both our SP and nominal QSO samples in order to maximize the source statistics. We first restricted our nominal and SP QSO samples to the LDR2 survey area using the Multi-Order Coverage map (MOC) provided on the LOFAR website.²⁰ We then restricted the sample to redshifts $0.5$ $<$ $z$ $<2.5$; the lower bound is to ensure minimal contamination from the host galaxy and the upper bound is due to a limitation with our dust fitting method (see § 3.3). This resulted in 690 and 33 753 QSOs in the SP and nominal samples, respectively. We then matched to the LDR2, adopting a $5''$ search radius which resulted in $\sim$ $19$ per cent (6389/33 753) and $\sim$ $31$ per cent (216/690) radio-detected QSOs in the nominal and SP samples, respectively (we further explore the radio detection fractions in § 4.2). Finally, we combined both the restricted nominal and SP QSO samples (removing duplicate QSOs that belong to both samples), resulting in our “combined QSO sample” which consists of 34 293 QSOs (see Figure 4).

To estimate the fraction of false matches and objects missed due to radio emission extended beyond $5''$, we utilized the LoTSS DR1 radio catalogue (LDR1; [21]), with optical counterparts from [89].²¹ The catalogue from [89] was created using a combination of likelihood-ratio matching and visual identification, resulting in optical and/or IR counterparts for $\sim$ 73 per cent of the LDR1 radio sources. For our false matching estimate, we will use the optical/IR positions of the LDR1 sources as the “true” positions. Restricting the nominal QSO sample to the much smaller region of LDR1 reduced the sample to only 3057 sources. Adopting the same $5''$ search radius as for LDR2, we matched this reduced sample to the radio positions in LDR1, resulting in 645 matches. We then repeated the matching, but instead adopting a 0\farcs75 search radius (following [30]) and matching to the positions of the optical counterparts from the value-added catalogue, resulting in 648 matches which represent our “truth” sample. From this, we were able to determine an estimate for the fraction of false matches by exploring how many of the 645 radio matches within $5''$ were not also included in our truth sample; from this we found a false match fraction of $\sim$ 1 per cent (8/645). Similarly, we were able to determine the fraction of missed matches by exploring how many of the truth sample were not included in the sample of 645; from this we found that $\sim$ 2 per cent (11/645) of sources are missed, most likely due to radio emission extended beyond $5''$. Although the false match fraction calculated here is greater than that found for FIRST (using a $10''$ matching radius: $\sim$ 0.2 per cent; [11]), which is at a comparable resolution to LoTSS (FIRST: $5''$; LoTSS: $6''$), the sensitivity of LoTSS is $\sim$ 9 times deeper than FIRST at equivalent frequencies, assuming a spectral index of $\alpha$ $=$ $0.7$. We would therefore expect a higher false match fraction from LoTSS due to the higher density of faint sources within a given region. Additionally, our current estimations are limited by source statistics. Our future DESI papers will benefit from the upcoming optical and IR cross-matched LDR2 catalogue for more accurate cross-matching and also a larger DESI SP QSO sample.

Results

Our results are based on the combined QSO sample (see Figure 4) that includes QSOs selected from our secondary target program (§ 2.1) and QSOs that were targeted by the nominal DESI QSO program, restricted to the LoTSS DR2 coverage area (to assess radio detection fractions). In the following sections, we explore the line-of-sight dust extinction towards the QSOs (§ 4.1) and the radio properties of the QSOs as a function of dust extinction (§ 4.2).

Quantifying the line-of-sight dust extinction

Under the assumption that dust is the main cause of reddening in reddened QSOs, we quantified the amount of dust extinction along the line-of-sight by fitting the spectra of the combined QSO sample with a dust-reddened VLT/X-shooter blue QSO composite (see § 3.3).

Overall, 58 per cent (20 031/34 293) of the combined QSO sample spectra are well fitted (MAD $<$ $4\times10^{-5}$; see § 3.3) by the dust-reddened composite. Figure 7 displays the best fitting solution for one of the most extreme cases, with a measured $E(B-V)$ $=$ $1.05$ mag ($A_V$ $\sim$ $4$ mag); for this source the Mg ii emission line is completely obscured, although broad H $\beta$ is still present indicating that it is a QSO. Inspecting the QSOs that are not well fit with a dust extinction curve, we found that they tend to have bluer optical colours (median $g-r=0.26$ and 0.21 mag for the QSOs with good and bad fits, respectively), fainter optical magnitudes (median $r=20.0$ and 20.6 mag, for the QSOs with good and bad fits, respectively), display strong absorption features, and/or reside at higher redshifts (median redshift $=1.22$ and $1.83$ for the QSOs with good and bad fits, respectively); as expected based on the lower effectiveness of our fitting procedure at $z$ $>$ $2$ due to the limited wavelength range of our chosen blue composite (see § 3.3). We therefore expect the QSOs with good fits to be slightly biased towards either the more luminous QSOs and/or those with higher dust extinctions. Additionally, from visually assessing the bluest QSOs ($E(B-V)$ $<$ $-0.1$ mag) with a MAD $<$ $4\times10^{-5}$, we found the fits to be unreliable (i.e., the extinction of the QSO compared to the blue composite appears to be consistent with zero, but the measured $E(B-V)$ is clearly underestimated), again highlighting the limitation of our fitting towards blue QSOs. However, due to the significant sample size of the bluest QSOs we expect this to have little effect on our results (see Appendix § 12). For the following analyses, we restrict to QSOs with good dust extinction fits and also with a measured $E(B-V)$ $>$ $-0.1$ mag. Examples of QSOs with good and bad dust extinction fits can be found in the online supplementary material.

Image not extracted

Figure 8 displays the distribution of dust-extinction values for the combined QSO sample. The 107 SP QSOs that were additionally identified in the SDSS are also highlighted. This demonstrates that DESI will not only accumulate a more statistically significant sample of reddened QSOs but also extends to QSOs with dust extinctions beyond those observed by the SDSS. Figure 9 displays the $g-r$ versus $r-z$ optical colour-colour space, with the colourbar representing the measured $E(B-V)$. We find that both the $g-r$ and $r-z$ colour scale with $E(B-V)$, demonstrating that optical colour is a good proxy for dust reddening. In order to maximize the highest $E(B-V)$ end of the QSO population, future QSO selections can utilize this correlation between optical colour and dust extinction in order to select the reddest QSOs.

Overall, DESI can identify QSOs with dust extinctions up to $E(B-V)$ $\sim$ 1.0 mag ($A_V$ $\sim$ $4$ mag; see Figure 7 for an example fit), comparable to the much smaller sample of NIR-selected heavily reddened QSOs from [4]. This demonstrates the capability of DESI to help bridge the gap between mildly obscured SDSS red QSOs (e.g., [35, 11, 51]), IR-selected red QSOs (e.g., [4, 2, 53]), and heavily obscured QSOs (e.g., [94]). See § 5.1 for a more detailed comparison with previous red QSO studies. The measured $E(B-V)$ values for the combined QSOs included in the DESI EDR can be found in the electronic table online.

Exploring the radio properties of dusty QSOs

In our previous work, we identified fundamental differences in the radio properties of red QSOs compared to typical blue QSOs based on the shallower SDSS QSO samples [11, 30, 29, 22]. In this paper we utilize the combined QSO sample, our spectral fitting results, and LDR2 (see § 3.4 and § 3.5) in order to explore the radio detection fraction as a function of reddening and extend our previous radio results to more heavily reddened QSOs.

$E(B-V)$	Combined QSO sample (all)	Combined QSO sample (L_{\rm6\,\upmu m}--$z$matched)
	# QSOs	Radio-detection frac. (#)	Radio-loud frac. (#)	# QSOs	Radio-detection frac. (#)	Radio-loud frac. (#)
-0.1--0.0	7412	17.2% (1275)	1.2% (87)	861	18.9% (163)	0.8% (7)
0.0--0.1	7720	20.1% (1552)	1.6% (127)	2529	20.3% (514)	1.5% (38)
0.1--0.2	2725	22.7% (618)	1.3% (36)	843	24.4% (206)	1.3% (11)
0.2--0.3	1068	29.7% (317)	0.9% (10)	739	27.6% (204)	0.9% (7)
0.3--0.4	431	29.7% (128)	0.2% (1)	297	30.3% (90)	0.3% (1)
0.4--1.1	329	38.0% (125)	0.9% (3)	297	37.0% (110)	1.0% (3)

In Figure 10 we show the LDR2 radio detection fraction for the combined QSO sample in bins of $E(B-V)$. On the basis of this figure it is clear that there is a strong positive correlation between the line-of-sight dust extinction and the radio detection fraction for DESI QSOs. This is consistent with our previous work in which we found that red QSOs in SDSS, selected by their $g-i$ optical colour, have a higher radio detection fraction compared to typical blue QSOs [11, 30, 29, 22]. To test that there are no redshift or luminosity biases that are driving this result, we matched the QSOs within each $E(B-V)$ bin in redshift and L_{\rm 6\,\upmu m} (using a tolerance of 0.05 and 0.2 dex for redshift and L_{\rm 6\,\upmu m}, respectively), using a 2:1 match ratio for the bluest bins to maximize source statistics. We show the LDR2 detection fraction for each matched bin in Figure 10 by the grey shaded region; although limited on source statistics, a strong positive trend is clearly still present, even when accounting for redshift and luminosity biases (see Appendix § 12 for additional tests). The number of QSOs and associated radio detection fractions for the full and L_{\rm 6\,\upmu m}–redshift matched samples can be found in Table 2. It is also worth noting that no radio pre-selection is incorporated in the DESI QSO target selection, in comparison to SDSS where some QSOs are selected based on a FIRST radio detection (our previous work addressed this pre-selection; see Section 2.2.1 in [11]). This result therefore corroborates our previous SDSS red QSO results, but utilizing a cleaner, less biased sample and expanding to QSOs with even higher line-of-sight dust extinction (SDSS red QSOs span the range $0.03$ $\lesssim$ $E(B-V)$ $\lesssim$ 0.3; [35, 11, 29, 31]).

To give some basic insight on the radio properties, in Figure 11 we show the L_{\rm 6\,\upmu m} versus $L_{\rm 144\,MHz}$ distribution for the radio-detected combined QSO sample in three bins of $E(B-V)$, with the radio-loud/radio-quiet threshold displayed by the dashed red line (see § 3.4); we find that 93 per cent of the radio-detected QSOs in the combined sample are radio-quiet. Assuming a $5\sigma$ upper limit for the LDR2 undetected QSOs (415 \upmuJy) we estimate a radio-quiet fraction for the overall combined QSO sample of $\sim$ 99 per cent. Therefore, the trend in radio detection fraction found in Figure 10 is driven by radio-quiet QSOs, consistent with our previous work. We do not find a significant difference in the median radio-loudness across increasing $E(B-V)$ bins ($\mathcal{R}$ $\sim$ $-5.9$ for each $E(B-V)$ bin). However, we do find a decrease in the number of radio-loud QSOs with increasing $E(B-V)$ ($\sim$ $8$, 5, and 3 per cent with increasing $E(B-V)$ bin); we explore this result more in § 5.2.

Discussion

We have presented the first results from our ongoing secondary target program in DESI that uses an optical and MIR colour selection to expand the optical colour space probed by the nominal DESI QSO selection. Our program complements the QSOs from the nominal DESI QSO selections, which are predominantly blue, and pushes to higher extinctions than those selected by both the SDSS and DESI QSO samples. Combining the QSOs from our secondary target program and the nominal DESI QSOs, the key result in this paper is the striking positive relationship between amount of line-of-sight dust extinction and the radio detection fraction for QSOs (Figure 10). This clearly demonstrates that the presence of dust (and/or gas, i.e., opacity) is an important factor for the production of radio emission in radio-quiet QSOs. In the following sections we discuss how these results complement and extend the previous red QSO studies, including the overlap between our sample and Extremely Red Quasars (ERQs; e.g., [84]; § 5.1), and explore the origin of the radio emission in radio-quiet dusty QSOs (§ 5.2).

A comparison to previous red and obscured QSO studies

Sample name	Key selection	$E(B-V)$	# sources	Reference
		[mag]
DESI combined QSOs	DESI QSO selections &$W123$	0.1--1.1	4553	This work
		0.25--1.1	1177	This work
SDSS red QSOs	DR7$(g-i)$90th percentile	$\lesssim0.2$	4899
	DR14$(g-i)$90th percentile		21 800
WISE red QSOs	$W1-W2$$>$0.7 &$W2-W3$$>$$2$&$W3-W4$$>$$1.9$&$E(B-V)$$>$$0.25$	0.25--1.2	21
W2M red QSOs	$0.5$$<$$W1-W2$$<2$&$2$$<$$W2-W3$$<$$4.5$&$E(B-V)$$>$$0.25$	0.25--1.2	63
HRQs	$J-K$$>$$1.6$&$W1-W2$$>$0.85	0.5--1.5	38
ERQs	$i-W3$$>4.6$	$\lesssim0.7$	97
Obscured	SED fitting	$\gtrsim2$	355
IR red QSOs		0.05--1	223

SDSS red QSOs

In our previous work exploring SDSS selected red and blue QSOs, we found that the majority of the red QSOs were red due to modest amounts of dust extinction ($E(B-V)\lesssim0.2$ mag; [31]; F22 [51]). For example, utilizing VLT/X-shooter data, F22 [51] found that the red QSO spectra displayed the characteristic frequency dependent spectral turnover at UV–optical wavelengths, indicating they were red due to line-of-sight dust extinction. To illustrate the improvements in the selection of dust-reddened QSOs by our DESI program, in Figure 12 we show a comparison between the X-shooter red and blue composites from F22 [51] and our combined DESI sample, in similar luminosity–redshift matched bins of measured $E(B-V)$ as in Figure 10 (see § 3.3). Firstly, we confirm that the DESI composites appear more dust reddened (i.e., a stronger UV spectral turnover) with increasing $E(B-V)$, providing confidence in our dust extinction fitting method. We also find that the reddest DESI composite has over $\sim$ 0.15 mag more dust extinction on average compared to the F22 [51] SDSS red composite, again demonstrating that with DESI we can obtain a statistically significant sample of luminous heavily reddened QSOs that are missed by the optically shallow SDSS survey.

We have also previously found that SDSS red QSOs display enhanced radio emission compared to typical blue QSOs ([11, 30, 29, 22]). Similarly, [53] used a WISE–2MASS colour selection to define a sample of red QSOs with dust extinctions $E(B-V)$ $>$ $0.25$ mag ($A_V$ $\gtrsim$ $0.8$ mag; by this definition, the SDSS QSOs defined as “red” in our previous work would be classified as “blue” in the [53] sample) and found a higher radio detection fraction for their more heavily reddened red QSOs compared to their blue sample (see Figure 18 therein), using both FIRST and VLASS. These results are in agreement with Figure 10, where we have now extended to more extreme red QSOs, which suggests that the radio enhancement in red QSOs is linked to a higher amount of line-of-sight dust extinction. Table 3 displays the number of combined QSOs with SDSS-like reddening ($E(B-V)$ $>$ $0.1$ mag) and the number that satisfy the [53] red QSO cut ($E(B-V)$ $>$ $0.25$ mag), resulting in 4533 and 1177 red QSOs, respectively. This demonstrates that the DESI survey will obtain a statistically significant sample of heavily reddened QSOs.

Obscured QSOs

Extending our previous red QSO work to fully obscured QSOs (e.g., Compton thick), [94] compared the VLA 1.4 and 3 GHz radio properties of a sample of obscured and unobscured QSOs in the COSMOS field (by their definition, even our extreme DESI red QSOs are unobscured) and found that the obscured QSOs have higher radio detection fractions compared to unobscured QSOs (see Figure 14 therein). However, interestingly they found that the obscured QSOs have consistent radio detection fractions to red QSOs. This could suggest that either 1) there is a “saturation” point in Figure 10, whereby increasing the amount of line-of-sight dust extinction towards a QSO will not result in a higher radio detection fraction, 2) the origin of the radio emission is different in obscured and unobscured QSOs, or 3) the differences in the radio detection fraction are more modest at the higher sensitivity radio data in COSMOS. Utilizing the increased source statistics in future DESI data releases we can populate the highest $E(B-V)$ end of Figure 10, testing whether the radio detection fraction increases or does indeed start to plateau, as would be the case in the former scenario; by the end of the five year survey we expect to add additional bins at $E(B-V)>0.4$ in Figure 10.

Extremely Red QSOs

Since our DESI study pushes to more reddened QSOs than the previous SDSS red QSO studies, it is interesting to investigate the overlap between our DESI red QSOs and the population of Extremely Red Quasars (ERQs; see Table 3 for a comparison of the red QSO selections discussed in this section), which are selected based on their optical–MIR colours, as compared to the pure optical or NIR selection approach often used to select red QSOs [1, 84]. ERQs have been found to host some of the most powerful outflows identified in QSOs (e.g., [32, 75, 83]), and therefore have also been suggested as a candidate transitional population in the evolution of QSOs (see discussion in § 5.2). However, due to the depth of SDSS, previous studies have only analyzed a small sample of the most luminous ERQs. In order to assess the capability of DESI (in particular, our SP) at selecting ERQs, we first defined an ERQ selection based on the optical–MIR colours of the SDSS ERQ sample from [84]. [84] defined an ERQ to have redshifts $>2.0$, $W3$ SNR $>3$, and $i-W3\geq4.6$ mag. For our ERQ classification, we based our selection on the full [84] ERQ sample that meet the $i-W3\geq4.6$ mag selection; 141 objects. Since the DR9 photometry does not include an $i$-band magnitude, we defined our DESI ERQ selection based on the $z-W3$ colours, qualitatively consistent with the ERQ selection of [1] and [84]; we classify a QSO as an ERQ if it has an optical–MIR colour that satisfies $z-W3>3.9$ mag. This resulted in $\sim$ $39$% (714/1852) of our SP QSO sample classified as ERQs, in comparison to only $\sim$ $7$% (6167/88 101) for the nominal QSO sample, demonstrating that our SP is also efficient at selecting ERQs; in future data releases we will build up a statistically significant sample of these rare systems. Our future papers will explore the dust extinction, emission-line, and radio properties of DESI ERQs and how they compare to the general red QSO population.

Mg II Absorption Line Systems

A number of studies have noted that a larger fraction of reddened QSOs host narrow absorption lines compared to typical blue QSOs [35, 118, 51]. The strongest and most distinctive of these absorption lines is the Mg ii$\lambda\lambda 2796,2803$ doublet, with a laboratory separation between the two lines of 7.1772 Å [121]. Mg ii absorption that is close to the Mg ii emission line (i.e., where $z_{\rm abs}$ $\sim$ $z_{\rm em}$) is thought to be associated with the host galaxy of the QSO, whereas blueshifted absorption (i.e., where $z_{\rm abs}$ $<<$ $z_{\rm em}$) is thought to be due to an intervening system. Previous studies have found that QSOs that host associated absorption lines (AALs) are, on average, redder, with stronger [O ii] emission, and are more likely to be radio detected, compared to both QSOs with intervening absorption lines and QSOs with no absorbers [119, 43, 120, 118]. The strength of the Mg ii absorption in systems with AALs has also been found to correlate with the dust extinction and radio properties, suggesting an intrinsic nature for the AALs [119, 118]. Despite the similarities in radio properties between QSOs with AALs and red QSOs, the measured dust extinctions are far more modest in the AAL QSOs, with an $E(B-V)$ $\sim$ $0.03$–0.06 mag [118]. Furthermore, [43] suggested that AAL QSOs may be a transitional population between heavily reddened QSOs and typical blue QSOs. In order to assess whether AAL QSOs are driving the relationship between reddening and the radio detection fraction found in Figure 10 and our previous SDSS work [11, 29, 30], we explored the properties of QSOs with Mg ii absorbers in our DESI combined QSO sample, utilizing the catalogue from [117].

[117] used an autonomous spectral pipeline to detect associated or intervening Mg ii absorption through an initial line-fitting process and a Markov Chain Monte Carlo (MCMC) sampler. To create the Mg ii absorber catalogue, they ran their detection pipeline on 83 207 DESI QSOs at $0.3$ $<$ $z$ $<$ $2.5$ which are part of the DESI EDR. They found 20.1 per cent (16 707/83 207) of the unique QSO spectra have Mg ii absorbers, with a large number that contain multiple absorption features.

Restricting the combined QSO sample to those included in the EDR significantly reduces the sample to 11 092 QSOs. Matching to the [117] Mg ii catalogue resulted in 4924 absorption systems in 30.7 per cent (3400/11 092) unique QSO spectra. Of these, 16.4 per cent (808/4924) are associated, defined by $-6000$ $<$ $v_{\rm off}$ $<$ $5000$ km s$^{-1}$ ($v_{\rm off}=c\times\frac{z_{\rm MgII}-z_{\rm QSO}}{1+z_{\rm QSO}}$; [117]), and the rest are intervening. Comparing the radio detection fractions, we found 24.3 per cent, 13.4 per cent, and 13.4 per cent of the AAL QSOs, intervening absorber QSOs, and QSOs without any absorption systems have LDR2 radio detections, respectively. We also found the AAL QSOs are considerably redder than the intervening and no absorber systems, with a median $E(B-V)$ $\sim$ $0.12$, 0.00, and 0.03 mag for the QSOs with associated, intervening, and no absorbers, respectively (see also Napolitano et al. in prep). Since the median $E(B-V)$ associated with the absorption feature has been found to be $\lesssim$ $0.03$ mag, this suggests that AALs are more likely to be associated with red QSOs. This is further shown in Figure 13, which displays the fraction of associated and intervening Mg ii absorbers as a function of $E(B-V)$; a larger fraction of associated absorbers are found in QSOs with higher amounts of dust extinction, consistent with previous studies [43, 118, 51]. This suggests the Mg ii absorption is associated with the dust/gas that is also responsible for reddening the spectrum. We found that the fraction of QSOs with intervening absorbers decreases with increasing $E(B-V)$; by definition, intervening absorbers will fall at the bluer end of the spectrum (relative to Mg ii), which will have a lower signal-to-noise for red QSOs due to dust extinction, and so this could reduce the effectiveness of the automated pipeline at detecting absorption features.

Figure 13 also displays the radio-detection fraction for the sample with and without both the intervening and associated absorber systems. We find that removing the absorber systems has little affect on the overall trend with $E(B-V)$, which is consistent with Figure 10 albeit with lower source statistics, suggesting that the amount of dust (and/or gas, i.e., opacity) in a QSO is the most important factor for determining whether a QSO is radio detected or not. However, since QSOs with more dust are more likely to host an AAL, QSO samples selected purely on the presence of Mg ii AALs are likely to be redder and more radio-detected compared to QSOs without absorbers.

Origin of the radio emission in dusty radio-quiet QSOs

The strong relationship between the line-of-sight dust extinction and radio detection fraction displayed in Figure 10 demonstrates an intrinsic link between dust and the production of radio emission in QSOs. However, it is currently unclear what is the dominant mechanism for the production of the radio emission. Calculating the fraction of radio-loud ($\mathcal{R}$ $>$ $-4.5$; see § 3.4) QSOs in the combined sample as a function of $E(B-V)$, we find that $\lesssim$ $2$ per cent of the QSOs (both radio-detected and undetected) are radio-loud across all bins (see Figure 10 and Table 2). This demonstrates that radio-loud QSOs, typically hosting a large-scale radio jet, are not driving the correlation between dust extinction and radio-detection fraction, in agreement with our previous red QSO studies [11, 29, 30].

The low frequency radio emission probed by LOFAR is dominated by synchrotron emission, which may be produced by star formation, low-powered/frustrated jets, and/or winds causing shocks in a dusty circumnuclear/ISM region (for more details on the production of radio emission in radio-quiet QSOs, see review by [17]). To perform a detailed analysis into the origin of the radio emission, multi-frequency and/or higher spatial resolution radio observations are required; our future radio study utilizing radio data from the upgraded Giant Meter Radio Telescope (uGMRT; PI: V. Fawcett), in addition to high spatial resolution e-MERLIN data [19], aims to do just this. In the following section we discuss some of the potential scenarios that could be driving the connection between the line-of-sight dust extinction and radio detection fraction in DESI QSOs.

In our previous work exploring the properties of SDSS red QSOs, we found no differences in the star-formation properties between the red and blue QSOs. This was based on deep IR data from the COSMOS region [29] and composite SEDs [31]. We therefore concluded that the radio emission in SDSS red QSOs was likely driven by AGN processes such as winds or jets. However, given that DESI pushes to redder systems, it is unclear whether this conclusion still holds for our combined QSO sample. In addition to comparing obscured and unobscured QSOs in COSMOS, [94] also compared the star-formation properties between a blue ($E(B-V)$ $<$ $0.05$ mag) and red ($0.05$ $<$ $E(B-V)$ $<$ $1.0$ mag) unobscured sub-sample. They found no difference in the star-formation properties between the two samples, despite finding a significantly higher radio-detection fraction for the red QSOs compared to the blue QSOs. Given that the red and blue QSOs from [94] cover a similar $E(B-V)$ range to Figure 10, this suggests that star formation is not a strong driver of the relation between dust extinction and radio detection fraction. Furthermore, Calistro Rivera et al. (in prep) perform comprehensive SED fitting of red and blue QSOs in the LOFAR deep fields and find no difference in the star-formation properties. However, individual star formation measurements are required to robustly demonstrate this for our sample. Calculating the $L_{\rm 144\,MHz}$ from star-forming galaxies with a star-formation rate of 100 M$_{\rm \odot}$yr$^{-1}$ (equivalent to local ultraluminous infrared galaxies; ULIRGs), based on the FIR-radio correlation from [114], we find that only $\sim$ 16 per cent of the radio-detected QSOs may have a substantial contribution from star formation to the radio emission, $\sim$ 83 per cent of which are at $z$ $<$ $1$. Therefore, the most likely mechanism that links the presence of dust in the QSOs to the radio emission is either frustrated jets and/or winds causing shocks in a dusty circumnuclear/ISM region.

Recent hydrodynamic simulations have predicted that radiation pressure on dust is an effective way to launch outflows (e.g., [95]). The optical and UV radiation from the accretion disc is absorbed and re-emitted by the dust as IR radiation, transferring momentum to the ambient gas [98, 99]. These AGN-driven outflows can have high velocities ($\sim$ $1000$–3000 km s$^{-1}$), extend up to tens of kpc [97], and are likely launched along the polar axis (e.g., [100]). When the outflow interacts with the surrounding dense ISM, they can cause shocks which accelerate electrons, producing radio emission [102, 101, 91]. Dusty winds are therefore a promising mechanism for the production of radio emission in radio-quiet red QSOs. This is supported by a number of studies that have observed powerful outflows in red QSOs [18, 74, 75, 31, 33, 90], although it is not always clear whether this is a result of luminosity biases (e.g., [45]; F22 [51]).

On the other hand, low-powered jet-ISM interactions have also been found to be important in radio-quiet QSOs [15, 96, 111, 115]. For example, jet-ISM simulations have found that an initially collimated radio jet can get spread into lower powered jet filaments, which can potentially drive shocks in the ISM [112, 113]. A significant fraction of the jet’s energy (up to $30$ per cent; [112]) can be imparted to the surrounding gas clouds, driving outflows. This is in agreement with observational studies that find correlations between the strength and/or spatial location of the ionized outflow with a young or weak radio jet [104, 107, 92, 103, 105]. The connection between these low-powered radio jets and the expanding ionized gas bubbles have been shown to result in both negative and positive feedback in the host galaxy [106, 15, 96], demonstrating the importance for studying the radio-ISM interactions in these radio-quiet systems.

With our current data it is difficult to distinguish between the jet versus wind scenarios (or whether it is some combination of both). However, in either case the resulting effect is a jet/wind shocking the surrounding ISM. This will likely drive out the surrounding gas/dust in a short-lived “blow-out” phase (e.g., [27]) that may eventually reveal a typical blue QSO [9]. Therefore, the relationship between dust and radio emission in QSOs may be crucial to our understanding of QSO evolution.

In future DESI papers we will robustly compare the emission line and radio properties of a statistically significant sample of QSOs as a function of $E(B-V)$, in order to understand to what extent differences in the observed outflow properties and/or accretion state drive the striking relationship in Figure 10. In addition to our future DESI studies, we have obtained four-band radio data observed with the uGMRT (PI: V. Fawcett) for a sample of 20 red and 20 blue radio-detected QSOs with high spatial resolution e-MERLIN 1.4 GHz data from [19]. With these data we can construct sensitive radio SEDs over 0.12–3 GHz (with the addition of VLASS), in order to search for spectral signatures of radio jets (e.g., [113]). For example, if the QSOs are similar to compact steep spectrum (CSS) and gigahertz-peaked spectrum (GPS) sources (see review by [47]), which are thought to be jet dominated, then we can expect a spectral turnover in the radio SED around $\sim$ $100$ MHz and $\sim$ $1$ GHz (for CSS and GPS sources, respectively). Therefore, we can explore 1) whether the radio SEDs of red and blue QSOs are different, and 2) if jets are the dominant mechanism for producing radio emission in radio-quiet dusty QSOs. Furthermore, future Very Long Baseline Interferometry (VLBI) observations will be crucial to resolve and confirm potential pc-scale radio jet morphologies in these systems.

Conclusions

In this paper, we present the first $\sim$ $8$ months of our DESI secondary target program that utilizes a MIR and optical selection to expand the DESI QSO colour space, in particular to reddened QSOs. Combining our SP QSO sample with the larger nominal DESI QSO sample, we have explored the dust extinction and radio properties. From our analyses we find that:

• Our secondary target program extends the nominal DESI QSO program to more extreme optical colours, the majority of which are red: From visually inspecting the spectra of the 3038 objects from our secondary target program, we find that $\sim$ $89$ per cent have high quality spectral and redshift classifications. Of these high quality sources, 81 per cent are QSOs (both Type 1 and Type 2), 13 per cent are galaxies, and 6 per cent are stars (Figure 3). Overall, our program provides $\sim$ $32$ per cent more QSOs that lie outside of the main colour space compared to the nominal QSO sample despite being $\sim$ $2$ per cent of the size (Figure 5), and also pushes to fainter optical magnitudes than the nominal QSO program (Figure 6). See § 2, § 3.1 and § 5.1.

• DESI is capable of detecting QSOs with up to $E(B-V)$ $\sim$ $1$ mag of dust extinction: Measuring the dust extinction by fitting the spectra of the combined QSOs with a dust-reddened blue QSO composite, we find up to $\sim$ 1 mag in $E(B-V)$ ($A_V$ $\sim$ 4 mag; Figure 8), comparable to NIR-selected heavily reddened QSOs. We also find that optical colour strongly correlates with dust extinction, demonstrating that optical colour is an effective way to select dust-reddened QSOs (Figure 9). Finally, splitting the combined QSO sample into $E(B-V)$ bins, we find that the reddest composite ($E(B-V)$ $>$ $0.2$ mag) has over $\sim$ 0.15 mag more dust extinction than the [51] SDSS red QSO composite, again demonstrating that with DESI we can push to much higher amounts of dust extinction than in SDSS (Figure 12). See § 4.1 and § 5.1.

• There is a strong link between the radio detection fraction and amount of line-of-sight dust extinction in QSOs: Utilizing the combined QSO sample to maximize the source statistics, we explored the 144 MHz radio properties of the QSOs using the LoTSS DR2 (LDR2) dataset. Comparing the LDR2 radio detection fraction of the combined QSO sample as a function of $E(B-V)$ we find a strong positive correlation (Figure 10), even when accounting for any redshift or luminosity biases, demonstrating that the amount of line-of-sight dust extinction is closely linked to the radio emission in a QSO. We find that $\lesssim$ $2$ per cent of the combined QSO sample are radio-loud (Figure 11), suggesting that large-scale radio-loud jets are not driving the relation between radio-detection fraction and dust extinction. Therefore, we suggest the origin of this predominantly radio-quiet emission is likely due to frustrated jets or outflows causing shocks in a dust/gas rich ISM, which may drive away the dust in a “blow-out” phase. See § 4.2 and § 5.2.

We have demonstrated that DESI will observe a statistically significant sample of dust-reddened QSOs, bridging the gap between the more modest SDSS red QSOs and the heavily reddened NIR-selected QSOs. From our results we have shown that red QSOs are likely undergoing a dusty blow-out phase, in which winds and/or jets clear out the surrounding dust and gas, eventually revealing a blue, typical QSO. Therefore, red QSOs are important objects to study in the context of AGN feedback and QSO evolution. In future papers we will use our growing DESI sample to investigate whether there are any differences in the outflow properties of QSOs as a function of dust extinction and to explore statistically significant sub-samples of rare/exotic QSOs (e.g., Extremely Red Quasars and Low-ionization Broad Absorption Line QSOs).

Acknowledgments

We would like to thank the anonymous referee for their constructive comments that greatly helped improved the clarity and structure of the paper.

VAF and CMH acknowledge funding from an United Kingdom Research and Innovation grant (code: MR/V022830/1) and (VAF) previously a quota studentship through grant code ST/S505365/1 funded by the Science and Technology Facility Council. DMA and ACE acknowledges the Science and Technology Facilities Council (through grant code ST/P000541/1). AB supported by the U.S. Department of Energy, Office of Science, Office of High Energy Physics, under Award No. DESC0009959. DR acknowledges the support of STFC grant NU-012097. ADM was supported by the U.S. Department of Energy, Office of Science, Office of High Energy Physics, under Award Number DE-SC0019022. We would like to thank E. Verner for providing us with the template presented in [40].

This work uses observations collected at the European Southern Observatory under ESO programme 0101.B-0739(A).

LOFAR data products were provided by the LOFAR Surveys Key Science project (LSKSP; https://lofar-surveys.org/) and were derived from observations with the International LOFAR Telescope (ILT). LOFAR (van Haarlem et al. 2013) is the Low Frequency Array designed and constructed by ASTRON. It has observing, data processing, and data storage facilities in several countries, which are owned by various parties (each with their own funding sources), and which are collectively operated by the ILT foundation under a joint scientific policy. The efforts of the LSKSP have benefited from funding from the European Research Council, NOVA, NWO, CNRS-INSU, the SURF Co-operative, the UK Science and Technology Funding Council and the Jülich Supercomputing Centre.

This material is based upon work supported by the U.S. Department of Energy (DOE), Office of Science, Office of High-Energy Physics, under Contract No. DE–AC02–05CH11231, and by the National Energy Research Scientific Computing Center, a DOE Office of Science User Facility under the same contract. Additional support for DESI was provided by the U.S. National Science Foundation (NSF), Division of Astronomical Sciences under Contract No. AST-0950945 to the NSF’s National Optical-Infrared Astronomy Research Laboratory; the Science and Technologies Facilities Council of the United Kingdom; the Gordon and Betty Moore Foundation; the Heising-Simons Foundation; the French Alternative Energies and Atomic Energy Commission (CEA); the National Council of Science and Technology of Mexico (CONACYT); the Ministry of Science and Innovation of Spain (MICINN), and by the DESI Member Institutions: https://www.desi.lbl.gov/collaborating-institutions. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the U. S. National Science Foundation, the U. S. Department of Energy, or any of the listed funding agencies.

The authors are honored to be permitted to conduct scientific research on Iolkam Du’ag (Kitt Peak), a mountain with particular significance to the Tohono O’odham Nation.

Data availability

The SV data and spectra underlying this article are publicly available in the DESI Early Data Release [124], with the remaining data and spectra obtained from the Main survey publicly available in a future data release (DESI collaboration et al. in preparation 2024).

The radio and VLT/X-shooter data used in this study are publicly available in online archives.

The data underlying the figures in this paper are available online: https://doi.org/10.5281/zenodo.8147342.

Electronic tables containing the DESI TARGETID for the full combined QSO and SP samples can be found online. Additional redshift and photometric information are provided for the sources included in the DESI Early Data Release [124].

Testing the performance of the modified pipeline

Comparing the redshifts and spectral types (“spectypes”) for the DESI SP sample obtained from the VI to the modified pipeline, we found that the modified pipeline correctly identifies $\sim$ $71$ per cent of all spectypes (see Table. tab:spec_tab_comp). We also found that $>$ $72$ per cent of our sample are high-quality QSOs compared to $\sim$ $47$ per cent classified by the modified pipeline; this is predominantly due to the number of reddened QSOs in our sample (see § 2.1 and [61]) and highlights the necessity of VI to reliably classify red QSOs in our sample.

A comparison between the resulting high-quality VI redshifts and those from the modified pipeline for our DESI SP sample is displayed in Figure 14. We found $\sim$ $66$ per cent (1795/2700) of all high-quality objects have correct redshifts ($\Delta z=|z_{\rm VI}-z_{\rm RR}|/(1+z_{\rm VI})<0.01$, which corresponds to a $\Delta v\leq3000$ kms$^{-1}$) as measured by the modified pipeline. In the SP QSO sample, $\sim$ $58$ per cent (1076/1852) have correct redshifts from the modified pipeline. Of the other $\sim$ $42$ per cent, $\sim$ $78$ per cent were also incorrectly classified as galaxies; these systems are all assigned a redshift $<$ $1.7$ by the modified pipeline, due to the maximum redshift that the modified pipeline can assign to a galaxy. The fraction of correct redshifts in our SP QSO sample is a lot lower than that found in the VI of the nominal DESI QSOs ($\sim$ $99$ per cent); this is again predominantly due to the high number of dust-reddened QSOs in our sample, which are substantially more likely to be miss-classified [61].

Population stacks

Figure 15 displays the composite of the 358 and 306 high-quality galaxies and Type 2 AGN in our DESI SP sample (see § 2.1), following the approach outlined in § 3.2. The Type 2 AGN display strong high ionization line emission lines, including [Ne v]\uplambda\uplambda3346,3426, that implies the presence of hard radiation from an AGN [116]. In the galaxy composite, the emission lines are narrower compared to the Type 2 composite and the [Ne v] lines are faint-to-non-existent, demonstrating little-to-no AGN activity in these objects. For a more detailed discussion of the galaxies in DESI, see [78, 64, 126, 79, 125, 80, 66]. The broader [O iii] in the Type 2 composite compared to the galaxy composite implies that the Type 2 AGN are more likely to host powerful outflows; future work will investigate the emission line properties of our SP sample.

Figure 16 displays the composite of the 164 high-quality stars in our DESI SP sample, split into “red” ($r-z$ $>$ $-0.4$) and “blue” ($r-z$ $<$ $-0.4$). The red stars show strong molecular features that are indicative of late-type cool stars (e.g., M-type), whereas the blue stars are mainly (DA) white dwarfs due to the broad and shallow Balmer features (see Figure 11 in [122]).²² For a more detailed discussion of the stars in DESI, see [81, 59].

SDSS and DESI comparison

	Redshift	Spectype
SDSS ID	SDSS	DESI VI	SDSS	DESI VI
131729.28+545720.6	5.07	1.65	QSO	QSO
153420.23+413007.6	0.92	1.36	QSO	QSO
155633.78+351757.3†	1.49	3.24	QSO	QSO
155921.70+435700.6	0.66	1.36	QSO	QSO
160113.56+422637.7	0.15	1.72	QSO	QSO
163148.98+495933.4	4.53	1.45	QSO	QSO
170016.63+350353.5	4.50	1.39	QSO	QSO
142843.33+120913.3	1.19	1.19	QSO	GALAXY
143827.41+370737.3	1.01	1.01	QSO	QSO (Type 2)
130100.88+320727.4	0.51	0.51	QSO	QSO (Type 2)

Matching our SP sample to the SDSS DR16 QSO catalogue [50] we found 117 matches, of which 110 are considered high-quality from the VI (107 are included in our SP QSO sample; see § 2.1). For these 110 matches, we compared our VI redshifts and spectypes to that from the SDSS DR16 catalogue. Overall, seven sources have significantly different redshifts in SDSS to those obtained by the DESI VI ($\Delta z=|z_{\rm VI}-z_{\rm SDSS}|/(1+z_{\rm VI})\geq0.05$; $\Delta v>15\,000$ km s$^{-1}$), two sources were classified as Type 2 AGN in the VI and one additional source was classified as a galaxy in the VI, although all three had correct redshifts (see Table 4). For these ten sources, we performed additional VI to confirm the redshift and spectype; the seven sources with redshift disagreements had clearly incorrect SDSS redshifts (including one extreme high redshift BALQSO; see Table 4). The sources that were VI’d as a galaxy or Type 2 AGN have weak broad emission-line features.

Testing the radio detection fraction plot

The key result of this paper is the striking positive correlation between the amount of line-of-sight dust extinction and the radio detection fraction (Figure 10). To check the robustness of this result we performed three tests: 1) varying the MAD cut we imposed to define a “good” spectral fit, 2) changing the redshift cut, and 3) exploring the 1.4 GHz FIRST and 3 GHz VLASS detection fractions (see Figure 17).

For the first test, we applied three different MAD cuts: $4\times10^{-5}$ (used for Figure 10), $7\times10^{-5}$, and $20\times10^{-5}$, which were satisfied by $\sim$ 58 per cent (20 031/34 293), $\sim$ 85 per cent (29 001/34 293), and $\sim$ 99 per cent (33 913/34 293) of the combined QSOs, respectively. We find that relaxing the MAD cut slightly decreases the overall radio detection fraction, but does not change the shape of the trend; this demonstrates that the choice of MAD cut does not significantly affect our basic result. For the second test, we used a redshift cut of $0.5$ $<$ $z$ $<$ $2.0$ and $1.0$ $<$ $z$ $<$ $2.5$. The first range was chosen to exclude the highest redshift sources ($z$ $>$ $2.0$) that may not have such reliable fits due to the limited rest-frame spectral coverage of our blue QSO composite used in the fitting and the second range excludes the lowest redshift sources ($z$ $<$ $1.0$) that may still have contribution from the host galaxy (§ 3.3). We find that the overall trend is the same for both redshift ranges.

Radio survey	Radio detection fraction
LDR2	6549/34 293 (19%)
FIRST	1001/18 816 (5%)
VLASS	631/34 293 (2%)

For the last test, we matched the combined QSO sample to the 1.4 GHz VLA FIRST and 3 GHz VLASS²³ surveys, adopting a $10''$ matching radius for both surveys which resulted in 5 per cent (1001/18 816; the sample was reduced to the FIRST survey area) and 2 per cent (631/34 293) matches, respectively. This is consistent with the sensitivity of FIRST (0.15 mJy) and VLASS (0.12 mJy); assuming a canonical spectral index of $\alpha=0.7$, extrapolating the LDR2 flux values to 1.4 and 3 GHz we expect 1079 and 882 QSOs to be detected in FIRST and VLASS, respectively. We find a positive trend in the FIRST radio detection fraction as a function of $E(B-V)$ out to $E(B-V)=0.3$ mag, after which the source statistics are too limited to conclude anything meaningful (Figure 17). This is despite FIRST being $\sim$ $9$ times less sensitive than LoTSS, again adding to the robustness of our result. The VLASS radio detection fraction displays a subtle increase with $E(B-V)$, but unfortunately the source statistics are too low to robustly conclude anything. Table 5 displays the radio detection fractions of the combined QSO sample for LDR2, FIRST, and VLASS.