Fundamental differences in the properties of red and blue quasars: measuring the reddening and accretion properties with X-shooter

We have recently found fundamental differences in the radio properties of red quasars when compared to typical blue quasars. In this paper we use X-shooter data, providing spectral coverage from $\sim 3000-25000$ Ang, of a sample of 40 red and blue luminous quasars at $1.45<z<1.65$ to explore the connections between the radio, emission-line, and accretion-disc properties. We fit various dust-extinction curves to the data and find that dust reddening can fully explain the observed colours for the majority of the red quasars in our sample, with moderate extinctions ranging from Av$\sim 0.06-0.7$ mags. We confront our spectra with a simple thin accretion-disc model and find this can describe the continua of both the blue and red quasars, once corrected for dust extinction; we also find no significant differences in the accretion properties. We detect ionized outflows in a number of red and blue quasars, but do not find any significant evidence that they are more prevalent in the red quasar population. Overall our findings imply that the radio emission is more closely connected to circumnuclear/ISM opacity rather than accretion disc or outflow differences.


INTRODUCTION
Quasi-stellar objects (QSOs), also known as quasars, are the most powerful class of Active Galactic Nuclei (AGN), with extremely high bolometric luminosities (up to 10 47−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 8 -10 9 M ). The majority of quasars are blue due to an unobscured view of the SMBH accretion disc which peaks in the ultra-violet (UV). However, there is a significant subset of QSOs with redder optical-infrared colours (referred to as "red QSOs"). Although red QSOs have been well studied in the literature (Webster et al. 1995;Serjeant & Rawlings 1996;Kim & Elvis 1999;Richards et al. 2003;Glikman et al. 2007;Georgakakis et al. 2009;Urrutia et al. 2009;Banerji et al. 2012;Glikman et al. 2012;Kim & Im 2018;Klindt et al. 2019;Fawcett et al. 2020;Rosario et al. 2020;Calistro Rivera et al. 2021;Rosario et al. 2021), the origin of the red colours is still debated. The most popular explanation is reddening due to dust extinction, which attenuates the blue light from the accretion disc (e.g., Webster et al. 1995;Glikman et al. 2007; Klindt et al. 2019). However, a red synchrotron component or stellar contamination from the host galaxy have also ★ E-mail: victoria.fawcett@durham.ac.uk been shown to contribute significantly to the colours of QSOs in strongly radio-loud/gamma-ray-loud and lower luminosity systems, respectively (e.g, Whiting et al. 2001;Glikman et al. 2007;Shen & Liu 2012;Klindt et al. 2019;Calistro Rivera et al. 2021).
To study reddening due to dust in QSOs, the two standard methods used are through fitting the continuum with extinction curves and measuring the Balmer decrement (the Hα to Hβ broad line flux ratio). For extinction-curve analyses, the extinction laws of the Small and Large Magellanic Clouds (SMC and LMC; e.g., Prevot et al. 1984;Clayton & Martin 1985) are used preferentially over a Milky Way standard extinction curve (MW; e.g., Cardelli et al. 1989). The SMC/LMC laws have a weak or absent 2175 Å graphite feature which, while prominent in the MW, has only been detected in a minority of QSOs (Jiang et al. 2011;Shi et al. 2020). Past studies have found that an SMC-like extinction curve best describes the dust extinction in red QSOs (Richards et al. 2003;Hopkins et al. 2004;Young et al. 2008;Urrutia et al. 2009). However, these analyses are usually limited to the SMC, LMC and MW-like curves, and based only on optical photometry or spectroscopy over a narrow wavelength range. The Balmer decrement provides an alternative constraint on the extinction towards the QSO through measuring the strength of different broad hydrogen emission lines (i.e., the Balmer and Paschen series). In previous studies, clear evidence has been found for larger Balmer decrements, and therefore larger amounts of dust extinction, in red QSOs compared to typical QSOs (Glikman et al. 2007;Kim & Im 2018).
However, even amongst the studies that favour dust as the cause of the reddening, there is further debate on the location and nature of this intervening dust. For example, some studies (Wilkes et al. 2002;Rose et al. 2013) postulate that the dust is associated with a grazing view through the putative AGN "dusty torus" (Antonucci 1993;Urry & Padovani 1995;Netzer 2015;Almeida & Ricci 2017). An alternative model (often referred to as the evolutionary scenario; Sanders et al. 1988) identifies red QSOs as a short-lived phase in which an obscured QSO rapidly reveals itself by blowing away the surrounding shroud of dust in the circumnuclear regions and interstellar medium (ISM) through the action of powerful winds (Hopkins et al. 2006(Hopkins et al. , 2008Farrah et al. 2012;Glikman et al. 2012;Alexander & Hickox 2012).
Recently, our group have found fundamental differences in the radio properties of red QSOs when compared to typical QSOs that cannot be explained just by the orientation of the dusty torus (Klindt et al. 2019;Rosario et al. 2020;Fawcett et al. 2020Fawcett et al. , 2021. The key results from these analyses is the identification of excess compact radio emission in red QSOs, which arises around the radio-quietradio-loud threshold (defined as the ratio of the 1.4 GHz luminosity to 6μm luminosity; see Klindt et al. 2019 for details). In Fawcett et al. (2020), we found that this radio excess decreases towards the extreme radio-quiet end and appears to be driven by AGN mechanisms rather than star-formation. In a more recent study, we have showed that this excess is compact even on extremely small scales (typically below 2-3 kpc; Rosario et al. 2021). These radio properties could be connected to the accretion engine, depending on parameters such as black-hole mass, accretion rate and black-hole spin, or connected to larger-scale processes such as winds or outflows (i.e., the commonly referred to "AGN feedback process"). However, star-formation cannot be conclusively ruled out for the low redshift/lower radio luminosity sources.
To test whether differences in the accretion properties between red and typical QSOs are driving the enhanced radio emission, accurate black-hole masses are required which can be then used for a characterization of the pure accretion disc (AD) emission through continuum fits (e.g., Capellupo et al. 2015). There are many methods to constrain black-hole mass, the most commonly used being the "single epoch virial mass determination" method. These methods combine a broad line full width at half-maximum (FWHM) with either a continuum or broad line luminosity to obtain the black-hole mass, building on calibrations between AD luminosity and the broad line region (BLR) size, derived from reverberation mapping studies (e.g., Kaspi et al. 2000;Bentz et al. 2013). In order to self-consistently constrain key accretion parameters, single-epoch spectroscopy with wide wavelength coverage is required: X-shooter (Vernet et al. 2011) at the Very Large Telescope (VLT) offers this capability, providing continuous spectral coverage across observed wavelengths ∼ 3000-25000 Å. In a series of papers (Capellupo et al. 2015;hereafter C15, Mejía-Restrepo et al. 2016;hereafter M16 and Capellupo et al. 2016; hereafter C16), X-shooter spectra were found to be very effective at constraining the accretion and emission-line properties in a sample of ∼ 1.5 QSOs. For example in C15, a thin AD model was found to provide a good fit to the majority of their QSO X-shooter spectra. The QSOs that could not be fit with an AD model were found to either have modest amounts of dust ( ∼ 0.1-0.45 mags) or required an additional disc-wind component which, once implemented, provided a good fit to the majority of the remaining spectra.
In this work we use new and archival X-shooter UV-NIR spectra for a sample of 40 red and typical QSOs at 1.45 < < 1.65, in order to determine whether the cause of the red colours is consis- The optical * − * colour (top) and L 6 μm (bottom) versus redshift distributions for our full rQSO (red stars) and cQSO (cyan stars) samples; the cQSOs included from C15 are outlined in black. The three QSOs that are excluded from our main analysis are shown as grey squares and the Fermidetected source is annotated. The DR7 QSOs from our parent sample are shown in grey. In the top plot the lower and upper * − * boundaries for our rQSO and cQSO selections, calculated from the DR7 QSO parent, sample are shown by the dashed red and cyan curves, respectively (see Section 2.1). tent with arising entirely from dust, or whether differences in the accretion properties could be driving both the colours and enhanced radio emission in red QSOs. In Section 2 we define our sample, describe the observations, the reduction process, and the modelling of the telluric absorption. In Sections 3, 4, and 5 we describe our dust extinction, emission line, and AD model fitting procedures, respectively. In Section 6.1 we quantify the amount of reddening present in our red QSOs through an extinction-curve comparison and analysis of the Balmer decrements, in Section 6.2 we explore the emissionline properties of our sample, and in Section 6.3 we explore whether there are any differences in the AD properties between red and typical QSOs. In Sections 7.1 and 7.2 we discuss the nature of dust in red QSOs and the presence of accretion-driven winds in red QSOs, respectively. Throughout our work we adopt a flat Λ-cosmology with 0 = 70 km s −1 Mpc −1 , Ω M = 0.3 and Ω Λ = 0.7.

Sample selection
Our sample is based on the Sloan Digital Sky Survey (SDSS) DR7 Quasar Catalogue (Schneider et al. 2010)  which restricts the sample to the luminous end of the QSO population that satisfy the main QSO selection algorithm (Richards et al. 2002), and removing the Faint Images of the Radio Sky at Twenty Centimeters (FIRST; Becker et al. 1995) pre-selected sources. We will hereafter refer to sources using a shorthand version of the SDSS identifier; e.g., SDSSJ160808.12+120702.4 is shortened to 1608+1207. We require mid-infrared (MIR) counterparts from the Wide-Field Infrared Survey Explorer (WISE; Wright et al. 2010), an all-sky survey which provides photometry in four bands (3.4, 4.6, 12 and 22 μm). It has been shown that MIR emission, arising from the hot dust heated by high-energy radiation from the accretion disc, is a useful, extinction-insensitive discriminant of QSOs (e.g., Stern et al. 2005;Lacy et al. 2005;Stern et al. 2012;Assef et al. 2013), as well as a measure of their intrinsic luminosity; for example the commonly used rest-frame 6 μm luminosity (L 6 μm ). We used the NASA/IPAC query engine to match SDSS DR7 QSOs to the All-Sky WISE Source Catalogue (ALL-WISE) adopting a 2. 7 search radius, which ensures a 95.5 per cent positional certainty (Lake et al. 2011), and required a detection with a signal-to-noise ratio (SNR) of greater than 2 in the WISE W1, W2 and W3 bands in order to derive an accurate estimation of L 6 μm . To compute the L 6 μm we used a log-linear interpolation or extrapolation of the fluxes in the W2 and W3 bands.
In this paper, we aim to extend our understanding of the fundamental differences between red and typical QSOs. To distinguish between these two populations, we follow the approach adopted in our previous papers (e.g., Klindt et al. 2019). We select the red QSO sample as the top 10 per cent of the redshift-dependent * − * distribution. Importantly, we then select a "control" sample in a consistent way to the red QSOs, as the median 50 per cent of the * − * distribution, in order to identify a set of typical QSOs that we can reliably compare to within our analyses 1 (see Klindt et al. 2019 for more details). From this colour-selected sample, we then targeted a sub-sample of red QSOs (rQSOs) and control QSOs (cQSOs), matched in redshift and 6 μm to be observed with VLT/X-shooter. We chose our redshift range and observing strategy based on the approach taken in C15 (Program ID 088.B-1034). C15 used X-shooter to observe 30 typical QSOs in a redshift range of 1.45 < < 1.65, which was chosen to ensure that the H α and the H β- [O ] complex lay in the H-band and J-band, respectively. In addition to our own dedicated observations (see Section 2.2), we included 15 QSOs from the C15 sample that were matched in 6 μm and * − * colour to our cQSO sample, and also satisfy our uniform flag and WISE constraints; these cQSOs from C15 are flagged in Table 1. The cQSOs from the C15 sample are also indicated in Fig. 1, but hereafter will be referred to as part of the full cQSO sample.
Due to scheduling restrictions, only 54 per cent of our dedicated observations were completed: the observed targets were left with a L 6 μm bias (see Fig. 1). In the relevant analyses (Section 6.3) we repeat our comparisons using a statistically-limited luminositymatched sub-sample in order to account for this bias. For the matching, we used an L 6 μm tolerance of 0.2 dex, following a similar approach to Klindt et al. (2019) (see Section 2.3.2 therein). Within our sample we found that a mistake in the amount of Galactic extinction was originally applied to two rQSOs (2109+0045 and 2239+1222) which, when corrected, moved the sources out of our * − * rQSO threshold. We do not consider these two sources in our main analyses and they will be referred to as "QSO" within the tables to distinguish them from the cQSOs and rQSOs.
The * − * colours and L 6 μm versus redshift for our full sample are shown in Fig. 1 (top and bottom panel, respectively). The radio-detected sources, and upper limits for the radio-undetected sources, are displayed in Fig. 2; the dashed line indicates the radioquiet threshold, defined here as the mechanical-to-radiative power of 0.1 per cent. This corresponds to R = −4.2, where R is the "radioloudness" parameter defined as the dimensionless quantity: We calculated L 1.4 GHz assuming a uniform radio spectral slope of α = −0.7, following a similar approach to Alexander et al. (2003). 2 This is different from the canonical definition for radio-loudness, which is defined as a ratio of the 5 GHz radio luminosity to 4400 Å optical luminosity (Kellermann et al. 1994); we do not use this definition since optical luminosities are heavily affected by dust extinction for the rQSOs. Our defined radio-quiet threshold of R = −4.2 is roughly consistent to the classical definition (see Klindt et al. 2019 for more details). In both figures, the QSOs excluded from further analyses are displayed as grey squares. Within our dedicated observed sample, we chose to include both FIRST radio-detected and undetected rQSOs and cQSOs in order to test whether the accretion and extinction properties are related to the radio properties. Our sample includes a flat-spectrum rQSO (2229-0832) which is also the only source in our sample that is detected in the fourth Fermi Large Area Telescope catalog (4FGL; Abdollahi et al. 2020). It is extremely radio-bright ( 1.4 GHz,int = 1.03 Jy) and has the largest radio luminosity in our sample; see Fig. 2 and Table 1. The Fermi detection suggests that the radio emission is relativistically beamed; also indicated by its flat spectrum. Our spectral analysis of 2229-0832 showed that the UV-near-IR emission is likely dominated by non-thermal processes (see Appendix A for a more detailed discussion), which we also confirm from our radiooptical analysis, following the method from Klindt et al. 2019 (see Appendix A therein). Due to the extreme nature of this source, we have excluded it from our analyses (indicated as 'Fermi-detected' in Figs. 1 and 2).

Observations
Our main analyses are based on observations obtained with the Xshooter instrument on the VLT (Vernet et al. 2011), which provides broad spectral coverage from ∼ 3000 to 25000 Å by observing three wavelength ranges simultaneously: the UV-blue (UVB), visible (VIS), and near-infrared (NIR). At the redshift of our sample, this corresponds to rest wavelengths of ∼ 1200 − 9800 Å.
Dedicated observations for this project were carried out in nodding mode (nodding length 5 ) between April-August 2018 (Program ID 0101.B-0739(A); PI: Klindt). We adopted slit widths of 1. 6, 1. 5 and 1. 2 in the UVB, VIS and NIR arms, respectively, corresponding to a resolving power of 3300, 5400 and 4300. For the brighter sources in our sample (J-band magnitude < 17.6 mags AB) we used 900 sec exposures, and for the fainter targets we used 1200-1400 sec exposures. For the majority of the sources, observations were taken under conditions where the seeing was ≤ 1. 6 and airmass ≤ 1.5. However, for two objects the seeing was highly variable throughout the night, reaching 1. 8 − 1. 9 at times: these QSOs are identified in Table 1. From visually inspecting the data for these two sources, we found no issues with the spectra that could arise from poor seeing during observations. Telluric standard stars were not provided by the observatory following the specification change from P101: instead, the ESO software tool molecfit was used to fit and correct telluric lines (see Section 2.3). To obtain the flux-calibrated spectra, spectrophotometric standards were taken once a night.
The full sample, including both our dedicated observations and the sub-sample from C15, consists of 12 rQSOs (6 radio-detected) and 28 cQSOs (12 radio-detected). We provide an overview of the observations for the full sample (including the three excluded QSOs) in Table 1.

Data reduction
We reduced the X-shooter data from our dedicated program using version 2.9.1 of the ESO X-shooter pipeline (Modigliani et al. 2010) within the ESO Reflex environment (Freudling et al. 2013). Data for each spectrographic arm were reduced individually, adopting the nodding mode setting in the pipeline. This provides an automated workflow, in which detector bias and dark current are subtracted followed by a wavelength and flux calibration, to produce the fluxcalibrated one-dimensional spectra. In nodding mode, a sequence of science frames at different nodding positions are combined into an image sequence, and then the pairs of combined frames are subtracted. Sky subtraction is carried out by taking the difference of the images from nodding position A and B. Due to a gradient in one of the four raw VIS images for the cQSO 1552+0939, we were unable to effectively subtract the sky from the overall source spectrum. For this source we therefore carried out the reduction in the VIS arm excluding the pair of raw frames corresponding to this gradient.
In the majority of our sources there is good agreement in the flux for the overlapping regions of the spectra from the three arms. However, in one case (cQSO 1513+1011) there was a clear offset between the visible and NIR spectra output by the pipeline due to an issue with the flux calibration of the standard star. For this source we used the publicly available reduced data for the NIR arm from the ESO archive, which we found to be consistent with the overlapping region of our VIS reduced spectrum.
To further test the validity of our flux calibration, we calculated synthetic , , and -band photometry by applying the SDSS photometric filters to our X-shooter spectra, and compared these data to both the SDSS DR7 and Pan-STARRS1 (PS1; Chambers et al. 2016) photometry (adopting a similar method to Selsing et al. 2016). We found that there was a small amount of scatter between all three samples, with an average standard deviation in the photometric offsets of ∼ 0.24 mag; consistent with the predicted variation due to QSO variability (0.26 mags;MacLeod et al. 2012). We did find a small mean offset of ∼ −0.15 mag between our synthetic photometry compared to both SDSS and PS1 which is not wavelength dependent. This suggests that some flux is lost from the X-shooter spectra which is likely due to extended emission beyond the QSO point source.
Ground-based observations suffer from atmospheric absorption (predominantly by H 2 O, O 2 , CO 2 and CH 4 ), which is particularly prominent in the NIR waveband. After running the X-shooter pipeline, we corrected the spectra for telluric absorption using a set of feature-free continuum windows (see Table 2; some windows were expanded for a few low SNR sources). This was achieved through the ESO software tool molecfit 3 (Smette, A. et al. 2015;Kausch, W. et al. 2015), which is based on synthetic modelling of the Earth's atmosphere. The model spectra are generated by the radiative transfer code Line-by-Line Radiative Transfer Model (LBLRTM), which takes an atmospheric profile created using local atmospheric information modelled from Global Data Assimilation System (GDAS 4 ) data and meteorological measurements from ESO Meteo Monitor (EMM 5 ). For all sources, the correction in the NIR arm did not remove all features, and so the regions heavily affected by absorption were masked in our final analyses.
For the C15 sources, we used reduced and telluric-corrected spectra obtained via private communication with the authors. We rebinned all the spectra to a consistent wavelength sampling and masked the same spectral regions in all of the QSOs.
Finally, we corrected all of the spectra for Galactic extinction using the Schlegel et al. (1998) map and the Fitzpatrick (1999) Milky Way extinction law. One cQSO (1521-0156) had an extreme level of Galactic extinction ( = 0.54 mags) and so for this source we masked the region around the 2175 Å feature, in order to minimize the uncertainties this introduced to our extinction and spectral fitting analyses (see Sections 3 and 4, respectively). The final spectra for the rQSOs, cQSOs and QSOs are shown in Figs. E1, E2, and E3, respectively.

DUST-EXTINCTION MEASUREMENTS
As previously mentioned, one rQSO (2229-0832) is Fermi-detected and has the largest radio luminosity in our sample; due to the extreme nature of the source, we exclude it from our analyses. In Appendix A, we conclude that this source is most likely reddened by a synchrotron component. Exploring whether the reddening for any of our other rQ-SOs can be explained with a synchrotron component, we find this to be the case only for 2229-0832. Moreover, in Fawcett et al. (in prep) we explore DR14 rQSO composites in bins of redshift and 6 μm , and find that sources with > 1 and log 6 μm > 45.0 erg s −1 (consistent with our X-shooter sample) display little to no host-galaxy absorption features (e.g., the Calcium H+K lines), and even the composite Table 1. Observation and sample details for our rQSO, cQSO and excluded QSO samples, including the 15 cQSOs from C15. The columns from left to right display the: (1) shortened source name utilized in this paper, (2) QSO sub-sample (rQSO, cQSO, QSO); the three sources labelled as QSO are excluded from our main analyses (see Section 2.1), (3) date of observations (those prior to 2013 are from Capellupo et al. (2015), while those post 2017 are new dedicated observations presented here), (4)-(5) RA and Dec as defined in the optical SDSS DR7 catalogue, (6) redshift obtained from SDSS DR7, (7) redshift defined via visual inspection, based on the observed wavelength of H α, (8) the optical SDSS * − * colour, (9) rest-frame 1.4 GHz FIRST radio luminosity (upper limits indicated by *) calculated assuming a radio spectral slope of α = −0.7, and (10) rest-frame 6 μm luminosity calculated using a log-linear interpolation or extrapolation of the fluxes in the W2 and W3 bands.
†QSOs from C15; C absorption identified; Offset between VIS and NIR; Fermi-detected; Poor seeing. from the lowest redshift and luminosity bin, in which we see strong host-galaxy features, requires an additional dust-reddening component in order to explain the continuum emission (consistent with the broad-band spectral energy distribution (SED) analysis from Calistro Rivera et al. 2021). Consequently, we expect the red colours of our X-shooter rQSOs to be predominantly due to dust extinction, as suggested by previous literature (e.g., Glikman et al. 2007;Klindt et al. 2019). Under this assumption, we predict the loss of flux at the rest-frame UV-VIS end of the spectrum to be frequency dependent, as seen in QSO dust-extinction curves. Modelling the dust extinc-tion in QSOs via the fitting of extinction curves has been undertaken in many previous studies. However, these analyses are usually only based on optical photometry or spectroscopy over a narrow wavelength range, which will not necessarily provide reliable measurements; the wide wavelength coverage of X-shooter can be used to robustly determine whether the red colours in the rQSO population are consistent with dust extinction. We adopt two methods widely utilized in the literature: fitting extinction curves to the continuum   (Richards et al. 2003; yellow), a MW law (Cardelli et al. 1989), a simple power-law (PL; purple), a flat "grey" law (Czerny et al. 2004; green) and a steeper law (Zafar et al. 2015;black). The different extinction curves correspond to variations in the dust grain composition and size. Only the MW curve displays the characteristic 2175 Å bump produced by graphite. A comparison of how these different extinction laws affect the QSO continuum is shown in Fig. D2. emission and measuring the Balmer decrement using the broad H α and H β emission-line fluxes. 6 We explored five different forms of extinction laws in this work: an SMC law (Richards et al. 2003; = 2.74) 7 , a MW law (Cardelli et al. 1989; = 3.1), a simple power-law (PL; ∝ λ −1 , = 4.0 chosen to give consistent values compared to the other curves), a flat "grey" law (Czerny et al. 2004; ∼ 2.0) and a steeper law (Zafar et al. 2015; = 2.2); a comparison of these extinction curves is displayed in Fig. 3. Only the MW-type curve displays the characteristic 2175 Å feature expected from small carbonaceous grains 6 The narrow lines are also sensitive to reddening but on larger scales than the nucleus. We use the broad lines to assess the reddening along the line-of-sight to the accretion disc and BLR. 7 ≡ / ( − ): Cardelli et al. (1989) present in our Galaxy. A steep extinction curve implies that the dust is predominately composed of smaller grains, which could be intrinsic or may be due to the destruction of larger grains by the intense AGN radiation field (Zafar et al. 2015). A flatter extinction curve therefore implies the dust is composed of larger grains that either formed far our from the central region, or are produced in an outflowing wind that escapes the central region quick enough to avoid the destruction of larger grains (Czerny et al. 2004;Gallerani et al. 2010).
To quantify the amount of dust extinction present in each source, we constructed a composite spectrum of our cQSOs which was then fit to both the rQSOs and cQSOs using a least-squares minimization code, which varied the normalization (fixed at the NIR end) and ( − ) parameters for each of the five different extinction curves. To create the composites, the Galactic extinction-corrected spectra were first shifted to rest-frame wavelengths (using VI ; see Section 4.1) and then adjusted to a common wavelength grid. The composite was then created by taking the geometric mean across all spectra; the geometric mean was used rather than the median or arithmetic mean since it preserves the spectral shape, extinction law, and mean extinction of the composite (see Appendix in Reichard et al. 2003). A comparison of our cQSO and rQSO composites is shown in Fig. 4 and a more detailed comparison of our cQSO composite to other QSO templates from the literature can be found in Appendix D. For the extinction-curve analysis, we then masked the emission-line regions and smoothed the spectrum of the cQSO composite with a Gaussian filter (σ = 3) to produce our unreddened QSO model. It should be noted that some of the cQSOs will also contain modest amounts of extinction, which will therefore be incorporated into the cQSO composite. However, we expect the extinction towards the cQSOs to be modest; indeed, splitting the cQSO composite into the bluest and reddest 50 per cent of sources, we estimate a spread in extinction values of ( − ) ∼ 0.013 mags ( ∼ 0.04 mags; consistent to the results found from fitting broad-band SEDs for statistical samples in Calistro Rivera et al. 2021).
The Balmer decrement provides a constraint on the amount of dust extinction along our line-of-sight towards the BLR, compared to the extinction-curve analysis which constrains the amount of dust extinction in the continuum region. However, due to the expected proximity of the accretion disc and BLR for the standard QSO model, to first order we would expect very similar levels of line-of-sight extinction. To calculate the dust extinction from the Balmer decrement, we used the equation from Calzetti et al. (1994): where k(λ) is the extinction curve evaluated at λ: in our calculation, we used a simple PL extinction curve which gives k(H α)= 3.35 and k(H β)= 4.53. The observed Balmer decrement [H α/H β] measured represents the broad H α to H β flux ratio, obtained from the emission-line fitting (see Section 4). The intrinsic Balmer decrement [H α/H β] intrinsic is usually set to the theoretical "Case B" value of 2.88 (Osterbrock 1989). However, this value is not representative of the typical broad line Balmer decrement in unreddened QSOs; for example, the Balmer decrement of our cQSO composite is 3.40 (consistent with the average or composite Balmer decrement from Dong et al. 2007;Glikman et al. 2007;Shen et al. 2011). In order for our Balmer reddening measurements to be consistent with our continuum-based reddening estimates (see Section 6.1.1), we adopted the Balmer decrement of our cQSO composite as the intrinsic ratio, taking a similar approach to Section 5.2 within Glikman et al. (2007). From Eq. 2, a BLR dust extinction of ( − ) = 0.2 mags, 2 × 10 3 3 × 10 3 4 × 10 3 6 × 10 3 Rest Wavelength (Å)  for the same measured H α flux, would correspond to a drop in the H β flux of around ∼ 20 per cent.

SPECTRAL FITTING
To characterize the continuum and emission-line properties of our QSO sample, we fitted the X-shooter spectra using the publicly available multi-component fitting code P QSOF 8 (Guo et al. 2018a). A detailed description of the code can be found in Guo et al. (2018b) and Shen et al. (2019). In this section we describe our approach to fitting the continuum emission (Section 4.1) and the emission lines (Section 4.2).

Continuum fitting: global approach
For each source, we first smoothed the spectra with a 10 pixel boxcar 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) and Fe continuum components, providing our own redshifts based on the H α broad line. Since all our sources are at a redshift > 1.16, P QSOF does not include a host-galaxy component; our sample is also luminous, and so we expect the effect of the host galaxy on the QSO spectra to be minimal (Calistro Rivera et al. 2021;Fawcett et al. in prep; see Section 3). The power-law continuum ( PL ) is defined as where λ 0 = 3000 Å is the reference wavelength. The parameters 0 and 1 are the normalization and power-law slope, respectively. The BC ( BC ; Dietrich et al. 2002) is defined as: where BE is the normalized flux density, λ ( ) is the Planck 8 https://github.com/legolason/PyQSOFit function at the electron temperature and λ is the optical depth at the Balmer edge λ BE = 3646 Å. Here BE , and λ are free parameters.
The Fe model ( Fe ) is defined as: For the QSOs with a measured dust extinction of ( − ) > 0.03 mags, 9 we also include a third-order polynomial to account for the intrinsic dust extinction, defined as: where are the polynomial coefficients. We set 'rej_abs = True', which removes any outliers falling below 3 of the continuum for wavelengths < 3500 Å to reduce the effect of absorption lines biasing the fitting.
From fitting the spectra, we found small disagreements between the SDSS redshift and the redshift inferred based on the observed wavelength of H α. The corrections were modest (average Δ ∼ 0.004) with the most extreme correction of Δ = 0.0275 applied to cQSO 1357-0307. Both the SDSS redshift ( SDSS ) and amended redshift ( VI ) are shown in Table 1. Throughout the rest of the paper we use the VI values.

Emission-line fitting: local continuum approach
To fit the emission lines we first subtracted from the spectra the Fe and BC models based on the global fit, before re-fitting the emission lines using small, localised power-law continuum windows (see shaded region in Fig. 4). The emission lines were modelled using Gaussian components: 1 Gaussian for the narrow lines and 3 Gaussians (1 narrow and 2 broad) for the broad lines. A local fitting approach for the emission lines was adopted, rather than a global approach, in order to reduce the systematic residuals from the telluric correction which can affect the performance of the global continuum fits. In M16, they compared both approaches and found that using a local continuum leads to a very small, but systematic, overestimation of the continuum emission, which will consequently lead to an underestimation of the FWHM. Depending on the emission line, we adopted additional constraints in the fitting process as described below.
For selected emission lines ([O ]λλ4959,5007, H β, [N ]λλ6549,6585 and [S ]λλ6718,6732) we tied either the flux ratio, line width or both (following a similar approach to Shang et al. 2007; see Table 4 therein). The two [O ]λλ4959,5007 emission lines were first fitted using a core and a wing component for each line; the fluxes of the core components were tied to the theoretical 1:3 ratio, respectively, and the FWHM of the wing components were limited to < 2000 km s −1 (similar to that in Calistro Rivera et al. 2021). The widths for all of the other narrow emission lines were set to the fitted [O ]λ5007 FWHM, if available, or were otherwise limited to < 1200 km s −1 . The maximum velocity offset for the narrow lines from the systemic redshift is also limited to 500 km s −1 . Within the H α complex, the flux of the [S ]λλ6718,6732 lines were fixed to a 1:1 ratio, and the flux ratio of the [N ]λ6549 and [N ]λ6585 lines were set at 1:3, respectively. The H α line was freely fit with three Gaussians to represent the broad and narrow emission-line components. The H β emission line was then fitted using H α as a template, fixing the peak and widths of the two broad components to be the same as that of H α. This approach was required due to the low signal-to-noise of the H β line in ∼ 28 per cent of the sources. For the broad emission lines, we limited the FWHM range between 1000 and 20000 km s −1 . We also allowed velocity offsets of up to 1000 km s −1 for all broad lines apart from C , for which we allowed offsets of up to 4000 km s −1 to accommodate for blueshifts associated with potential outflows (Rankine et al. 2020). Our chosen offsets and line width limits are based on previous literature (Trakhtenbrot Guo et al. 2018b). As a further check, after fitting the lines we visually inspected each spectral fit to verify that our approach provided a good characterization of the emission lines and did not miss any components with broader emission or larger velocity offsets.
The uncertainties in the spectral quantities were estimated using a Monte Carlo approach within the least squares fitting code kmpfit, using 50 iterations. An example of our continuum and emission-line fitting is show in Fig. 5.

ACCRETION-DISC FITTING AND BLACK-HOLE MASS MEASUREMENTS
To assess whether rQSOs and cQSOs have different accretion properties (e.g., as may be predicted if the rQSOs are an earlier black-hole growth phase), we compared the X-shooter spectra to standard AD models. For this purpose, we adopted a geometrically thin, radiative AD model from Slone & Netzer (2012), and fit these models to the individual spectra of both the rQSOs and cQSOs, correcting for dust extinction. Three parameters completely describe these models: black-hole mass (6.3 < M BH < 10.3 log[M ]), black-hole spin (−1.0 < < 0.998) and mass accretion rate (−3.0 < < 1.0 log M yr −1 ); these limits are based on standard QSO models and are similar to that used in C15.
To reduce the number of free parameters in the model, we can constrain M BH using the viral relation from single-epoch spectra (as done in e.g., C15). In our work we use the FWHM of the broad H α line to calculate M BH , since this is less affected by dust extinction than the more widely used Mg line, consistent with previous red QSO black-hole mass estimates. 10 We adopted the M BH calibration from Shen & Liu (2012): where = 1.390, = 0.555 and = 1.873 for FWHM H α and 5100 . We also estimated M BH from the Mg emission line for completeness and to provide a consistent comparison to the M BH values estimated by C16, using = 1.816, = 0.584 and = 1.712 for FWHM Mg and 3000 .
To fit the AD model and determine the best-fitting parameters we used the Python χ 2 minimization code lmfit 11 , adopting a leastsquares algorithm. We set M BH to the values obtained from this approximation and allowed this to vary within the errors. We ran the fitting twice: first fixing the value of to the non-negative bestfitting values obtained from our dust-extinction fitting (see Sections 3 and 6.1.1), and then allowing this to vary as a free parameter in order to test whether a dust component is required in order to produce the best fits for the rQSOs. In both cases, the mass accretion rate and black-hole spin were left as free parameters.

RESULTS
Using our X-shooter spectra for 12 rQSOs and 28 cQSOs at 1.45 < < 1.65, we have fitted the continuum and emission-line components with physically motivated models (see Sections 3, 4, and 5) to constrain a range of key parameters; in the following sections we present our results.

Continuum reddening: extinction-curve comparison
Under the assumption that dust is the main cause of reddening in the rQSOs (see Section 3), we quantified the amount of dust extinction along the line-of-sight by fitting the QSOs with a dust-reddened cQSO composite, adopting five different extinction curves: SMC, 10 Since these are luminous QSOs we expect any contamination to the broad H α emission line from star-formation to be negligible.  Table 3. Table of the best-fitting dust-extinction parameters and the corresponding reduced χ 2 (χ 2 r ) obtained from fitting the SMC, MW, PL, flat and steep extinction curves to our rQSOs and cQSOs. Overall, the ( − ) values found for the five different curves are very similar. For each rQSO, the ( − ) and χ 2 r values of all extinction curves both within 0.9 < χ 2 r < 3 and 0.3 of the best-fitting χ 2 r value are highlighted in bold. It is worth noting that there will be additional uncertainties introduced by the spread of extinction values within our composite ( ( − ) ∼ 0.013 mags) that would need to be taken into account in order to convert to a fully unobscured QSO (see Section 3). The reddening estimated from the Balmer decrements are shown for comparison. Associated Mg absorption.  MW, PL, steep and flat (see Fig. 3). We also investigated whether there was a preference between the five different extinction curves, which would provide an insight into the nature of the dust; i.e., the size and location of the dust grains. We varied the normalization and amount of dust extinction applied to our cQSO composite using each of the five extinction-curve models, and then fit this to our QSO spectra using a least-squares minimization code. The bestfitting ( − ) constraints and the corresponding reduced χ 2 (χ 2 r ) values are displayed in Table 3; we note that the ( − ) values were derived with respect to the cQSO composite which will contain some amount of dust extinction, although this is generally negligible (see Section 3 and Calistro Rivera et al. (2021)). A movie that demonstrates how increasing amounts of dust extinction affects the shape of the cQSO composite in comparison to the rQSO spectra is shown in Fig. 8 of the online Supplementary material.
Overall, for any given QSO, the ( − ) values obtained from the five curves are in good agreement; the standard deviation in the maximum difference between the best-fitting ( − ) values obtained from five extinction curves is ∼ 0.014 mags. We define a fit to be good if the χ 2 r value is both within the range 0.9 < χ 2 r < 3 and within 0.3 of the χ 2 r value for the best-fitting extinction curve (the ( − ) and χ 2 r values for the good-fitting extinction curves are highlighted in bold in Table 3). By this definition, only one of the rQSOs (2241-1006) could not be well fitted using any of the curves. The poorly fitted rQSO is on the borderline of our rQSO definition and potentially has additional extinction due to an intervening absorber. The other rQSOs are well fitted by at least one of the extinction curves, which suggests that their red colours can be fully explained by dust extinction. For the rQSOs, we find extinction values in the range of ( − ) ∼ 0.02-0.23 mags ( ∼ 0.06-0.7 mags; with the additional uncertainty of ( − ) ∼ 0.013 mags introduced from Av (mag) Figure 6. Δ( * − * ) versus the ( − ) values obtained from the best fitting extinction curves for the rQSOs (red) and cQSOs (cyan) (see Table 3). Faded points indicate poor fits (χ 2 r > 3 or χ 2 r < 0.9). There is a clear relationship between optical colour and dust reddening. our cQSO composite; see Section 3), consistent with Richards et al. (2003) and Figure 4 within Klindt et al. (2019). We find that many of the rQSOs have good fits from more than one extinction curve, preventing us from identifying a clear best-fitting extinction curve. However, our analyses do allow us to rule out certain extinction curves; we find that for all but one rQSO, it is possible to robustly rule out both the MW and flat extinction curves 12 . This suggests that the characteristic 2175 Å bump is rare in reddened QSOs, in agreement with previous studies (e.g., Zafar et al. 2015;Temple et al. 2021). The sharper rise of the PL, SMC, and steep extinction curves suggests that the dust is composed of smaller grains, formed within the vicinity of the SMBH where nuclear activity destroys the larger grains (Li & Draine 2001;Zafar et al. 2015; see Section 7.1 for a more detailed discussion). However, a larger sample of rQSOs with well-constrained ( − ) measurements is required to place robust constraints on the nature of the dust grains. Due to the lack of dust (by definition), it is not possible to distinguish between different extinction curves for the cQSOs. Excluding the poorly fitted sources, we find that there is a strong correlation between Δ( * − * ) and dust extinction (see Fig. 6) which shows our method of selecting dustreddened rQSOs is robust. We use the ( − ) values obtained from the best fitted extinction curve per source for further analyses.
In our sample, 4 rQSOs and 3 cQSOs display narrow Mg absorption features within their spectra, with < indicating a potential intervening non-associated absorber. For these sources there may be an additional extinction component associated with the absorbing system along the line-of-sight, and so the intrinsic QSO extinction could be slightly lower than that measured from our extinction-curve analysis; these sources are indicated within Table 3.
Fitting a dust-reddened cQSO composite to our rQSO composite (see Fig. 4), we find that it is consistent with an average dust extinction of ( − ) = 0.106 ± 0.001 mags, again indicating that our rQSO sample is only moderately reddened. This is consistent with Calistro Rivera et al. (2021), who find a median value of ( − ) = 0.12 +0.21 −0.08 mags for a sample of ∼ 1000 SDSS red QSOs, obtained from broad-band SED fitting. Interestingly, the emissionline profiles of the dust-reddened cQSO composite are similar to that in the rQSO composite, apart from C which is more suppressed in the rQSO composite. This difference in the C profile could be intrinsic, linked to the higher incidence of both broad and narrow absorption lines in the rQSO sample (see Section 7.2 for a more detailed discussion), or a result of our small rQSO sample. The potential differences in the C emission line between rQSOs and cQSOs will be explored using a much larger sample in our future DR14 spectroscopic study (Fawcett et al. in prep).

Broad line reddening: Balmer decrements
Using the measured broad H α and H β emission-line fluxes, we calculated the Balmer decrements for each rQSO and cQSO (see Section 3). A comparison of the ( − ) values calculated from the continuum fitting with those calculated from the measured Balmer decrements is shown in Fig. 7. The heavily reddened FIRST-2MASS (F2M) QSOs from Glikman et al. (2007) are shown for comparison in addition to the expected intrinsic scatter within the DR7 QSO Balmer decrements from Shen et al. (2011). The rQSOs display, on average, larger H α/H β ratios, including a significantly larger fraction of upper limits that indicate the H β line was too weak to fit (5/12: 42 per cent of the rQSOs, compared to just 3/28: 11 per cent of the cQSOs), consistent with a larger amounts of dust extinction along the line-of-sight towards the BLR. Excluding the upper limits, the majority of the cQSOs and rQSOs lie within the intrinsic Balmer decrement scatter as a function of ( − ). However, at the higher extinction values ( ( − ) > 0.55 mags) probed by the F2M sample, there is a discrepancy between the reddening derived from the continuum and that from the Balmer decrements, which could be due to poor H β fits in noisy spectra, an intrinsically steep spectrum or a difference in the BLR-AD geometry between moderately and heavily reddened QSOs (Glikman et al. 2007;Kim & Im 2018). In future work (Fawcett et al. in prep) we will have the improved source statistics from our dedicated program using the Dark Energy Spec-troscopic Instrument (DESI; DESI Collaboration et al. 2016), to push to more dust-reddened systems than those observed with SDSS and test whether there are any differences between the continuum and BLR dust reddening at higher ( − ) values.

Continuum and emission-line properties
From our spectral fitting (see Section 4), we calculated the continuum luminosities, corrected for the continuum-measured dust extinction using the best-fitting ( − ) values, at rest-frame 1350 Å, 3000 Å, 5100 Å and 6200 Å, in addition to the emission-line properties for C , C ], Mg , H β, [O ]λ5007, and H α. The emission-line fitting parameters can be found in Table 6.2. To assess the quality of our spectral fitting, we compared the Mg FWHM and luminosity to that from the Rakshit et al. (2020) catalogue, which contains spectral properties for the SDSS DR14 QSOs, using P QSOF to estimate the continuum and major emission-line properties. We find a median offset of −1392 km s −1 and 0.25 dex for the FWHM and luminosity measurements of the rQSOs, respectively, and −941 km s −1 and 0.01 dex for the FWHM and luminosity measurements of the cQSOs, respectively. The larger difference in the luminosity measurements for the rQSOs is a consequence of the lack of dust-extinction correction applied in the Rakshit et al. (2020) catalogue; therefore, the luminosities for dust-reddened QSOs will be underestimated.
A comparison between the extinction-corrected 5100 Å luminosity ( 5100 ) and the extinction-corrected H α line luminosity and 6 μm is displayed in Fig. 8. Both the H α and 6 μm luminosities show a strong correlation with 5100 , with the best fitting power-law displayed as dashed and dotted lines. After correcting for extinction, the rQSOs have luminosities within the scatter of the cQSOs for the same 5100 . However, the rQSOs are clearly biased towards higher luminosities as a result of our incomplete observations (see Section 2.2), which we account for in the following analyses by comparing results obtained with the full sample to those obtained using the statistically limited 6 μm -matched sub-samples (see Section 2.1). Fig. 9 shows a comparison between the M BH from Rakshit et al. (2020), calculated from the Mg line, compared to the M BH estimated in this work, calculated from the Mg and H α broad line FWHM and continuum luminosities using the relations from Shen & Liu (2012) (see Section 5). In both plots, we have calculated our M BH values from dust-extinction corrected luminosities using the best-fit values (see Section 6.1); this will have a greater effect on the Mg estimation compared to H α, since shorter wavelengths are more affected by dust extinction than longer wavelengths. 13 We applied the same level of dust-extinction correction for the Rakshit et al. (2020) M BH values, to those applied for each of our sources. In Fig. 9, the squares and crosses display the median M BH calculated from our X-shooter emission line Mg (left) and H α (right) fits, with and without a dust-extinction correction, respectively. A comparison between the Mg -and H α-derived M BH from our X-shooter data is displayed in Fig. B1.
The median extinction-corrected M BH for the rQSOs and cQSOs are found to be 9.24 ± 0.17 and 9.12 ± 0.05 log M calculated from the H α line, and 9.24 ± 0.15 and 9.23 ± 0.09 log M calculated from the Mg line, respectively. The median M BH of the rQSOs and cQSOs estimated from both of these estimations are consistent within 1σ uncertainties; this is in agreement with the black-hole masses reported by Calistro Rivera et al. (2021)  values of ∼ 9 log M for both the rQSOs and cQSOs. Further to this, after matching in 6μm luminosity (M BH : 9.05 ± 0.26 and 9.12 ± 0.10 log M for the rQSOs and cQSOs, respectively) there appears to be no significant differences between the black-hole masses of rQSOs and cQSOs. It should be noted that individual estimates of M BH calculated via the virial approach will have an intrinsic scatter due to uncertainties on the inclination angle and geometry of the BLR (see Yong et al. 2016 for a more detailed discussion), and also between different estimators (see Shen & Liu 2012 for a more detailed discussion). This effect should be minimal when averaging over large samples, but due to our small samples this could be affecting the median black-hole mass estimates.

Accretion-disc properties
Following Section 5, we fitted a thin AD model to the spectra, fixing to the best-fit values from our extinction-curve analysis (see Section 6.1.1), in order to compute the mass accretion rate, blackhole spin and Eddington ratio parameters (see Fig. 10 for an example of the best-fitting AD model to a rQSO). For 83 per cent (10/12) of the rQSOs and 86 per cent (24/28) of the cQSOs, the thin AD model provided a good fit to the QSO continuum (χ 2 r < 5); the best fitting parameters and χ 2 r values are displayed in Table 5. Since both the majority of the rQSOs and cQSOs are well fitted with a thin AD model, this suggests that there are no significant differences between the accretion discs of cQSOs and rQSOs, once the effects of dust extinction are taken into account. We also refitted the thin AD models to the spectra, this time leaving as a free parameter, and found consistent levels of dust extinction for the majority of the rQSOs. For a few moderately dust-reddened sources we found AD solutions . An example of the best-fitting AD model for one rQSO. The de-reddened spectrum is shown in grey. The solid magenta curve gives the best fitted AD model solution to the original spectrum, and the dashed magenta line gives the best-fitted model to the de-reddened spectrum. Note: the upturn in the spectrum at λ < 1500 Å with respect to the AD model is due to the Ly α emission line. not requiring dust, but overall we obtained consistent AD parameters to those found adopting a fixed (see Fig. 11 for a comparison of the , with and without fixing in the AD fitting). This also suggests that dust-reddening is a consistent explanation for the red colours in the rQSOs, even if differences in accretion processes are considered. It is worth noting that we did not include an additional disc-wind component in our model, which has been found to affect the shape of QSO spectra in a similar way to dust-reddening (C15); more complex AD models will be tested in a future, larger X-shooter study. We find good agreement to C16 between our best fitting AD parameters for the C15 cQSOs (see Appendix C for a more detailed discussion).
We find no significant differences in the accretion rate or blackhole spin for our rQSOs compared to the cQSOs (see Table 6 for the median AD properties). The best-fitting black-hole spin values span the whole range (−1 to 0.998), with an average value of −0.15 ± 0.47 and 0.46±0.23 for the rQSOs and cQSOs, respectively. If differences in the post-merger spin evolution of black-holes were responsible for the radio excess in rQSOs (Garofalo & Bishop 2020), then we might expect systematic differences between the black-hole spin distributions of the two samples. However, in the luminosity-matched samples we find the two black-hole spin distributions to be consistent; applying the two-sided Kolmogorov-Smirnov (K-S) test this is statistically not significant ( -value = 0.88). The median mass accretion rates are 0.70 ± 0.07 and 0.45 ± 0.14 log [M yr −1 ] for the rQSOs and cQSOs, respectively; these differences are statistically significant ( -value = 0.09), but after matching the two samples in 6 μm the distributions are consistent ( -value = 0.42), with average values of 0.74 ± 0.17 and 0.81 ± 0.09 log [M yr −1 ] for the rQSOs and cQ-SOs, respectively. Splitting out the radio-detected QSOs within our Table 5. Table of the best-fitting AD parameters and virial black-hole mass estimates for our rQSOs and cQSOs (see Section 5). Columns from left to right display the: (1) shortened source name utilized in this paper, (2) QSO sub-sample (rQSO or cQSO), (3) virial black-hole mass (M BH ) estimate using the broad H α line (see Section 5), (4) black-hole spin ( ), (5) mass-accretion rate ( ), (6) bolometric luminosity ( bol ), (7) Eddington ratio (λ edd ), defined as Edd / bol and (8) the reduced χ 2 of the fit (χ 2 r ). The errors correspond to 1σ standard errors.  Table 6. Table of the median AD properties for our rQSO, cQSO, 6 μmmatched, radio-detected (rad det) and radio-undetected (rad undet) samples. See Table 5  samples we explored whether there were any differences in accretion properties between the radio-detected and undetected sources (see Table 6). Significant differences would support previous studies that suggest radio-loud quasars have higher black-hole spin values as per the "spin paradigm", whereby relativistic jets can be formed from the rotational energy of the spinning SMBH (Wilson & Colbert 1995;Schulze et al. 2017). However, we do not find any significant differences, although our sample is too small to draw robust conclusions. Fig. 12 shows the distribution of M BH as a function of Eddington ratio. There are small, tentative ( -value = 0.04) differences in the Eddington ratios, with rQSOs displaying on average lower Eddington ratios compared to cQSOs (median values of −0.95±0.10 and −0.68± 0.04 for the rQSOs and cQSOs, respectively); tentative differences remain after matching in 6μm ( -value = 0.03), although the sample size is now very small (7 rQSOs and 7 cQSOs). This is contrary to what was found in previous red QSO studies (e.g., Urrutia et al. 2012;Kim et al. 2015;Banerji et al. 2015;Kim & Im 2018; plotted on Fig. 12) who claim that red QSOs have higher Eddington ratios compared to typical QSOs which could be evidence for a phase  Figure 11. Best fitting mass accretion rate obtained from the AD model without fixing versus that obtained with fixing to the best-fitting extinction values from our extinction-curve analysis (see Section 6.1). Overall, we find similar mass accretion rate values obtained from both fitting methods. of rapid black-hole growth. However, the radio properties of our samples are not representative of the overall QSO population; for example, we have many more radio-loud cQSOs compared to rQSOs. It should also be noted that no difference in M BH or Eddington ratio (inferred from the FWHM of broad emission lines) between SDSS selected rQSOs and cQSOs was found in our previous studies (Klindt et al. 2019;Calistro Rivera et al. 2021), using the Shen & Liu (2012) and Rakshit et al. (2020) samples, respectively. Consequently, a larger X-shooter sample with fewer biases is required to robustly test whether there are significant differences in the Eddington ratio between rQSOs and cQSOs.

DISCUSSION
We have used the X-shooter spectra of a sample of 40 red and blue/typical QSOs at 1.45 < < 1.65, in order to constrain their line and continuum-emission, dust-extinction, and AD properties. From these data we can assess whether the enhanced radio emission we have previously found in red QSOs are due to differences in the accreting engine or outflow properties.
From fitting a dust-reddened cQSO composite to our QSO samples (see Section 6.1.1), we find that the red colours in 11/12 (∼ 92 per cent) of our rQSOs can be fully explained by dust, and rule out both the MW and flat extinction curves for all but one of the rQSOs. The amount of dust reddening in our rQSO sample is modest, with ( − ) values ranging from ∼ 0.02-0.23 mags ( ∼ 0.06-0.7 mags; see Table 3). We also find a strong correlation between Δ( * − * ) and ( − ) which suggests our method of selecting dust-reddened QSOs is robust (see Fig. 6). On the basis of our emission-line fitting approach, we also measured the broad H α/H β Balmer decrements and find on average good agreement within the intrinsic scatter for the rQSOs and cQSOs as a function of ( − ). We also find higher Balmer decrements, including upper limits, for the rQSOs compared to the cQSOs which indicates higher levels of dust along the line-of-sight towards the BLR (see Section 6.1.2).
Using the black-hole masses calculated from the broad H α emission (Section 6.2) and the dust-extinction measurements from our extinction-curve fitting (Section 6.1.1), we fitted a simple thin AD model to our rQSOs and cQSOs to test whether thin disc models adequately explain the nature of the accretion-disc emission in these QSOs (see Section 6.3). We find good fits for 83 per cent (10/12) of the rQSOs and 86 per cent (24/28) of the cQSOs, with no significant differences in the black-hole spin, mass accretion rate or black-hole mass parameters when accounting for luminosity. We find tentative evidence that rQSOs have lower Eddington ratios compared to cQ-SOs (see Table 6); this difference remains after matching in 6 μm , although we are extremely limited in source statistics (only 7 rQSOs and 7 cQSOs are matched in 6 μm ). Leaving as a free parameter we also recover consistent levels of dust extinction to our continuum measurements, and find similar accretion properties to those obtained using a fixed (see Fig. 11). This also suggests that the red colours in rQSOs are fully consistent with modest levels of dust obscuration of an otherwise normal QSO.
On the basis of our analyses using X-shooter spectra the only significant differences between rQSOs and cQSOs appears to be the presence of larger amounts of dust along the line-of-sight towards the rQSOs. One scenario that can connect this dust to the fundamental differences in the radio properties between rQSOs and cQSOs is winds and outflows interacting with a dust and gas rich ISM/circumnuclear environment. Therefore, in the following subsections we focus on constraining the location and composition of the dust and search for differences in the presence and properties of outflows between our rQSO and cQSO samples.

The nature of dust in red QSOs
The shape of a dust-extinction law, determined by the total-toselective extinction (where ≡ / ( − ); Cardelli et al. 1989), is primarily due to the composition of the dust and distribution of grain sizes. In general, a lower value of corresponds to a steeper extinction curve and smaller dust grains, and vice versa. Consequently, we can potentially learn about the physical properties of the dust by identifying the best-fitting dust-extinction laws, which are defined by different values of , for each of the QSOs. Typically, studies have used an SMC-like extinction curve ( = 2.74) to explain dust reddening in QSOs, although the effectiveness of this for describing the dust in all QSOs has been debated. For example, Czerny et al. (2004) construct a QSO extinction curve based on the blue and red composite spectra from Richards et al. (2003) and find a flatter "grey" curve that deviates from an SMC-like curve at shorter wavelengths, which has also been found in other studies (Gaskell et al. 2004;Gaskell & Benker 2007;Mehdipour & Costantini 2018;Temple et al. 2021). They found their derived grey extinction curve is best described by dust containing no graphite, with a fairly large grain size of > 0.016 μm, and suggested that this lack of graphite might support a scenario where the dust is formed in an outflowing wind (e.g., Elvis et al. 2002). On the other hand, Zafar et al. (2015) analysed high extinction QSOs for which the SMC law could not provide a good fit, deriving an average extinction curve that is steeper than the SMC law, with a best-fit value of = 2.2 ± 0.2. They suggest that the steepness of this curve implies a lack of large dust grains, which may have been destroyed by the activity of the AGN.
Across the overall rQSO sample, we do not find a significant preference for one extinction curve over another (see Section 6.1.1 and Fig. 6); however, we rule out the flat and MW extinction curve for all but one rQSO. This could imply that rQSOs have dust composed of smaller grains, described by a steeper curve, suggesting a higher level of nuclear activity or that the dust is formed within the vicinity of the central nucleus (Zafar et al. 2015). This is consistent with previous work which have found the characteristic 2175 Å bump (expected from small carbonaceous grains present in our Galaxy)  Kim & Im (2018), and the heavily reddened QSOs from Banerji et al. (2012Banerji et al. ( , 2015 (maroon squares); both of these samples favour higher Eddington ratios, albeit with very different selections. Our parent SDSS DR7 uniformly-selected cQSOs are plotted in grey. The corresponding black-hole mass and Eddington ratio distributions for our rQSOs, cQSOs and the SDSS DR7 sample are shown in the adjacent histograms as a solid red, solid cyan and dashed grey lines, respectively. The uncertainties for our rQSO and cQSO samples are calculated using MCMC. The error bars for the F2M QSOs are taken from Urrutia et al. (2012) and correspond to 0.1 dex estimates. Banerji et al. (2015) do not provide uncertainties for the Eddington ratios and so we assume a 0.1 dex estimate. absent in rQSOs (Richards et al. 2003;Czerny et al. 2004;Zafar et al. 2015;Temple et al. 2021). This result is also in agreement with Willott (2005), who claim that flatter extinction curves are driven by redshift selection effects incorporated in composites that are constructed from large QSO samples. Therefore, distinguishing between different extinction curves requires the analysis of individual spectra or composites constructed over a narrow redshift range.
We can also potentially gain insight into the dust composition from empirically determining the value of in each spectrum. To do this, we normalize the spectra relative to the observed-frame flux of the QSOs in the WISE W1 band (3.4 μm), and then calculate and ( − ) in order to determine ≡ / ( − ). The W1 band is likely dominated by thermal emission from dust grains at = 1.5 (e.g., Stern et al. 2012;Assef et al. 2013), and therefore normalizing the spectra relative to this band can allow us to estimate the intrinsic and ( − ) directly from the spectra. From our analysis we do not find a significant difference in the value between the rQSOs and cQSOs, although we do find tentative evidence for relatively low values (∼ 2) for the rQSOs. These low values are equivalent to a steeper dust-extinction curve, consistent with the results from our extinction-curve analysis. However, we are restricted to moderate levels of due to the shallowness of SDSS, as well as a small number of sources. Our future study analysing the dust-extinction properties of DESI rQSOs, which will push to much more extinguished systems than those identified in the SDSS, will robustly constrain the nature of dust in rQSOs.

Accretion-driven winds in red QSOs?
Overall, we find no significant differences in the accretion properties between rQSOs and cQSOs, and therefore we consider other mechanisms that could be driving the differences in the radio emission. One possibility is that rQSOs have more prominent winds and outflows than cQSOs, which could lead to enhanced radio emission via interactions and shocks with the ISM gas. We can search for the signatures of outflows by analysing the emission-line profiles: strong line asymmetries in the C broad line or very high velocity dispersions in the forbidden lines (e.g., [O ]λ5007) are clear indicators of powerful ionized outflows (Carniani et al. 2015;Zakamska et al. 2016;Rankine et al. 2020;Calistro Rivera et al. 2021;Jarvis et al. 2021).
In Calistro Rivera et al. (2021), they found a higher incidence of high velocity C blueshifts (> 1000 km s −1 ) in rQSOs, compared to typical QSOs, for a statistical sample of ∼ 1800 QSOs using the DR14 catalogue (Rakshit et al. 2020). From measuring the C blueshifts of our sample, we find tentative evidence for larger blueshifts in the rQSOs compared to the cQSOs; 50 per cent (6/12) of rQSOs display a C blueshift with > 1000 km s −1 , compared to 29 per cent (8/28) of the cQSOs. However, after matching in 6 μm , this difference is not significant ( -value = 0.43), and analysing a larger luminositymatched sample obtained from Rankine et al. (2020), we also do not find a difference (Fawcett et al. in prep).
The [O ]λ5007 emission line traces outflows on larger scales (outside of the BLR), compared to UV emission-line tracers. Fitting the [O ]λ5007 emission line with a core (narrow) and wing (broad) component, Calistro Rivera et al. (2021) found evidence for larger FWHMs in the wing components of a sample of SDSS red QSOs compared to typical QSOs, at < 1 with 6 μm > 44.5 ergs −1 , suggesting that outflows are more prevalent within the red QSO population. In our X-shooter sample, the [O ]λ5007 emission line is only detected in 82 per cent and 46 per cent of the rQSOs and cQSOs, respectively, and so any analyses using this line will be limited. Overall, we find evidence for broad [O ]λ5007 profiles in both the red and control sample; however, we do not find any significant differences in the FWHM of either the [O ]λ5007 core or wing components ( -values = 0.96 and 0.97 for the [O ]λ5007 core and wing FWHM, respectively). We find tentative evidence for larger blueshifts in the [O ]λ5007 wing components of rQSOs; 42 per cent (5/12) of rQSOs display an [O ]λ5007 wing with a blueshift of > 400 km s −1 compared to 21 per cent (6/28) of the cQSOs. However, matching in 6 μm this difference is not significant ( -value = 0.32), although there are too few sources for this to be conclusive. The lack of any significant differences in the [O ]λ5007 line could be due to a number of factors such as the small number of sources or low number of [O ]λ5007 detections in our sample. The [O ]λ5007 and C profiles and fits for the cQSOs and rQSOs are shown in Figs. E4, E5, E6, and E7. A future larger X-shooter sample is required to robustly tie down any differences in the outflow properties between rQSOs and cQSOs.
In Alexander et al. (in prep) we find a link between rQSOs and Low-ionization Broad Absorption Line QSOs (LoBALs; BALQSOs that display additional low-ionization species such as Mg and Al , e.g., Trump et al. 2006); LoBALs tend to have redder optical colours on average and have enhanced radio emission compared to nonBALs (also identified in other studies, e.g., Morabito et al. 2019). However, after removing the BALQSOs from the parent rQSO sample there is still a radio enhancement compared to the typical QSOs, suggesting that the presence of dust is more crucial to the fundamental differences in the radio properties of rQSOs than the presence of a BAL wind (Alexander et al. in prep). This, together with the lack of evidence for stronger winds in rQSOs, could suggest that the radio emission in rQSOs is more closely related to circumnuclear/ISM opacity rather than outflow power. For example, the radio emission could be produced via shocks (causing synchrotron emission) between the wind/outflow and the circumnuclear/ISM gas and dust (e.g., Liu et al. 2014;Zakamska et al. 2016).

CONCLUSIONS
We have studied a sample of rQSOs and cQSOs at 1.45 < < 1.65 using wide-band X-shooter optical-NIR spectra in order to measure the dust-extinction, emission and accretion-disc properties. From our analyses we find that: • The colours of red QSOs are fully consistent with dust extinction (see Fig. 6 and 7): From our extinction-curve fitting analysis, we find that the optical colours in the majority (11/12) of our rQSOs are consistent with dust, with moderate amounts of extinction ( − ) ∼ 0.02-0.23 mags ( ∼ 0.06-0.7 mags). Exploring the nature of the dust in rQSOs, we rule out the flat and MW extinction curves and find similarly good fits for both the steep and PL extinction curves. Analysing the Balmer decrements we find that the rQSOs have larger Balmer decrements, including the upper limits, compared to the cQSOs; this is evidence for more dust along the line-of-sight to the BLR in the rQSOs. See Sections 6.1.1, 6.1.2 and 7.1.
• There are no significant differences in the accretion properties of red and typical (control) QSOs after correcting for dust extinction (see Table 6 and Fig. 12): We find that a simple thin accretion disc can describe the accretion engine of 83 per cent of the rQSOs and 86 per cent of the cQSOs. Refitting the AD models with as a free parameter, we recover consistent accretion parameters and find similar levels of dust extinction to those found in the extinction-curve analysis, for the majority of the QSOs. After correcting for any bias in luminosity between the two samples, we find no significant differences in the accretion properties of rQSOs such as black-hole spin, mass accretion rate or black-hole mass. We tentatively find lower Eddington ratios in the rQSOs compared to the cQSOs, but a larger unbiased sample is required to robustly confirm or refute this result. See Sections 6.2 and 6.3.
• There are no significant differences in the outflow properties of red and typical QSOs: From fitting the emission lines, we find a larger fraction of rQSOs display strong blueshifts in both the C ( > 1000 km s −1 ) and [O ]λ5007 wing ( > 400 km s −1 ) profiles, although these differences are not significant. We also find no significant differences between rQSOs and cQSOs in terms of the FWHM of either the core or wing [O ] Our results suggest that the main difference between red and typical QSOs is the presence of dust rather than significant differences in the accretion or outflow properties. On the basis of these results, a potential self-consistent scenario that links the enhanced radio emission to the dust in red QSOs is a more dust and gas-rich environment, in which the radio emission is due to shocks produced by either an outflow or a jet interacting with a higher opacity ISM/circumnuclear environment in the red QSOs.
This paper demonstrates the capability and quality of X-shooter spectra in determining the nature of red QSOs, and future larger Xshooter samples will robustly tie down and refine many of these results. Our future work will also utilize the Dark Energy Spectroscopic Instrument (DESI; DESI Collaboration et al. 2016) optical spectra (observed wavelength of 360-980 nm) of ∼ 2-3 million quasars, which will push to much fainter and more obscured systems ( up to ∼2-3 mags) than those targeted in the SDSS. This sample, combined with NIR spectroscopic follow-up, can potentially provide us with both with the impressive source statistics and the broader range in required to robustly constrain many of the results presented in this paper. These data will also provide the potential to measure differences in the dust composition with increasing and to, ultimately, conclusively determine the connection between the presence of dust and the enhanced compact radio emission in red QSOs.

APPENDIX A: 2229-0832: A SYNCHROTRON-DOMINATED RED QSO?
The QSO 2229-0832 is the most radio bright QSO in our sample, as well as being the only Fermi-detected source (see Section 2.1). From fitting a dust-reddened cQSO composite to this source, we find that none of the dust-extinction curves provide a good fit; this could suggest that other factors are dominating the red colours, such as synchrotron emission. In a previous study we explored whether synchrotron emission can explain the optical colours in a sample of SDSS DR7 rQSOs (Klindt et al. 2019). To do this, we fitted a synchrotron emission model to the radio-optical data of 3FGL J0045.2-3704 to generate a synchrotron-dominated template (see Appendix A in Klindt et al. 2019, andKlindt et al. 2017 for more details). This template was then scaled from the 1.4 GHz flux of the QSO sample to estimate the expected -band flux from a synchrotron component, comparing to that from the DR7 catalogue. Only ∼ 8 per cent of the radio-detected rQSOs were found to have optical colours that were potentially contaminated by synchrotron emission. Following this approach, we find that the predicted -band magnitude estimated from a synchrotron component is comparable (< 0.36 mag difference) to the SDSS -band magnitude for 2229-0832, suggesting that a synchrotron component may be dominant in this source. For the other sources in our sample, the predicted -band magnitude estimated from a synchrotron component are inconsistent (> 0.5 mags difference) with the SDSS -band magnitude.
One clear signature of synchrotron emission is the characteristic double-bump structure in the SED, usually found in Flat Spectrum Radio Quasars (FSRQs) and BL Lacertae Objects (BL Lacs). The lowest energy bump is due to synchrotron emission in the radio-Xray and the higher energy bump is due to a synchrotron self Compton mechanism in the X-ray-gamma-rays. Therefore finding a gammaray counterpart in addition to a double-bumped SED in a QSO is a good indicator of a dominant synchrotron component. The QSO 2229-0832 displays two clear bumps in the SED which is further evidence that the red colours of this source are likely dominated by synchrotron emission rather than dust.

APPENDIX B: BLACK-HOLE MASS COMPARISON
Fig. B1 displays the extinction-corrected black-hole masses calculated using the broad Mg versus H α from our X-shooter spectra. The squares show the median values for the rQSOs and cQSOs (corrected for dust extinction); the median rQSO black-hole mass estimated from the two emission lines is in good agreement, but the median cQSO H α black-hole mass is slightly lower than that from the Mg emission line. This difference could be due to small number Median rQSO (AV corr) Median cQSO (AV corr) Median rQSO (no AV corr) Median cQSO (no AV corr) Figure B1. The extinction-corrected black-hole masses calculated from the Mg line versus the H α line for our X-shooter spectra. See Fig. 9 for a description of the data, markers and colours plotted.
statistics or the intrinsic scatter between different black-hole mass estimators (Shen & Liu 2012).

APPENDIX C: COMPARISON TO C15
In this section we compare some of the emission-line and black-hole properties for the 15 C15 cQSOs calculated in this work to those from C16 and M16. In Fig. C1 (left) we compare M BH calculated from the H α broad line (using relations from Shen & Liu 2012; see Section 5) to those obtained from AD fitting in C16. Overall, we find good agreement between the two methods. Fig. C1 (middle) shows a comparison of the mass accretion rate obtained from our accretion-disc fitting (see Section 5), compared to that in C16. Fig. C1 (right) shows a comparison of the total Hα emission line luminosities (not corrected for dust extinction) based on our emission line fitting (see Section 4), compared to that from M16.

APPENDIX D: CONTROL COMPOSITE COMPARISON
Our measured continuum extinction values are calculated with respect to our cQSO composite, under the assumption that this is a good representation of the typical QSO population. To test this assumption, we compared our cQSO composite to existing QSO templates: Vanden Berk et al. (2001) (hereafter V01) and Selsing et al. (2016) (hereafter S16).
V01 is the most widely used QSO template in the literature. It is composed of over 2200 SDSS spectra that span a redshift range of 0.044 < < 4.789. Due to the large range in redshifts of the QSOs that comprise this composite and the narrow wavelength range covered by SDSS, a variety of different types of systems will be contributing to the composite at different wavelengths. For example, the lowest redshift QSOs, which will be more host-galaxy dominated, will be predominantly contributing towards the reddest end of the composite, whereas the highest redshift QSOs, and therefore the more luminous by selection biases, will be contributing predominantly towards the bluest end. S16, on the other hand, use only seven QSOs in a redshift range of 1 < < 2, observed by X-shooter, to construct their composite. This has the advantage of the broader wavelength range of X-shooter, combined with a narrower range in redshift, which will result in similar types of systems contributing at all wavelengths. However, the seven QSOs are chosen to be very bright ( < 17 mags), and therefore may not be representative of the overall QSO population. Our X-shooter sample consists of 28 and 12 luminous (16 < < 18 mags) cQSOs and rQSOs. Fig. D1 displays the comparison between our cQSO composite, and those from V01 and S16, all normalized to 1450 Å. Below ∼ 4000 Å, the three composites are fairly similar. Above 4000 Å, the V01 template flattens out compared to both S16 and our composites; this is expected due to the strong host-galaxy contamination from the low redshift QSOs included in the V01 template. The S16 composite appears to be redder relative to our cQSO composite. Comparing the * − * colours of the seven X-shooter QSOs included in the S16 composite to our red and control QSO colour selections, we find that six have colours consistent with our cQSOs and one is slightly redder (although not red enough to be selected as an rQSO); this could be driving the differences in the composite. Equally, the low source statistics in both of the composites, the difference in luminosities between the samples, or the level of host-galaxy contamination could also be affecting the shape.

APPENDIX E: SPECTRA, [O ] AND C PROFILES
Figs. E1, E2, and E3 display the spectra for the rQSOs, cQSOs and three QSOs not included in our final sample (see Section 2.1), respectively. Figs. E4 and E5 display the C and H β-[O ] thumbnails for the cQSOs, respectively, and Figs. E6 and E7 display the C and H β-[O ] thumbnails for the rQSOs, respectively. The red and green lines indicate the narrow and broad Gaussian fits, respectively, and the blue line indicates the overall fitted profile. There are a number of rQSOs and cQSOs that display both broad [O ] and extreme C blueshifts; however, there is are no significant differences between the rQSOs and cQSOs. This suggests that, although both rQSOs and cQSOs can host powerful outflows, it is not clear that they are more prevalent in the rQSO population. This paper has been typeset from a T E X/L A T E X file prepared by the author.  cQSO composite: this work Selsing+2015 Vanden Berk+2001 Figure D1. Comparison of our cQSO composite (cyan; see Section 3) to that in V01 (black) and S16 (blue). All three composite are similar until ∼ 4000 Å; the V01 composite flattens out due to host-galaxy contamination by low redshift QSOs. The S16 composite is slightly redder than our composite, which could be due to a difference in colour in the underlying sample, a difference in luminosities, greater host-galaxy contribution or just small number variations.
2 × 10 3 3 × 10 3 4 × 10 3 6 × 10 3 Rest Wavelength (Å)  Figure D2. Comparison of how the five different extinction laws used in this study affect the continuum shape of the cQSO composite (solid cyan curve) with ( − ) = 0.11 mags. The rQSO composite is also displayed (solid red curve); reddening the cQSO composite with a simple PL gives the best fit to the rQSO composite. The dashed grey lines correspond to the major emission lines. cQSOs Flux (scaled) Figure E4. C profiles for the cQSOs. The y-axis is scaled to 1.4 times the maximum flux. The red and green lines indicate the narrow and broad Gaussian fits, respectively, and the blue line shows the overall fit profile. The source name and rest wavelength of C is also displayed. cQSOs Flux (scaled) Figure E5. H β- [O ] profiles for the cQSOs. The y-axis is scaled to 1.3 times the maximum flux. The red and green lines indicate the narrow core and broad wing fits, respectively, and the blue line shows the overall fit profile. The source name and whether it is radio-quiet (RQ) or radio-loud (RL) is also displayed. Significant [O ] detections are indicated by a star next to the source name.