Quasars were first discovered in the radio band (Schmidt 1963) and soon recognized as luminous extragalactic sources in multiple bands from radio to X-ray. The tremendous energy of quasars originates from the accretion of their central supermassive black holes (SMBHs). Due to their high luminosities, quasars are powerful tools to probe the distant universe. They are often used to study SMBHs, their host galaxies, the intergalactic medium, etc. Surveys over the past 20 yr, e.g., the 2dF QSO Redshift Survey (Boyle et al. 2000), the 6dF QSO Redshift Survey (Croom et al. 2004), the Sloan Digital Sky Survey (SDSS) Quasar Survey (Richards et al. 2006), and the SkyMapper Southern Survey (Wolf et al. 2020;Onken et al. 2022), have searched roughly 1 million quasars (Flesch 2021). Among them, SDSS contributed the majority of the known quasars (Lyke et al. 2020). However, among such a large sample, most of the quasars locate at low redshift. Then the spatial density of high-redshift quasars becomes a critical question.
The quasar luminosity function (QLF) has been widely used to measure how the spatial density of quasars evolves with luminosity and redshift. It has also been used to constrain the quasar contribution to the cosmic X-ray and infrared background (e.g., Hauser & Dwek 2001;Hopkins et al. 2007;Shen et al. 2020) and the contribution of quasar UV photons to the cosmic H and He reionization (e.g., Worseck et al. 2011;Jiang et al. 2016). The QLF at z < 3.5 in the UV/optical band has been well studied. A pure luminosity evolution (PLE) model with a double power-law shape can efficiently describe the QLF at z = 0-2 (e.g., Boyle et al. 1988Boyle et al. , 2000;;Croom et al. 2004), suggesting that the characteristic luminosity in the QLF evolves with redshift while the faint-and bright-end slopes remain unchanged in this redshift range. The PLE scenario is not enough to describe the QLF at higher redshifts. Therefore, a luminosity evolution and density evolution (LEDE) model was proposed to fit the QLF at z ∼ 2-3.5 (e.g., Ross et al. 2013).
At z > 3.5, the bright-end QLF has been measured reasonably well, but the faint end has not been well determined. A full QLF fit usually relies on the combination of large-scale surveys (e.g., SDSS) and small pencil-beam surveys (e.g., Glikman et al. 2010Glikman et al. , 2011)). Using early SDSS data, Fan et al. (2001) studied 39 luminous quasars and suggested that the bright-end shape of the QLF evolves with redshift at z > 3. Glikman et al. (2010Glikman et al. ( , 2011) ) studied the faint end of the QLF at z ∼ 4 using 24 quasars and found a shallow slope α = -1.6 that is consistent with previous studies. For QLFs at z ∼ 5, McGreer et al. (2013) constructed a well-defined sample of 52 quasars from SDSS and measured the QLF at 4.7 < z < 5.1. Yang et al. (2016) extended the bright end of the QLF at z ∼ 5. Table 1 lists some recent studies of UV/optical QLFs that cover a redshift range of z ∼ 3-5. Recently, QLFs at z ∼ 6-7 have also been established (e.g., Jiang et al. 2016;Matsuoka et al. 2018;Wang et al. 2019).
As seen above, significant progress has been made in determining QLFs in different redshift and luminosity ranges. However, the evolution of the quasar population in a wide redshift and luminosity range has not been well characterized. Some studies have tried to analyze such an evolution based on the combination of different quasar samples from the literature (Manti et al. 2017;Kulkarni et al. 2019;Shen et al. 2020;Kim & Im 2021). For example, Kim & Im (2021) selected and combined some binned QLFs at z ∼ 2.4, 3.9, 5.0, and 6.1 from the literature. They found that a pure density evolution (PDE) model is enough to describe the QLFs at 2 < z < 6. In such studies, one has to assume that there are no systematic effects among the individual measurements of QLFs, which is not always the case.
In this paper, we use SDSS multiepoch imaging and followup spectroscopy to construct a well-designed sample of more than 1000 quasars at 3.5 < z < 5, which is roughly 1 mag fainter than the SDSS main quasar sample (Richards et al. 2006). This sample allows us to derive reliable QLFs within a wide luminosity range at high redshift. The layout of the paper is as follows. In Section 2, we introduce the target selection and spectroscopic observations of our quasar candidates. In Section 3, we present our quasar sample, calculate its area coverage, estimate sample incompleteness, and derive QLFs. In Section 4, we compare our result with previous studies and discuss the evolution of the QLF at high redshift. We summarize the paper in Section 5. Throughout this paper, we use point-spread function (PSF) magnitudes, and magnitudes are expressed in the AB system (i.e., SDSS magnitudes are converted to AB magnitudes). We adopt a Λ-dominated flat cosmology with H 0 = 70 km s -1 Mpc -1 , Ω m = 0.3, and Ω Λ = 0.7.
In this section, we will briefly introduce the SDSS imaging survey and then present the details of the quasar candidate selection from the SDSS imaging data. Our program was one of the SDSS-III Baryon Oscillation Spectroscopic Survey (BOSS) ancillary programs, so at the end of the section, we will provide a summary of the BOSS spectroscopic observations of the quasar candidates.
The SDSS is an imaging and spectroscopic survey using a dedicated wide-field 2.5 m telescope (Gunn et al. 2006) with five broad bands, ugriz, at Apache Point Observatory. An SDSS imaging run consists of six parallel scan lines, and two interleaving runs slightly overlap, leading to duplicate observations in a small area. The imaging survey was along great circles and had two common poles, so regions near the survey poles overlap substantially. In addition, if a run (or part of a run) did not satisfy the SDSS quality criteria, the relevant region was reobserved, yielding duplicate observations in this region. Due to the above survey strategy and geometry, SDSS has a large number of duplicate observations (referred to as overlap regions in this paper). The total area of the overlap regions is more than one-third of the SDSS footprint. Detailed information about these overlap regions can be found in Jiang et al. (2015Jiang et al. ( , 2016)). We selected quasars in part of the overlap regions in this paper. These overlap regions provide a unique data set that allows us to select quasars fainter than those found from the SDSS single-epoch data.
In addition to the single-epoch main survey, SDSS conducted a deep imaging survey of ∼300 deg 2 (Stripe 82) on the celestial equator in the South Galactic Cap (Annis et al. 2014;Jiang et al. 2014). Stripe 82 roughly spans 20 h < R. A. < 4 h and -1°.26 < decl. < 1°.26 and was scanned around 70-90 times. The combined data are 1.5-2 mag deeper than the singleepoch data.
There are a variety of quasar selection methods using optical imaging data. Early searches of type 1 quasars largely rely on pointlike morphology and blue UV continuum colors. This color selection is efficient for quasars at relatively low redshift, as quasars and stars have different loci in color-color diagrams (Fan 1999). Later, more methods and more sophisticated techniques were developed. For example, transfer learning has been a useful tool (Fu et al. 2021) for sky regions with large dust extinction and large contamination like the Galactic plane. Other methods, including the likelihood approach (Kirkpatrick et al. 2011), the neutral network approach (Yèche et al. 2010), and extreme deconvolution (Bovy et al. 2011), have been applied to recent quasar surveys such as the SDSS-III BOSS (Ross et al. 2012).
Our goal was to select quasars at 3.6 < z < 5.5. Searches of higher-redshift quasars in SDSS have been carried out by other programs (e.g., Fan et al. 2006;Jiang et al. 2016). We chose to use the traditional color-color diagrams to select our targets. Specifically, we selected quasars at 3.6 < z < 4.5 and 4.5 < z < 5.5 using the gri and riz colors, respectively (see Figure 1 and Table 2). At z > 3.6 (z > 4.5), the Lyα emission line enters the r (i) band, and the Lyman forest absorption makes the quasars much fainter in bluer bands. Therefore, the color-color diagrams are efficient for the selection of quasars in these two redshift ranges (Fan et al. 1999;Richards et al. 2002).
Our targets (like other targets for the SDSS BOSS) were selected from the SDSS Data Release 7 (DR7) imaging data. We first present our quasar selection in the overlap regions. We use 1 to denote the SDSS primary detection and 2 to denote the SDSS secondary detection. For example, i 1 (i 2 ) is the i-band magnitude for the primary (secondary) detection. The SDSS primary detections generally have slightly higher signal-to-noise ratios than the secondary detections. The selection procedure consists of two major steps. In the first major step, we retrieved a preliminary candidate list from the SDSS Query CasJobs online server. We searched the following area at high Galactic latitude: 100°< R.A. < 300°, decl. > -5°, and Galactic latitude b > 40°. The object type of both the primary and secondary detections is "star"; i.e., they were classified as point sources. Although distant quasars are pointlike objects in ground-based images, faint point sources can be misclassified as extended objects. We will correct this effect in Section 3. The positional separation between a primary detection and its secondary detection was required to be smaller than 0 5, which ensures that they are the same object. We excluded objects with the SDSS processing flags "BRIGHT," "EDGE," "SATUR," and "BLENDED." We then imposed initial color cuts to reduce the number of objects in the preliminary target list. The initial cuts are similar to (but much looser than) the final color cuts addressed in Appendix A. We do not expect to lose real quasars in this Figure 1. The ri vs. gr (left) and iz vs. ri (right) color-color diagrams for the illustration of our quasar candidate selection. The black dots are randomly selected point sources that define stellar loci. Each panel shows about 10,000 point sources. In the top panels, the blue and red dots represent randomly selected quasars with i < 20.2 mag from the SDSS DR5 quasar catalog (Schneider et al. 2007). The top left panel includes about 1000 and 500 quasars in the two redshift ranges, respectively. The top right panel includes about 300 and 50 quasars in the two redshift ranges, respectively. The point sources and quasars are not from the same area in the sky. The black lines indicate our quasar selection criteria. These criteria are slightly different in different magnitude ranges. step. In addition, objects with previous spectroscopic observations were excluded using specObjID = 0 (meaning no spectroscopic observations). This was to reduce the number of targets for follow-up spectroscopy.
After we obtained the preliminary list of targets, we combined the primary and secondary detections. For each object, we first converted its primary and secondary magnitudes to flux and calculated the weighted mean flux. Errors were added in quadrature. We then converted the combined flux and errors to AB magnitudes and errors. For example, i (i err ) denotes the combined i-band magnitude (error). The combined magnitudes and errors will be used in the following color selection.
Our color selection criteria were primarily based on the criteria for SDSS I and II from Richards et al. (2002). The selection of candidates at 3.6 < z < 4.5 (gri candidates) was based on the ri versus gr diagram (left column of Figure 1). In Figure 1, the black dots represent randomly selected point sources from the SDSS Query CasJobs online server, the blue and red dots represent a sample of randomly selected quasars from SDSS DR5 (Schneider et al. 2007), and the solid lines indicate our selection criteria. Note that there were no known quasars fainter than i = 20.2 mag here. In order to reduce the number of contaminants, we used slightly different criteria in three different magnitude ranges, 19.0 < i < 20.2, 20.2 < i < 20.8, and 20.8 < i < 21.3 mag. All selection criteria are provided in Appendix A.
The selection of quasar candidates at 4.5 < z < 5.5 (riz candidates) was based on the iz versus ri diagram (right column of Figure 1). We also used slightly different criteria in the three magnitude ranges, and the criteria are shown in Appendix A. These criteria are very similar to those used in the literature (e.g., Richards et al. 2002;McGreer et al. 2013;Wang et al. 2016;Yang et al. 2016).
The target selection in Stripe 82 is straightforward. The combined images and photometric catalogs were available in the archive, and the data were much deeper (Annis et al. 2014). We selected candidates down to i = 21.5 mag and included objects brighter than i = 19.0 mag from the SDSS Query CasJobs online server (using Run = 106 or 206 and specObjID = 0). The overall selection criteria are very similar to Equations (A1) and (A4). They are shown in Appendix A. The search area is 22 h < R.A. < 4 h . We did not use the region of R.A. < 22 h , as the Galactic latitude becomes lower.
Our targets were observed by the BOSS spectrograph in SDSS-III (Eisenstein et al. 2011;Dawson et al. 2013). The BOSS main survey covered ∼10,000 deg 2 in the north and south galactic caps and was completed in 2014 (Alam et al. 2015). The BOSS quasar survey mainly focused on quasars at 2.2 < z < 3.5. Our program was selected as one of the ancillary programs to fill spare fibers. The SDSS bitmasks used in SDSS targeting can be found on the websitefoot_1 of our program. From the selection procedure above, we obtained 4374 quasar candidates, including 3454 candidates from the overlap regions and 920 candidates from Stripe 82. A total of 3406 candidates were spectroscopically observed. The mean fraction of targets with spectroscopic observations reaches about 78%, and we will correct this incompleteness in Section 3.3.
From the spectroscopic observations, we obtained 887 quasars. They have been included in the SDSS DR16 quasar catalog (Lyke et al. 2020). Their redshift distribution is shown in Figure 2. This sample consists of 35 quasars at z < 3 and 852 quasars at z > 3. The z > 3 sample (hereafter the new sample) includes 652 quasars in the overlap regions and 200 quasars in Stripe 82.
As we mentioned earlier, we did not observe the objects that had already been spectroscopically observed in SDSS I and II. Some of them also satisfy our target selection criteria. We recovered this quasar sample (hereafter the archival sample) as follows. For quasars in the overlap regions, we changed one criterion (using specObjID! = 0) and repeated the selection procedure. For quasars in Stripe 82, we directly used the criteria to match quasars in the DR7 quasar catalog (hereafter DR7Q; Schneider et al. 2010). We recovered a total of 346 quasars. Our final high-redshift sample is the combination of the archival and new samples and consists of 1198 quasars at z > 3 (Table 2).
We measure the continuum properties of the quasars assuming a power-law shape f l µ l a l . We fit this power law to the spectral regions with little line emission. The resultant slope α λ distribution is shown in Figure 3. The mean α λ value is about -1.1, similar to previous measurements of highredshift quasars (e.g., Fan et al. 2001;Schneider et al. 2001) but much softer than the results from low-redshift works due to the short wavelength coverage for high-redshift quasars (Vanden Berk et al. 2001). Continuum luminosity/magnitude M 1450 is also calculated in this step. Figure 4 shows the redshift and M 1450 distributions of the archival and new samples. The median value of M 1450 is around -25.5 mag. The new sample is about 1 mag deeper, on average, so the QLF calculated in this paper will reach a lower luminosity compared to that from the SDSS single-epoch data.
The area coverage of the SDSS overlap regions is complex. We use Hierarchical Equal Area isoLatitude Pixelization (HEALPix; Górski et al. 2005) to estimate the effective area of the overlap regions. The basic idea is to pixelize the sky sphere into a mesh of quadrilateral pixels, and the effective area is calculated by adding up all pixels that cover our data points.
For a given data set, the starting resolution level of HEALPix is important for the area calculation. We follow Jiang et al. (2016) and adopt HEALPix level 10 (i.e., 11.8 arcmin 2 pixel -1 ) as the best starting level for the overlap regions (see Figure 5 in Jiang et al. 2016). We classify all pixels into three categories. The first category consists of empty pixels that do not cover any objects, so they do not contribute to the effective coverage. The close neighbors to empty pixels are boundary pixels that will result in the uncertainty of the area calculation. The remaining pixels are all in the third category (hereafter nonboundary pixels). All nonboundary pixels at level 10 contribute to the total effective area. We then gradually increase the resolution level for the boundary pixels. The resultant new pixels are once again classified into two categories, boundary and nonboundary pixels. The new nonboundary pixels are added to the total area, and the new boundary pixels are refined again by increasing the resolution level. This procedure stops when the resolution roughly matches the average surface density of the data points. In this paper, we reach the best resolution at level 12.
From the above calculation, the total area of the overlap regions is 1292 ± 266 deg 2 . There are two main types of overlap regions. One is that a large field of the sky or a whole SDSS run was observed twice or more. In this case, the relative area uncertainty is very small. The extreme case is Stripe 82. The other type is that only (part of) the narrow scan lines in an SDSS field or run are overlap regions. In this case, there is a significant fraction of boundary pixels in the area calculation that produce a large uncertainty. Many overlap regions in this work belong to the second type, so the uncertainty of the calculated total area is nonnegligible. This uncertainty will be included in the measurement of our QLF later.
The calculation of the Stripe 82 area is very straightforward, since it is one rectangular piece of the sky. Its coverage is 225 deg 2 with a negligible uncertainty.
In this subsection, we estimate our sample incompleteness, which is critical to derive QLFs. The first incompleteness is from the fact that BOSS did not observe all of our targets. For example, about 82% (71%) of the targets in the gri (riz) sample at i < 20.2 mag were observed in the overlap regions. When we correct this incompleteness, we assume that the quasar fraction in the unobserved candidates is the same as that in the observed sample. The incompleteness slightly varies with the i-band magnitude. This variation is considered as the uncertainty of this incompleteness, and the results (1%-3%) are negligible.
The second incompleteness arises from the morphological bias. The quasar candidates that we observed are point sources, but faint point sources with low signal-to-noise ratios can be misclassified as extended sources by the SDSS photometric pipeline. We correct this bias for the targets in the overlap regions (i.e., the single-epoch data). The Stripe 82 imaging data are much deeper, and we assume that the targets in this region do not suffer from the morphological bias. We will see below that this assumption is reasonable. To estimate the bias for the single-epoch data, we use 27,593 point sources classified in Stripe 82. We divide the data into narrow magnitude bins. For each magnitude bin, we calculate the fraction of the objects that are misclassified as extended sources in the single-epoch data.
The resultant fractions are from 0.04 at 19.0 mag < i < 19.6 mag to 0.55 at 21.2 mag < i < 21.3 mag. There is a clear relation between the fraction and brightness. The fraction in the brightest range is nearly zero, suggesting little bias for the targets in Stripe 82. These fractions are considered as sample incompletenesses in individual bins and will be included when we calculate QLFs. In order to estimate the uncertainty of the incompleteness, we resample the data 1000 times. For each time, we randomly select 500 point sources per bin to estimate the incompleteness and finally regard the standard deviation as uncertainty. The resultant uncertainties (1%-2%) are negligible.
The next incompleteness comes from our color selection criteria, i.e., the color cuts introduced in Section 2. This incompleteness is described by a selection function, the probability that a quasar with a given magnitude (M 1450 ), redshift (z), and intrinsic spectral energy distribution (SED) meets the color selection criteria. We calculate the average selection probability p s (M 1450 , z) by assuming that intrinsic SEDs have certain distributions. Following the procedure in Fan (1999) We follow the procedure in Jiang et al. (2016) to add photometric errors. We withdraw a large representative sample of point sources in the overlap regions and Stripe 82 and derive error distributions as a function of magnitude in the u, g, r, i, and z bands. Finally, we construct a grid of 2 million mock quasars in a redshift range of 3.5 < z < 5.5 and a luminosity range of -27.5 < M 1450 < -23.5, with step sizes of ΔM = 0.02 and Δz = 0.02. Then we calculate the selection function p s (M 1450 , z), the fraction of simulated quasars that meet our selection criteria. Figure 5 shows the selection functions for the overlap regions and Stripe 82. We estimate the uncertainty of the selection function using the same method as we did for the morphological incompleteness, and the result is around 1%. As we will see, the uncertainties of the incompleteness corrections are negligibly small compared to other uncertainties, so they are not included in the following QLF calculations.
We use a traditional 1/V a method (Avni & Bahcall 1980) to derive the binned differential QLFs. The available volume for a quasar with absolute magnitude M and redshift z in a magnitude bin ΔM and redshift bin Δz is
where p(M, z) is the final selection function that includes all incompleteness corrections discussed above.
In general, the binned QLF and its statistical uncertainty can be expressed as
where the sum is over all quasars in each bin. When a density approaches the Poisson limit, its uncertainty is corrected using Equation (7) in Gehrels (1986).
We divide our sample into several luminosity and redshift bins. We focus on three redshift ranges (3.6 < z < 4.0, 4.0 < z < 4.5, and 4.5 < z < 4.9) that include a subsample of 1106 quasars. This subsample is used for our QLF measurement. The magnitude limits in each redshift range are determined by the faintest and/or brightest quasars and the selection functions. The binned QLF results are listed in Table 4 and displayed in Figure 6 as the blue (overlap regions) and red (Stripe 82) circles. The horizontal locations of the symbols are at the centers of each magnitude bin, and the horizontal bars indicate the magnitude coverage ranges. The binned QLFs calculated for the two data sets are consistent within 1σ.
We combine the two data sets from the overlap regions and Stripe 82 and derive a parametric QLF using the maximumlikelihood method (Marshall et al. 1983). This method aims to minimize the function S, which is equal to -L 2 ln , where L is the likelihood function,
where p(M, z) includes all of the incompleteness corrections discussed above. The first term is the sum over all observed quasars in the sample. The second term is integrated over the whole magnitude and redshift range of the sample. It represents the total number of expected quasars for a given luminosity function. The confidence intervals are determined from the logarithmic-likelihood function using a χ 2 distribution of ΔS ( S S min = -) (Lampton et al. 1976). We choose a double power-law form (Boyle et al. 2000) as the parametric QLF model,
where α and β are the faint-and bright-end slopes, M * is the characteristic magnitude (or break magnitude), and Φ * is the density normalization. We assume that these parameters do not change in small redshift ranges, such as the ranges considered here. We will discuss the QLF evolution in Section 4.2. We perform a grid search to determine the best-fit results and the confidence intervals. The grid resolutions of log Φ * , M * , and α are 0.05, 0.05, and 0.1, respectively. There is a strong degeneracy between M * and α, so we set a bright limit of -28.0 mag for M * . The best-fit results are listed in Table 5.
Figure 6 shows the results in three redshift ranges. The open circles denote the data points (some of the faintest bins) that have very low completeness and deviate significantly from the general trend. It is unclear what causes this deviation. This has frequently been seen in previous studies and is likely due to some unknown selection effects. We did not use these data points in the above calculation. Our sample covers a limited range of luminosity, so it is not able to constrain both slopes α and β. Therefore, we use three fixed values for β in each redshift range (see Figure 6) and derive the other three parameters. The best-fit α values are about -1.8 at 3.6 < z < 4.9, indicating that the results for different β values and redshift ranges are not significantly different. In addition, most of the best-fit M * values for three redshift ranges are lower than -27 (see Table 5), making the double power-law model degenerate into a single power-law model. These results suggest that our sample alone is not enough to constrain all parameters in the above QLF model. In the next section, we will combine our binned QLFs with some results in the literature.
In Figure 7, we show a collection of previous QLF measurements at 3.6 < z < 4.9 (Richards et al. 2006;McGreer et al. 2013McGreer et al. , 2018;;Yang et al. 2016;Boutsia et al. 2018Boutsia et al. , 2021;;Kulkarni et al. 2019;Schindler et al. 2019;Kim et al. 2020) Notes. The magnitudes for the overlap regions are combined magnitudes, as introduced before. The magnitudes for Stripe 82 are from the SDSS Query CasJobs online server or DR7Q, depending on whether they have previous spectroscopic observations. All magnitudes are expressed in the AB system and have been corrected for the extinctions. a Overlap: the overlap regions; S82: Stripe 82.
(This table is available in its entirety in machine-readable form.) QLFs and study quasar evolution from z = 7.5 to 0. In the bottom panel of Figure 7, we particularly compared our results with previous QLF measurements at z ∼ 5 (e.g., McGreer et al. 2013McGreer et al. , 2018;;Yang et al. 2016;Kim et al. 2020). We did not include samples with no or very few spectroscopic observations (e.g., Akiyama et al. 2018;Niida et al. 2020).
Figure 7 shows that our results are generally consistent with previous measurements. In the top panel, our luminosity coverage for z = 3.8 partly overlaps with the luminosities covered by Richards et al. (2006), Boutsia et al. (2018), and Kulkarni et al. (2019). In this overlap range, the binned QLFs from different studies roughly agree with each other. Our binned QLFs at -26.5 < M 1450 < -25.5 are about 1.5 times the results of Richards et al. (2006). It is unclear whether their results were underestimated or our results were overestimated. In the middle panel, for z = 4.25, our result is well consistent with the previous results from Richards et al. (2006) and Kulkarni et al. (2019) except for two high data points from Kulkarni et al. (2019). The bottom panel shows several studies of the QLF at z = 4.7, and most of these results are consistent with ours within a 1σ level. It is worth noting that the discrepancies are relatively larger at the faint and bright ends, where the uncertainties are also significantly large.
In order to better constrain the shape of QLFs, we combine our QLF measurements with the results from some previous Notes. a Φ is in units of Mpc -3 mag -1 . b ΔΦ is in units of 10 -9 Mpc -3 mag -1 . studies (Richards et al. 2006;McGreer et al. 2013McGreer et al. , 2018;;Yang et al. 2016;Boutsia et al. 2018Boutsia et al. , 2021;;Kulkarni et al. 2019;Schindler et al. 2019;Kim et al. 2020). Detailed quasar selection functions in these studies are not all publicly available, so we choose to use their binned QLFs. We note that results from binned data may be subject to biases (e.g., La Franca & Cristiani 1997;Page & Carrera 2000). In addition, previous samples are not fully independent, as they included some common samples of quasars. To reduce the impact from potential biases, we only use our results in the luminosity ranges that our sample covers. In addition, we exclude those points with very low completeness. The points that are not used in the fitting are shown as the open symbols in Figure 7, and the purple line in each panel denotes the best-fit QLF. We fit the observed QLF data points (Φ obs ) using the maximumlikelihood estimation. We use the logarithmic-likelihood function ,
where σ obs is the 1σ uncertainty of Φ obs from the literature. We use the emcee Python packagefoot_2 (Foreman-Mackey et al. 2013) for the Markov Chain Monte Carlo sampling of the QLF parameters. The chosen priors on the parameters are shown in Appendix B. The best-fit results and uncertainties are estimated based on the 16th, 50th, and 84th percentiles of the samples in the marginalized distributions. Figure 7 shows the best-fit QLFs, and the parameters are listed in Table 5. ) with reduced errors. The binned QLF data points used in Boutsia et al. (2021) are higher than our measurements, resulting in a little steeper faint-end slope α and a higher-density Φ * than our results. The QLF at z = 4.25 is also well constrained by our results. As mentioned earlier, the faintest binned QLF calculated by Kulkarni et al. (2019) is very high, so they obtained a steep faint-end slope 2.20 0.14 0.16 a = --+
. Our slope α is slightly flatter.
For the QLF at z = 4.7, we adopt the results with 1.8 0. for the bright end. In this high-redshift range, the current sample size is still small; thus, the above measurements are associated with relatively large uncertainties.
Based on the measurements from the above subsection, we explore quasar evolution at 3.5 < z < 5. We use a linear function to describe the evolution of the four QLF parameters in this small redshift range,
Here we consider four cases. 1.8 0.2 0.5 We fit these four models to the observed QLF data points using the maximum-likelihood estimation introduced earlier. We calculate the reduced χ 2 as
where n is the number of data points, and ν is the number of free parameters. The results are listed in Table 5. These four models have a similar fitting performance in terms of 2 c n . From cases 1-4, the performance is only mildly improved. Besides, the PDE model in case 1 has the smallest number of free parameters and is enough to describe the evolution of the QLF at z ∼ 3.5-5. This is consistent with the conclusion of Kim & Im (2021).
We compare the four models in Figure 8. The case 1 (purple solid line) and case 4 (purple dashed line) results have a similar fitting performance, so it is difficult to find the best-fit model based on current data. These two results are both higher than the result of Kulkarni et al. (2019) in the bright end. This discrepancy can be attributed to the different data that the two studies used. For example, Kulkarni et al. (2019) did not use the data from Boutsia et al. (2018) and Schindler et al. (2019). These data may increase the measurement of the bright-end QLF. In addition, the redshift coverage of Kulkarni et al. (2019) is much larger than that of this work. The data in the range of z < 3.5 and z > 5 will affect the result in 3.5 < z < 5 when computing the QLF evolution.
Compared with Kim & Im (2021), our QLF slopes are steeper. Such steep slopes in the faint and bright ends have been found in recent measurements (e.g., Giallongo et al. 2019;Grazian et al. 2022). The measurements of the slopes sensitively depend on the data points at the brightest and faintest ends that usually suffer from large incompleteness and uncertainties. In addition, the combination of different samples introduces extra uncertainties that are often difficult to characterize. A large sample with full coverage of both ends is needed to improve the measurement of the QLF.
Finally, we explore the cumulative spatial density evolution of quasars at high redshift. The spatial density of quasars brighter than a given magnitude M is calculated by integrating the QLF,
where we use the PDE model (case 1) of this work as Φ(M, z).
Figure 9 shows the cumulative density as a function of redshift for different magnitude ranges using our model and previous results (Richards et al. 2006;McGreer et al. 2013;Ross et al. 2013;Jiang et al. 2016;Yang et al. 2016;Wang et al. 2019;Kim & Im 2021). It indicates a rapid density decline from z ∼ 3.5 to 7, consistent with previous studies. For example, Fan et al. (2001) fit an exponential decline, ρ(<M, z) ∝ 10 kz , and found k = -0.47 at high redshift. McGreer et al. (2013) found that the slope at 4 < z < 5 was k = -0.38 ± 0.07. Previous results also suggested that the spatial density of quasars drops faster with increasing redshift at z > 3.5. As shown in Figure 9, the slopes at 5 < z < 6 and 6 < z < 7 are k = -0.72 ± 0.11 and -0.78 ± 0.18, respectively (Jiang et al. 2016;Wang et al. 2019). The PDE model from Kim & Im (2021; dashed lines) also supports this scenario, as shown in Figure 9.
From our sample, we measured k = -0.7 ± 0.1 at 3.5 < z < 5, which is slightly steeper than previous measurements for the same redshift range but similar to those measured at z = 5-7. The main reason for the steeper slope is that our binned QLF within -26.5 < M 1450 < -25.5 at z = 3.8 is about 1.5 times the previous results (see gray squares in Figure 7 and blue crosses in Figure 9). However, the cumulative densities between our model and the observed results are consistent within 1σ. To confirm our results, we need a larger and more complete sample, such as a quasar sample from the Chinese Space Station Telescope widearea slitless spectroscopic survey (Zhan 2021). If confirmed, the quasar density at 3.5 < z < 5 declines faster than previous measurements and as fast as the density evolution at z > 5.
In this paper, we have built a sample of more than 1000 quasars at z > 3, including 974 quasars in 1292 deg 2 of the SDSS overlap regions and 224 quasars in 225 deg 2 of Stripe 82. The spectroscopic observations were conducted by the SDSS-III BOSS. The sample spans an absolute magnitude range of -27.5 mag < M 1450 < -24.0 mag. This is roughly 1.5 mag fainter than the SDSS main quasar sample selected from the single-epoch data.
We have constructed QLFs at 3.5 < z < 5 based on this sample and studied quasar evolution from z = 5 to 3.5. We first corrected sample incompleteness caused by the misclassification of the object morphology, the color selection of the candidates, and the incomplete spectroscopy of the candidates. We then derived the binned QLFs at 3.6 < z < 4.0, 4.0 < z < 4.5, and 4.5 < z < 4.9 and modeled the QLFs using a double power-law form. The luminosity coverage of our sample is not large enough to constrain all parameters in the double power-law model, so we fixed the bright-end slope β. We found that the faint-end slopes for the three redshift ranges are α ∼ -1.8, with moderate to large uncertainties from 0.1 to 0.3 to >0.5. The relatively large uncertainties are mainly due to the relatively small sample size and the fact that our sample does not reach a very low luminosity.
We have made use of some studies from the literature and improved the measurement of the QLFs. We combined their binned QLFs with ours and characterized the QLFs in a larger luminosity range of -29 mag < M 1450 < -23 mag. We found that the faint-end slopes of the QLFs are around -1.7, and the bright-end slopes are from -4.0 to -3.5. Finally, we investigated the evolution of the QLFs from z ∼ 5 to 3.5 and found that a simple PDE model can efficiently describe the QLF evolution in this redshift range. This is consistent with some recent results. We also found that the quasar density at 3.5 < z < 5 declines faster than previously thought, and its evolution parameter k is similar to that at z > 5.
Figure 9. Cumulative density evolution of quasars at 0 < z < 7. The green, cyan, and orange lines represent magnitude ranges M 1450 < -25, -26, and -27 mag, respectively. The bold solid lines are the densities at 3.5 < z < 5 calculated from the PDE model in this paper. The cyan shaded region indicates the 1σ uncertainty of our PDE model. The solid lines at the ranges of 0 < z < 3.5, 5 < z < 6, and 6 < z < 7 are from Ross et al. (2013), Jiang et al. (2016), and Wang et al. (2019), respectively. The dashed lines show the PDE model at 2 < z < 6 from Kim & Im (2021). The blue symbols denote observational cumulative densities at z = 3.75, 4.25, 4.9, 5.05, 6.1, and 6.7 measured by Richards et al. (2006;crosses)
The criteria for the gri candidates in overlap regions at 20.2 mag < i < 20.8 mag are as follows:
The criteria for the gri candidates in overlap regions at 20.8 mag < i < 21.
The criteria for the riz candidates in overlap regions at 20.8 mag < i < 21.3 mag are as follows:
The Astrophysical Journal, 928:172 (14pp), 2022 April 1 Pan et al.
https://www.sdss.org/dr12/algorithms/ancillary/boss/highz/
https://emcee.readthedocs.io/en/stable/