ON THE RADIO AND OPTICAL LUMINOSITY EVOLUTION OF QUASARS

The optical emission of quasars, or active galactic nuclei (AGNs), is dominated by the radiation of the plasma accreting onto supermassive black holes, while the radio emission is dominated by the plasma outflowing from the black hole/accretion disk systems. Hence different but complementary information can be gathered in both photon energy ranges regarding the cosmological evolution of AGNs and its relation to structure formation in the universe. It is therefore important to analyze in detail redshift distributions of quasars in both frequency regimes, investigating carefully any possible differences between these two.

The rapid evolution of active galaxies identified in radio catalogs as "quasi-stellar radio sources" (or QSRs) in the redshift range z ≲ 2 was established soon after their discovery (e.g., Schmidt 1968). Subsequent optical discoveries of similar sources, most of which had no detectable radio emission, led to the emergence of the class of "radio-quiet quasars" (or "quasi-stellar objects," QSOs for short; e.g., Osterbrock 1989). These optically selected sources also showed similar strong evolutionary trends, similar to the radio-selected ones. These evolutions are modeled as density evolution, luminosity evolution, or a combination of the two in numerous works (e.g., Schmidt 1968; Petrosian 1973; Marshall et al. 1983; Dunlop & Peacock 1990; Maloney & Petrosian 1999; Willott et al. 2001) and can be designated as the evolution of the luminosity function (LF, for short).

By now the evolution of the LF has been described not only for optical and radio luminosities but also for X-ray, infrared, and bolometric luminosities (e.g., Ueda et al. 2003; Richards et al. 2006; Matute et al. 2006; Hopkins et al. 2007; Croom et al. 2009). Most of these studies have treated the evolution with a bi-variate function Ψ_i(L_i, z), where L_i is the luminosity (or, in this case, luminosity spectral density) in some photon energy range, e.g., L_i = L_opt or L_rad. The shape of the LF and its evolution are usually obtained from a flux-limited sample f_i > f_{m, i} with L_i = 4 π d²_L(z) (1/K_i(z)) f_i, where d_L is the luminosity distance and K_i(z) stands for the K-correction. For a power-law emission spectrum of index ε_i defined as $f_i \propto \nu ^{-\varepsilon _i}$, one has $K_i(z) = (1+z)^{1-\varepsilon _i}$.

However, because no matter how a quasar is discovered, optical observations are required for determination of redshift, then the flux limit of optical observations (f_{m, opt}) and the optical luminosity enter the picture, so that one now must consider the joint LF and its evolution Ψ(L_opt, L_i, z), a tri-variate function with L_i = L_rad or L_X, for example.

In general, the first step required for investigation of a multi-variate distribution is the determination of whether the variables of the distributions are correlated or are statistically independent. For example, in the case of a single LF, the correlation between L and z is what we call luminosity evolution, and independence of these variables would imply absence of such evolution. Mathematically, independence means that the function is separable, Ψ_i(L_i, z) = ψ_i(L_i) × ρ(z), in which case one is left with the determination of a single variable LF ψ_i(L_i) and the density evolution ρ(z). As shown by Petrosian (1992), the most exact nonparametric method for this task from a flux-limited (or a more generally truncated) sample is the Lynden-Bell (1971) method. However, this simple and elegant method cannot be used for cases when variables are correlated (e.g., when there is luminosity evolution). Efron & Petrosian (1992, 1999) (EP for short) developed new methods for determination of the existence of correlation or independence of the variables from a flux-limited and more generally truncated data set, and prescribed how to remove the correlation by defining new and independent variables (say L'_i ≡ L_i/g_i(z) and z, where the function g_i(z) describes the luminosity evolution) and then how to determine the mono-variate functions ψ_i(L'_i) and ρ(z). Thus, one can determine both the luminosity and density evolutions g_i(z) and ρ(z), as well as the LF at any redshift.⁶

In the case of quasars with the optical and some other band luminosity, we have at least a tri-variate function. In this case, one must determine not only the correlations between the redshift and individual luminosities (i.e., the two luminosity evolutions) but also the possible correlation between the two luminosities, before individual distributions can be determined. Knowledge of these correlations and distributions is essential for not only constraining robustly the cosmological evolution of active galaxies, but also for interpretation of related observations, such as the extragalactic background radiation (e.g., Singal et al. 2010; Hopkins et al. 2010).

Another related aspect of this subject, which has attracted considerable attention over the years, is the distribution of the"radio-loudness parameter" for the quasar population, and the distinction between the so-called radio-loud (RL for short) and radio-quiet (RQ for short) quasars. The question of whether there are two distinct populations was addressed soon after the discovery of quasars using small samples with radio flux limits greater than 1 Jy. Initially it was found that the distribution of the ratio of radio to optical luminosity R ≡ L_rad/L_opt, the so-called radio-loudness parameter, was a fairly broad power law with index β_R ∼ −2.3 in the range $2.8<{\rm \log }R<5.2$ (Schmidt 1972; Petrosian 1973). At that time this ratio was defined for the radio luminosity at ν_rad = 0.5 GHz and optical luminosity at 2500 Å  (or the frequency ν_opt = 1.2 × 10¹⁵ Hz). Nowadays it is defined with radio luminosity at 5 GHz so that for a mean radio spectral index ε_rad ∼ 0.6 the old data would be in the range 2.2–4.6 of the modern definition of log R.⁷ Later, however, the survey limits were extended to much lower fluxes (especially in the radio domain), and this has resulted in a much wider range of the ratio that extends to values well below one, namely $-3< {\rm \log } R<5$. Within this broader range, weak hints of the bi-modality described by Kellerman et al. (1989) suggested that log R = 1 could be chosen as the RL/RQ demarcation value. Using this value for the division between RL and RQ quasars, the differences between the two classes have been investigated, including the possibility of distinct cosmological evolution of the RL and RL populations (e.g., Miller et al. 1990; Goldschmidt et al. 1999; Jiang et al. 2007). Still, the more recent analyses of different samples of objects reported in the literature so far gave rather inconclusive results on whether any bi-modiality in the distribution of the radio-loudness parameter for quasars is inherent in the population (see Ivezić et al. 2002, 2004; Cirasuolo et al. 2003).

There have been many papers dealing with this ratio and RL versus RQ issue, as well as luminosity ratios at other wavelengths, e.g., IR/radio, optical/X-ray, etc. However, none of these works have dealt with the intrinsic distribution (and/or evolution) of the ratio, which is related to the tri-variate LF Ψ(L_opt, L_rad, z) by⁸

$$\begin{eqnarray} G_R(R,z) &=& \int _0^{\infty } { \Psi (L_{{\rm opt}}, R \,\, L_{{\rm opt}}, z) \, L_{{\rm opt}} \, dL_{{\rm opt}} } \nonumber \\ &=& \int _0^{\infty } { \Psi \left({L_{{\rm rad}} \over R}, L_{{\rm rad}}, z \right) \, L_{{\rm rad}} \, {{ dL_{{\rm rad}} } \over {R^2}} }. \end{eqnarray} \tag{ 1 }$$

These works did not take the observational selection effects properly into consideration nor did they address the correlations between the radio and optical luminosities. Neglecting these effects when attempting to determine the distribution of radio loudness is usually given the justification that the ratio is essentially independent of cosmological model and redshift (as long as the K-corrections are the same). The broad distribution of observed ratios obtained in this way (see Figure 14 below) deviates from a simple power law and may even have a hint of bi-modality, seemingly justifying at face value the choice of log R ∼ 1 as the separation point between RL and RQ sources. However, as shown in Appendix A, even in the simplest cases the observed distribution (and its moments) could be very different from the intrinsic ones. Thus, for determination of the true distribution of R the data truncations must be determined and the correlations between all variables must be properly evaluated.

Our aim in this paper is to take all these effects into account in determination of the evolution of optical and radio luminosities and their ratio and to find their distributions, with particular emphasis on the radio loudness question. In Section 2, we describe the data we use. In Section 3, we provide an overview of the procedure used. In Section 4, we present our results on the correlations and evolutions of the LFs. In Section 5, we describe the density evolution and the luminosity density evolutions, while in Section 6 we calculate the LF corrected for luminosity evolution, which we call the "local" LF. Finally, in Section 7 we evaluate the distribution of radio loudness, R. This work assumes the standard ΛCDM cosmology throughout, with H₀ = 71 km s⁻¹ Mpc⁻¹, Ω_Λ = 0.7, and Ω_m = 0.3.

In order to evaluate the luminosity evolution in both radio and optical, and to separate and compare these effects, we require a data set that has both radio and optical fluxes to reasonable limits and across a range of redshifts, that contains a significant number of both RL and RQ objects. The overlap of the FIRST (Faint Images of the Radio Sky at Twenty centimeters) bright quasar radio survey with the Automatic Plate Measuring Facility catalog of the Palomar Observatory Sky Survey (POSS-I), as presented by White et al. (2000), is such a data set. It contains 636 objects with optical R-band optical magnitudes, 1.4 GHz total and peak pixel fluxes, and spectroscopic redshifts. The survey has a limiting R-band magnitude of 17.8 or f_{m, R} = 0.22 mJy, a limiting peak pixel 1.4 GHz flux of 1 mJy, and redshifts that range from 0.02 to 3.425. Figures 1 and 2 show the radio and optical luminosities versus redshifts of the quasars in the survey, assuming the standard K-corrections for power laws with optical and radio indices ε_opt = 0.5 and ε_rad = 0.6, respectively (where $f_i \propto \nu ^{-\varepsilon _i}$). Figure 3 shows the radio-loudness parameter R versus redshift for the data set.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

This sample spans a very wide range of luminosities (5 dex in optical and 7 dex in radio) with a significant number of sources in the range 0.1 < R < 10⁴ (with 336 RL and 300 RQ). Therefore, it is well suited for our analysis here. We have examined some other combined radio and optical survey data sets and found them to be not as well suited for this analysis. For example, the combined FIRST radio survey with the 2dF optical survey as reported in Cirasuolo et al. (2003) features only 12 RQ objects (of 113 total), and the combined FIRST with the Large Bright Quasar Survey as reported in Hewett et al. (2001) has only 77 objects and different optical flux limits for the various fields, making the method employed here cumbersome. Of course a much larger sample could be achieved combining optical data from the Sloan Digital Sky Survey (SDSS) with FIRST radio data, as in Jiang et al. (2007). Forming such a data set, however, would necessitate care in associating optical sources with radio ones, given the very large number of optical sources. Issues such as how large a radius to allow for an association and what to do about multiple matches are important considerations without immediate answers. It is our intention to demonstrate the techniques employed here with a well-established smaller data set before moving on to a more comprehensive one involving SDSS optical data.

The LF gives the number of objects per unit comoving volume V per unit source luminosity, so that the number density is dN/dV = ∫dL_iΨ_i(L_i, z). To examine luminosity evolution, without loss of generality, we can write an LF in some waveband i as

$$\begin{equation} \Psi _{i}(L_i,z) = \rho (z)\,\psi _i\big(L_{i}\big/g_{i}(z), \eta _i^j\big)\big/g_i(z), \end{equation} \tag{ 2 }$$

where g_i(z) and ρ(z) describe the luminosity evolution and comoving density evolution with redshift, respectively, and η^j_i stands for parameters that describe the shape (e.g., power-law indices and break values) of the i-band LF (we use the normalization ∫^∞₀ψ_i(L_i)dL_i = 1).⁹ In what follows we assume a non-evolving shape for the LF (i.e., η^j_i = constant, independent of L and z), which is a good approximation for determining the global evolutions. Once these are determined, this hypothesis can be tested and results amended. However, for more complicated functional forms with variable η^j_i, e.g., for luminosity-dependent density evolution (LDDE), the determination of the variations will require a large sample with significant numbers of objects in reasonably narrow redshift and luminosity bins.

Given this assumption then once the luminosity evolution g_i(z) is calculated, the density evolution ρ(z) and local LF ψ_i(L'_i) ≡ ψ_i(L_i/g_i(z))/g_i(z) can be determined.¹⁰ We consider this form of the LF for luminosities in different bands, allowing for separate (optical and radio) luminosity evolution.

• 1.  
  As is often done, one might naively assume that the joint LF Ψ(L_opt, L_rad, z) is separable into two forms like Equation (2) with a common density evolution. However, as discussed in Section 1, because the optical and radio luminosities of the quasars are, in general, highly correlated, the simultaneous determination of the LFs of both requires care. The first step in this procedure should be to determine the degree and form of the correlation between the optical and radio luminosities. As described below, the EP method allows us to determine whether any pair of variables are independent or correlated. Once it is determined that they are correlated one should seek a coordinate transformation to define a new pair of variables which are independent. This requires a parametric form for the transformation. One can define a new luminosity which is a combination of the two; we can define a "correlation-reduced radio luminosity" L_crr = L_rad/F(L_opt/L_fid), where the function F describes the correlation between L_rad and L_opt and L_fid is a fiducial luminosity taken here to be 10²⁸ erg s⁻¹ Hz⁻¹. This is a convenient choice for L_fid as it is lower than the lowest 2500 Å luminosity considered in our sample, but results do not depend on the particular choice of numerical value. For the correlation function we will assume a simple power law

  $$\begin{equation} L_{{\rm crr}} = {{L_{{\rm rad}}} \over {(L_{{\rm opt}}/L_{{\rm fid}})^{\alpha }}}, \end{equation} \tag{ 3 }$$

  where α is a bulk power-law correlation index to be determined by a fit to the data. This is essentially a coordinate rotation in the log–log luminosity space. As shown in Section 4 below, EP also prescribe a method to determine a best-fit value for the index α which orthogonalizes the new luminosities. Given the correlation function we can then transform the data (and its truncation) into the new independent pair of luminosities (L_opt and L_crr), whose distribution can be represented as

  $$\begin{eqnarray} \Psi (L_{{\rm opt}},L_{{\rm crr}},z) &=& \rho (z) \, \times \, \psi _{{\rm opt}}\big(L_{{\rm opt}}\big/g_{{\rm opt}}, \eta _{{\rm opt}}^j\big)\big/g_{{\rm opt}} \nonumber \\ && \times \psi _{{\rm crr}}\big(L_{{\rm crr}}\big/g_{{\rm crr}}, \eta _{{\rm crr}}^j\big)\big/g_{{\rm crr}}. \end{eqnarray} \tag{ 4 }$$
• 2.  
  The next step is determination of the two independent luminosity–redshift correlation functions g_opt and g_crr which describe the luminosity evolutions. The procedure for determination of these functions is similar to the ones for removing the correlations between the luminosities except now we make coordinate transformations in the L_opt − z and L_crr − z spaces. We assume simple forms

  $$\begin{equation} g_i(z)=(1+z)^{k_i} \end{equation} \tag{ 5 }$$

  so that L'_i = L_i/g_i(z) refer to the local (z = 0) luminosities.¹¹ The full procedure is detailed in Section 4.
• 3.  
  The density evolution function ρ(z) is determined by the method shown in EP (see Section 5 below). Once all correlations are removed we end up with a local separable LF as in Equation (4).
• 4.  
  The local LFs of uncorrelated luminosities L'_opt and L'_crr can then be used to recover the local radio LF by a straight-forward integration over L'_crr and the true local optical LF as

  $$\begin{eqnarray} \psi _{{\rm rad}}(L_{{\rm rad}}^{\prime }) &=& \int _0^{\infty } { \psi _{{\rm opt}}(L_{{\rm opt}}^{\prime }) \, \psi _{{\rm crr}}\left({{ L_{{\rm rad}}^{\prime } } \over ({{L_{{\rm opt}}^{\prime }/L_{{\rm fid}}})^{\alpha } } } \right)} \nonumber\\ &&\times\, {{{dL_{{\rm opt}}^{\prime }} \over ({{L_{{\rm opt}}^{\prime }/L_{{\rm fid}}})^{\alpha }} } \,}. \end{eqnarray} \tag{ 6 }$$

  As stated above this procedure can be used for the determination of the radio LF at any redshift, from which one can deduce that the radio luminosities also undergo luminosity evolution with

  $$\begin{equation} g_{{\rm rad}}(z) = g_{{\rm crr}}(z) \, \times \, [g_{{\rm opt}}(z)]^{\alpha } \end{equation} \tag{ 7 }$$

  (cf. Equation (3)).
• 5.  
  Similarly we can determine the local distribution of the radio to optical luminosity ratio, R' = L'_rad/L'_opt = L'_crr × L'_opt^{α − 1} × L_fid^−α, as

  $$\begin{equation} G_{R^{\prime }} = \int _0^{\infty } { \psi _{{\rm opt}}(L_{{\rm opt}}^{\prime }) \, \psi _{{\rm crr}}\left({{ R^{\prime } \, L_{{\rm fid}}} \over ({{L_{{\rm opt}}^{\prime }/L_{{\rm fid}}})^{\alpha -1} } } \right) \, {{dL_{{\rm opt}}^{\prime }} \over {{L_{{\rm opt}}^{\prime }}^{\alpha - 1} \, L_{{\rm fid}}} } \,} \end{equation} \tag{ 8 }$$

  and its evolution

  $$\begin{equation} g_{R}(z) = g_{{\rm crr}}(z) \, \times \, [g_{{\rm opt}}(z)]^{\alpha -1} = {g_{{\rm rad}} \over g_{{\rm opt}}}. \end{equation} \tag{ 9 }$$

We now describe results obtained from the use of the procedures described in Section 3 on the data described in Section 2. Here we first give a brief summary of the algebra involved in the EP method. We follow closely the steps described in Maloney & Petrosian (1999). This method uses the Spearman rank test to determine the best-fit values of parameters describing the correlation functions using the test statistic

$$\begin{equation} \tau = {{\sum _{j}{(\mathcal {R}_j-\mathcal {E}_j)}} \over {\sqrt{\sum _j{\mathcal {V}_j}}}} \end{equation} \tag{ 10 }$$

to test the independence of two variables in a data set, say (x_j, y_j) for j = 1, ..., n. Here $\mathcal {R}_{j}$ is the y rank of the data point j in a set associated with it. For untruncated data (i.e., data truncated parallel to the axes), the associated set of point j includes all of the x_k < x_j. If the data are truncated one must form the associated set consisting only of those points of lower x value that would have been observed if they were at the x value of point j given the truncation. As an example, if we have one-sided truncations as in Figures 1 and 2, then the associated set A_j = { k: y_k > y_j, y⁻_k < y_j }, where y⁻_k is the limiting y value of data point i (see EP for a full discussion of this method).

If (x_j, y_j) were independent then the rank $\mathcal {R}_j$ should be distributed uniformly between 0 and 1 with the expectation value and variance $\mathcal {E}_j=(1/2)(j+1)$ and $\mathcal {V}_j=(1/12)(j^{2}+1)$, respectively. Independence is rejected at the n σ level if | τ | > n. To find the best-fit correlation the y data are then adjusted by defining y'_j = y_j/F(x_j) and the rank test is repeated, with different values of parameters of the function F.

4.1. Radio–Optical Luminosity Correlation

The radio and optical luminosities are obtained from radio and optical fluxes from a two-flux-limited sample so that the data points in the two-dimensional flux space are truncated parallel to the axes which we consider to be untruncated. Since the two luminosities have essentially the same relationship with their respective fluxes, except for a minor difference in the K-correction terms, we can consider the luminosity data to also be untruncated. In that case, as mentioned above, the determination of the associated set is trivial and one is dealing with the standard Spearman rank test. Assuming the correlation function between the luminosities F(x) = x^α, we calculate the test statistic τ as a function of α. Figure 4 shows the absolute value of τ versus α, from which we get the best-fit value of α = 1.3 with 1σ range ±0.2. As expected, α is near unity and the value α = 1 is not ruled out with a high significance. As discussed in Appendix B, this correlation may be inherent in the population or may be an artifact of the flux limits and wide range of redshifts, although this particular point is not important for the analysis going forward, as the rotation to L_crr is a technique to achieve independent variables (L_opt and L_crr) in the context of the data present to recover the inherent redshift evolutions.

Zoom In Zoom Out Reset image size

4.2. Luminosity–Redshift Correlations

We now describe our results on determination of the luminosity evolution, i.e., the luminosity–redshift correlation functions g_i(z), which according to Equation (5) reduces to determination of the values of the indices k_i. The basic method for determining the best-fit k_i is the same as above but in this case the procedure is more complicated for several reasons. First, as evident from Figures 1 and 2, the L_i − z data are heavily truncated due to the flux limits. Second, we now are dealing with a three-dimensional distribution (L_crr, L_opt, z) and two correlation functions (g_crr(z) and g_opt(z)).

Specifically, since we have two criteria for truncation, the associated set for each object includes only those objects that are sufficiently luminous in both bands to exceed both flux minima for inclusion in the survey if they were located at the redshift of the object in question. Consequently, we have a two-dimensional minimization problem, because both the optical and correlation-reduced radio evolution factors, $g_{{\rm opt}}(z)=(1+z)^{k_{{\rm opt}}}$ and $g_{{\rm crr}}(z)=(1+z)^{k_{{\rm crr}}}$, come into play, as the luminosity cutoff limits for a given redshift are adjusted by powers of k_opt and k_rad too.

We form a test statistic $\tau _{{\rm comb}} = \sqrt{\vphantom{A^A}\smash{\hbox{${\tau _{{\rm opt}}^2 + \tau _{{\rm crr}}^2}$}}}$, where τ_opt and τ_crr are those evaluated considering the objects' optical and correlation-reduced radio luminosities, respectively. The favored values of k_opt and k_crr are those that simultaneously give the lowest τ_comb and, again, we take the 1σ limits as those in which τ_comb < 1. For visualization, Figure 5 shows a surface plot of τ_comb

Zoom In Zoom Out Reset image size

We have verified this method with a simulated Monte Carlo (MC) data set in which objects are distributed in redshift and given randomized luminosities in accordance with set optical and radio evolutions. The algorithm can recover the evolutions correctly provided that they are not wildly different, i.e., one very positive and the other very negative.

Figure 6 shows the best-fit values of k_opt and k_crr and taking the 1σ and 2σ contours. Results are shown for the entire data set taken as a whole and also with the data split into the RL and RQ subsets. The radio luminosity evolution itself can be recovered by Equation (7).

Zoom In Zoom Out Reset image size

Given the tight constraints achieved when the data set is considered as a whole, and the sharp bifurcation when the set is split into the RL and RQ populations, it is evident that splitting the population before determination of the luminosity evolutions introduces a bias into the determinations. This is expected because differing evolutions will have a strong effect on the likelihood that objects at a given redshift will fall on the RL or RQ side of the division according to the standard definition used here, and the data for each set will be artificially truncated along R = 10 as a function of the evolutions.

We see that positive evolution in both radio and optical wavebands is favored. The minimum value of τ_comb favors an optical evolution of k_opt = 3.0 and a radio evolution of k_rad = 5.4, but uncertainty at the 1σ level allows the range of k_opt from 2.5 to 3.25, and k_rad from 5.3 to 5.75. Therefore, we conclude that quasars have undergone a significantly greater radio evolution relative to optical evolution with redshift. In the above analysis, we have assumed sharp truncation boundaries and that the data are complete above the boundaries. As discussed in Section 8, this may not be the case for the FIRST radio data. If an estimate of the uncertainty from the consideration of possible radio incompleteness at faint fluxes is included, the favored optical range enlarges to k_opt from 1.25 to 3.75, with a slightly lower best-fit value of k_opt = 2.0. Due to combined effects on k_opt and k_crr, the value of k_rad is not much affected by allowing for possible radio incompleteness at faint fluxes, perhaps counterintuitively.

Next we determine the density evolution ρ(z). One can define the cumulative density function

$$\begin{equation} \sigma (z) = \int _0^z {\rho (z) \, dz}, \end{equation} \tag{ 11 }$$

which, following Petrosian (1992) based on Lynden-Bell (1971), can be calculated by

$$\begin{equation} \sigma (z) = \prod _{j}{\left(1 + {1 \over m(j)}\right)}, \end{equation} \tag{ 12 }$$

where j runs over all objects with a redshift lower than or equal to z and m(j) is the number of objects with a redshift lower than the redshift of object j which are in object j's associated set. In this case, the associated set is again those objects with sufficient optical and radio luminosity that they would be seen if they were at object j's redshift. The use of only the associated set for each object removes the biases introduced by the data truncation. Then the density evolution ρ(z) is

$$\begin{equation} \rho (z) = {d \sigma (z) \over dz}. \end{equation} \tag{ 13 }$$

However, to determine the density evolution, the previously determined (in Section 4) luminosity evolution must be taken out. Thus, the objects' optical and radio luminosities, as well as the optical and radio luminosity limits for inclusion in the associated set for given redshifts, are scaled by taking out factors of $g_{{\rm opt}}(z)=(1+z)^{k_{{\rm opt}}}$ and $g_{{\rm rad}}(z)=(1+z)^{k_{{\rm rad}}}$, with k_opt and k_rad determined as above.

Figures 7 and 8 show σ(z) and ρ(z) for the objects in the data set. We evaluate and display the density evolution separately for the RL and RQ objects and for the data set as a whole to compare them. It is seen that the two groups, divided in this way, exhibit very similar density evolution. The number density of quasars seems to peak at between redshifts 1 and 1.5, a little earlier than generally thought for the most luminous quasars (e.g., Shaver et al. 1996), and earlier than that found in Richards et al. (2006), but similar to the peak found for less luminous quasars by Hopkins et al. (2007), and in agreement with Maloney & Petrosian (1999).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Knowing both the luminosity evolutions g_i(z), and the density evolution ρ(z), one can form the luminosity density functions $\pound_{i}(z)$, which are the total rate of production of energy of quasars as a function of redshift. We show this for the radio luminosity density $\pound_{{\rm rad}}(z)$. As evident, the two populations of RL and RQ quasars have very similarly shaped radio luminosity density functions (Figure 9).

Zoom In Zoom Out Reset image size

6.1. General Considerations

In a parallel procedure, we can use the "local" (redshift evolution taken out or "de-evolved") luminosity L'_i distributions (and de-evolved luminosity thresholds) to determine the "local" LFs ψ_i(L_i'), where again the i represents the waveband, and the prime indicates that the luminosity evolution has been taken out. We first obtain the cumulative LF

$$\begin{equation} \Phi _i(L_{i}^{\prime }) = \int _{L_{i}^{\prime }}^{\infty } {\psi _i(L_{i}^{\prime \prime }) \, dL_{i}^{\prime \prime }} \end{equation} \tag{ 14 }$$

and, following Petrosian (1992), Φ_i(L'_i) can be calculated by

$$\begin{equation} \Phi _i(L_{i}^{\prime }) = \prod _{k}{\left(1 + {1 \over n(k)}\right)}, \end{equation} \tag{ 15 }$$

where k runs over all objects with a luminosity greater than or equal to L_i and n(k) is the number of objects with a luminosity higher than the luminosity of object k which are in object k's associated set, determined in the same manner as above. The LF ψ_i(L'_i) is

$$\begin{equation} \psi _i(L_{i}^{\prime }) = - {d \Phi _i(L_{i}^{\prime }) \over dL_{i}^{\prime }}. \end{equation} \tag{ 16 }$$

In Section 4, we have determined the luminosity evolution for two independent functions, the optical luminosity L_opt and the correlation-reduced radio luminosity L_crr. We can form the local optical ψ_opt(L'_opt) and correlation-reduced radio ψ_crr(L'_crr) LFs straightforwardly, by taking the evolutions out. As before, the objects' luminosities, as well as the luminosity limits for inclusion in the associated set for given redshifts, are scaled by taking out factors of $g_{{\rm crr}}(z)=(1+z)^{k_{{\rm crr}}}$ and $g_{{\rm opt}}(z)=(1+z)^{k_{{\rm opt}}}$, with k_crr and k_opt determined in Section 4. We use the notation L → L' ≡ L/g(z).

6.2. Local Optical Luminosity Function

Figures 10 and 11 show the local cumulative Φ_opt(L'_opt) and differential ψ_opt(L'_opt) local optical LFs of the quasars in the White et al. (2000) data set, while Figure 12 shows the local correlation-reduced radio LF, ψ_crr(L'_crr).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The optical LF shows evidence of a break at 2 × 10³⁰ erg s⁻¹ Hz⁻¹, which was present already in data used in Petrosian (1973). Fitting a broken power law yields slope values of −2.0  ±  0.2 and −3.2  ±  0.2 below and above the break, respectively. If we allow for the possibility of additional uncertainty resulting from the consideration of possible radio incompleteness at faint fluxes (see discussion in Section 8), the range on the power law above the break increases to −2.8 ± 0.4. As the optical LF has been studied extensively in various AGN surveys, we can compare the slope of ψ_opt(L'_opt) obtained here to values reported in the literature. For example, Boyle et al. (2000), using the 2dF optical data set (but with no radio overlap criteria) use a customary broken power-law form for the LF, with values ranging from −1.39 to −3.95 for different realizations, showing reasonable agreement.¹²

6.3. Local Radio Luminosity Function

With ψ_opt(L'_opt) and ψ_crr(L'_crr), we can determine the local radio LF ψ_rad(L'_rad) with Equation (6). Figure 13 shows the local radio LF ψ_rad(L'_rad) reconstructed in this way. It is seen that the local radio LF contains a possible break around 10³¹ erg s⁻¹ Hz⁻¹, with a power-law slope of −1.7 ± 0.1 below the break and −2.4 ± 0.1 above it. These ranges for the power law above the break are increased slightly to −2.2 ± 0.3 if the effects of possible radio incompleteness are included, as in Section 8. The slope above the break seen here is similar to earlier results of Schmidt (1972) and Petrosian (1973) which probed only those luminosities. A more complete comparison can be made with Mauch & Sadler (2007), who form radio LFs of local sources in the Second Incremental Data Release of the 6 Degree Field Galaxy Survey (6dFGS) radio catalog. For the sources they identify as AGNs, they find a break at 3.1 × 10³¹ erg s⁻¹ Hz⁻¹, with slopes of −2.27  ±  0.18 and −1.49  ±  0.04 above and below the break (converting to luminosity units).

Zoom In Zoom Out Reset image size

As stated in the introduction, naively one may expect that because the ratio R is independent of cosmological model and nearly independent of redshift, the raw observed distribution would provide a good representation of the true distribution of this ratio. In Figure 14 we show this raw distribution by the triangles, arrived at by using the raw values of R from the data and forming a distribution in the manner of Equations (15) and (16) with no data truncations. It appears that this naive approach shows a hint of possible bi-modality with ${\rm \log } R\sim 1$ as the dividing value.¹³

Zoom In Zoom Out Reset image size

As discussed in Section 3, we can reconstruct the local distribution of $G_{R^{\prime }}(R^{\prime })$, as in Equation (8), which provides for a more proper accounting of the biases and truncations. The results of this calculation are also shown in Figure 14. The distribution calculated in this way clearly is different to the raw distribution and does not show any apparent bi-modality. There is still a possible feature in the same region (${\rm \log } R\sim 1$) where the raw distribution shows a dip. This feature is of marginal significance and results from the similarly shaped feature in ψ_crr(L'_crr) centered around L'_crr = 10⁻⁸ erg s⁻¹ Hz⁻¹. Even if significant, this change in slope cannot be taken as evidence for two physically distinct populations, but could be a useful point to make an arbitrary division into RL and RQ objects.

We also know that the best-fit redshift evolution of the ratio, given Equation (9), is g_R(z) = (1 + z)^2.6. The change in the distribution of R with increasing redshift is also shown in Figure 14.¹⁴ Another way to look at this is that we have found that the radio luminosity evolves at a different rate than the optical luminosity, with the consequence that their ratio is a function of redshift. The radio loudness of the population increases by a factor of 5 by redshift 1 and by a factor of 28 by redshift 3. This is in disagreement with the result presented by Jiang et al. (2007) who show a decrease in fraction of RL sources with increasing redshift, which could be the case if the radio luminosity were to evolve more slowly than the optical luminosity. They however do not determine individual evolutions or LFs. On the other hand, Miller et al. (1990) have noted that the fraction of RL quasars may increase with redshift, which they attribute to a difference in the evolutions of the two populations (RL and RQ). Donoso et al. (2009) compute radio and optical LFs at different redshifts and reach the same conclusion. Cirasuolo et al. (2006) also find that the RL fraction may modestly increase at high redshift. Although not directly comparable, LaFranca et al. (2010) show a similar evolution for R_x, the ratio of radio to X-ray luminosity, as we show here for R. We note that our results favor one population, in the sense that the distribution of G(R), recovered from considering the data truncations inherent in the survey and correlations between the luminosities, is continuous.

Power-law parameterization. One may raise the concern that the simple power-law parameterization used for the redshift luminosity evolutions (Equation (5)) may not be the most ideal one. In particular, it may not accurately represent the evolutions at the highest redshifts considered here. To check this, we repeat the analysis with a different parameterization for the luminosity evolution which allows for a flattening at higher redshifts,

$$\begin{equation} g_i(z)={ {(1+z)^{k_i}} \over {1+\left({{1+z} \over {4}}\right)^{k_i}} }, \end{equation} \tag{ 17 }$$

where i again represents the optical or correlation-reduced radio luminosity. In this parameterization, the functional form for the radio luminosity evolution g_rad(z) and the evolution of the radio-loudness parameter g_R(z) are lengthier expressions involving k_opt k_crr and α, given Equations (7) and (9).

This alternate parameterization for the evolutions does not appreciably affect the results. The best-fit evolution factors as a function of redshift under the alternate parameterization differ very little from those in the simple parameterization to redshift 3.5 (the highest object in the sample). In a data set with higher redshift objects, the form of the parameterization will be more consequential.

Luminosity-dependent density evolution. Another concern may be that LDDE, which is not considered in the functional forms for the LF used here, may more accurately represent the evolution of the LF. As a check of this effect, we divide the data into high- and low-luminosity halves (cutting on optical luminosity), and check the similarity of the computed density evolutions for the two sets versus that computed assuming the absence of LDDE. Given that an artificial difference is already introduced in the two halves because there are a lack of low-luminosity objects in the high-redshift sample and a lack of high-luminosity redshifts in the low-redshift sample (see Figures 1 and 2), we conclude that the density evolutions determined in this way are sufficiently similar to justify neglecting LDDE.

Optical measurement errors. There is the possibility that errors in the optical magnitudes could lead to a bias which could affect the results. The bias introduced by these errors would be negligible if the log N–log S was flat (i.e., N(> S)∝S⁰ and dN/dS∝S⁻¹). Since the number density of sources increases with decreasing flux, it is more likely that a source will be included than excluded. However, the magnitude of this effect will depend on the faint-end source counts slope and the shallower the slope the smaller the effect. The magnitude of this effect is proportional to [(1 + m_below) δ]², where δ is the fractional error in flux and m_below is the faint-end differential source counts slope (i.e., $dN/dS \propto S^{m_{{\rm below}}}$). For the reported POSS-I magnitude errors of 0.2, for the faint-end magnitudes of ∼17.8, δ is less than 0.2. For m_below ∼ −2 and δ ∼ 0.2 the bias will be less than 4% averaged over the faintest fluxes, which would be manifest in raising the faint-end cumulative source counts slope by an amount still considerably smaller than this, as it is fitted over a larger range of fluxes. Previously, Caditz & Petrosian (1993) investigated this effect for a flux-limited data set and showed that the difference was minor between using a Gaussian distribution of fluxes for each source instead of assuming a well-defined flux, for shallow faint-end source counts such as here.

Radio incompleteness. Lastly, the selection function for the FIRST objects in the White et al. sample used here might not be a sharp Heavyside function at a peak pixel flux density of 1 mJy, but rather smeared out. According to Figure 1 of Jiang et al. (2007), the selection function of FIRST for SDSS optically identified quasars is such that at an integrated flux density of 1 mJy only about 55% of sources are seen, and this number rises to 75% at 1.5 mJy and about 85% at 2 mJy.¹⁵ This particular selection function would likely not be identical for the POSS-I optically identified quasars of the sample we use here. Also, we have considered the sample to be limited by the peak pixel flux (i.e., surface brightness limited) in the radio rather than being limited by the integrated flux, in accordance with the criteria set forward in White et al. So it is difficult to directly compare the potential radio selection function here with the one in Jiang et al.

The way to test the effects on our analysis is to repeat the analysis limiting the sample to a higher radio flux, where the sample would presumably be more complete, and determine the extent to which the calculated parameters change in a systematic way. We have done so with a lower radio flux limit of 2 mJy (486 objects), as opposed to the original 1 mJy (636 objects). The effect, propagated through the analysis, is primarily to extend the 1σ uncertainties on k_opt and k_crr in the direction of lower k_opt (1.25 on the low end) and higher k_crr (3.75 on the high end), and to move the best-fit values to 2 and 2.5, respectively. There is no discernable effect on the value of the correlation parameter α. The modified best-fit values and increased 1σ uncertainty for k_opt and k_crr only slightly alter the 1σ range and the best-fit value of k_rad, since the particular error ellipse shape means that lower k_opt values accompany higher k_crr values. The main effect on physical parameters then is to shift the best-fit and low-end 1σ values for k_opt downward. We also find that there is a negligible effect on the density evolutions, there is a small flattening effect on the high-end optical and radio LF slopes, which are reported in Section 6, and a negligible effect on the shape of the G_R(R) distribution. To the extent that faint flux radio incompleteness is present in the sample considered here, it does not seem to have a large systematic effect on the determination of the parameters in this analysis.

We have used a general and robust method to determine the radio and optical luminosity evolutions simultaneously for the quasars in the White et al. (2000) data set, which combines 1.4 GHz radio and R-band optical data for 636 quasars ranging in redshifts from 0.02 to 3.425 and over seven orders of magnitude in radio loudness. We find that the quasars exhibit more substantial radio evolution than optical evolution with redshift (Section 4.2 and Figure 6). We also show that when divided into RL and RQ sets accordingly to the standard definition (divided by the value of the radio-loudness parameter R = 10), the two sub-populations exhibit similar density evolution. The local optical and radio LFs that we obtain are consistent with previous determinations.

Differences are noted with previous determinations of the radio luminosity evolution of quasars. Willott et al. (2001) also used a power-law parameterization of the radio evolution with redshift for the radio-bright sources they consider. Our result for the radio luminosity evolution, when evaluated for the data set as a whole, is not consistent to within uncertainty with their results (power laws ranging from 3.1 to 3.6). Strazzullo et al. (2010) have recently obtained results with a radio survey to low (13.5 μJy) flux limits, quoting k_rad = 2.7 ± 0.3 with the same parameterization used here for the sub-population that they identify as AGNs. This is also not in agreement with the radio evolution determined in our analysis. In a future paper, we intend to carry out the same analysis using the much larger SDSS sample of quasars used by Jiang et al. (2007).

There has been much discussion as to whether RL and RQ quasars, defined solely by means of the radio-loudness parameter as explained above, constitute a true continuum or two populations that can be said to be distinct in some way (e.g., Kellerman et al. 1989; Ivezić et al. 2002, 2004; Cirasuolo et al. 2003). Our analysis favors the former. First, we found that the division of the quasar population into the two aforementioned classes introduces strong biases into the simultaneous determination of the radio and optical luminosity evolutions (Section 4.2; see also Figure 6). More importantly, as shown in Section 7, forming a distribution in the raw values of the radio-loudness parameter R, without taking into account the biases introduced by the truncations in the data and the correlated luminosity evolutions, results in a shape for the distribution which is very different from the true distribution; e.g., it shows a possible dip in G_R(R) at R ∼ 10, while the true distribution is rather smooth and shows at most a modest feature. Even if the feature is real, it does not suggest any bi-modality, but rather a continuous range of physical properties in a single population. We also find that at higher redshifts and optical luminosities the radio loudness become more pronounced. This is opposite to the trends presented by Jiang et al. (2007), but in agreement with others such as Donoso et al. (2009).

Accessing the unbiased distribution of the radio-loudness parameter for quasar sources is crucial not only for understanding the cosmological evolution of this class of active galaxies, but also for understanding jet launching processes in the vicinities of supermassive black holes. In this context, we note that the observed optical fluxes of quasars are dominated by the emission from an accretion disk accreting at relatively high rates, around 1%–100% of Eddington, and therefore radiating with ≃ 10% efficiency. Hence, the optical luminosity is a very good measure of the total accretion power in quasar sources. On the other hand, the observed radio fluxes of the class of objects discussed are expected to originate in the outflowing magnetized plasma. In particular, the radio emission of quasars is produced predominantly via the synchrotron emission of relativistic well-collimated jets (in the case of very RL sources) or via the cyclotron and/or free–free emission of at most mildly relativistic disk winds (in the case of very RQ nuclei). In both cases, the observed radio luminosities should be considered as proxies for the kinetic luminosities of the outflowing matter. Therefore, the radio-loudness parameter R characterizes the efficiency of the production of jets/outflows for a given accretion power.

The lack of any clear bi-modality centered around R = 10 in the distribution of the radio-loudness parameter for quasars, as advocated here, implies then that there is no critical change in the parameters of the central engine between the RL quasars (those producing extremely powerful relativistic jets) and the RQ ones (those producing only mildly relativistic and uncollimated disk winds). This is a crucial piece of evidence for understanding still debated mechanisms for jet launching in AGNs. Note, for example, that our finding is hardly consistent with the idea that RL quasars possess counterrotating (with respect to the black hole spin) accretion disks, as opposed to RQ quasars with co-rotating disks only. Instead, the nuclei of jetted and non-jetted quasar sources seem to be intrinsically very similar, differing only smoothly and continuously in some particular respects. But it has to be emphasized that here we do not discuss the whole population of AGNs (including, e.g., Seyfert galaxies) but only strictly the quasar population (see the related recent discussion in Sikora et al. 2007, and references therein). We note that our analysis does not address the question of whether there may be a large population of quasars with values of R beyond the range of the present sample; in particular we cannot rule out a significant population with very low values of R (log R < 0.01).

Another (possibly related) result we find is that the radio loudness may increase with increasing optical and radio luminosities, as the best-fit value for the correlation parameter α is 1.3, although a strictly linear correlation is not ruled out with much significance and, furthermore, one may dispute that the correlation is inherent in the population (see the discussion in Appendix B). However, if this super-linear correlation were indeed intrinsically the case, it would imply the existence of some connection between the efficiency of formation of relativistic jets and accretion power, which may in turn depend on the combination of the evolving accretion rate and black hole spin (see in this context Tchekhovskoy et al. 2010), given the observed increase in radio loudness with redshift. Note, however, that even though we have used simple one-parameter functions to describe the emerging correlations, it is possible that some of them are more complex. For example, the correlation index α between the radio and optical luminosities may be close to unity only for low-luminosity objects, but much larger than that for more luminous (and therefore RL) quasars. More data and further analysis are needed to address this and similar issues, which may provide further constraints on theoretical models.

Another application of the analysis presented is related to the understanding of the origin of the cosmic background radiation in the radio frequency regime. In Singal et al. (2010) we estimated the fractional contribution of quasars to the cosmic radio background, assuming the level reported by Fixsen et al. (2011). In general, the flux of the objects fitting the definition of RL is well characterized by current interferometric radio surveys so that their contribution to the total radio background intensity can be estimated to be 15%–25% of the observed value. In the earlier work we also estimated the total contribution to the background of the RQ objects, and found it to be between 1% and 2% for favored models of quasar luminosity evolution. This estimate was based on integrating values of the quasar bolometric LF, as reported in the literature, over redshift, applying a mapping between optical and radio luminosity, and assuming that the optical and radio luminosities had identical redshift evolutions. We also noted there that the contribution we estimated was dependent on the latter assumption and would be revised in the case of differing radio and optical luminosity evolutions. As we see here that the population of quasars has greater radio evolution relative to optical, the contribution of RQ quasars to the radio background will be somewhat larger than the value reported in our previous work. We will present a quantitative determination in a forthcoming paper.

J.S. thanks S. Kahn and R. Schindler for their encouragement and support. Ł.S. was supported by the Polish Ministry of Science and Higher Education through the project N N203 380336.

Distribution of radio loudness, R. As described in the introduction, one of the aims of this paper is the determination of the distribution of the radio to optical luminosity ratio R = L_rad/L_opt from an observed sample of radio and optical fluxes and redshifts. This requires a proper accounting for correlations and evolutions of optical and radio luminosities and the observational selection effects. Here we describe the how these factors affect the observed distribution of R.

The true or intrinsic distribution of R values is related to the radio and optical LFs Ψ(L_opt, L_rad, z) as

$$\begin{eqnarray} G_T(R,z)&=&\int ^\infty _0\Psi (L_{{\rm opt}}, R \,\, L_{{\rm opt}}, z) \, L_{{\rm opt}} \, dL_{{\rm opt}}\nonumber\\ &=&\int ^\infty _0\Psi \left({L_{{\rm rad}} \over R}, L_{{\rm rad}}, z \right) \, L_{{\rm rad}} \, {{ dL_{{\rm rad}} } \over {R^2}} \end{eqnarray} \tag{ A1 }$$

(compare to the separable, local form in Equation (8)).

The observed distribution, on the other hand, is different because the observational selection effects truncate the data. For example, for a sample with well-defined flux limits f_rad ⩾ f_{m, rad}, f_opt ⩾ f_{m, opt}, and R_obs ≡ f_{m, rad}/f_{m, opt} the observed distributions are

$$\begin{eqnarray} &&G_{\rm obs}(R,z) = \int ^\infty _{L_{\rm min,{\rm opt}}(z)}\Psi (L_{{\rm opt}}, R \,\, L_{{\rm opt}}, z) \, L_{{\rm opt}} \, dL_{{\rm opt}}\qquad\\ &&\hspace{-0.9pc} {\rm for}\,\,\,\,\, R>R_{\rm obs}\nonumber\;{\rm and} \end{eqnarray} \tag{ A2 }$$

$$\begin{eqnarray} G_{\rm obs}(R,z)& = & \int ^\infty _{L_{\rm min,rad}(z)}\Psi \left({L_{{\rm rad}} \over R}, L_{{\rm rad}}, z \right) \, L_{{\rm rad}} \, {{ dL_{{\rm rad}} } \over {R^2}} \nonumber \\ & = & \int ^\infty _{L_{\rm min,{\rm opt}}(z)R_{\rm obs}/{\hat{R}}}\Psi (L_{{\rm opt}}, R \,\, L_{{\rm opt}}, z) \, L_{{\rm opt}} \, dL_{{\rm opt}}\qquad\\ &&\hspace{-4.7pc}{\rm for}\,\,\,\,\, R<R_{\rm obs},\nonumber \end{eqnarray} \tag{ A3 }$$

where

$$\begin{eqnarray} &&{\hat{R}}=(K_{{\rm opt}}/K_{{\rm rad}})R,\nonumber\\ && L_{\rm min,{\rm opt}}(z)=4\pi d_L^2(z)K_{{\rm opt}} f_{m,{\rm opt}},\\ && L_{\rm min,rad}(z)=4\pi d_L^2(z)K_{{\rm rad}} f_{m,{\rm rad}},\nonumber \end{eqnarray} \tag{ A4 }$$

and K_i is the K-correction factor for waveband i. All of these are obtained from the observed distributions of the fluxes and redshifts. If we approximate the observed distribution of fluxes and redshifts by a continuous function n_obs(f_opt, f_rad, z) then

$$\begin{eqnarray} \Psi (L_{{\rm opt}}, L_{{\rm rad}}, z)&=&n_{\rm obs}\left({L_{{\rm opt}}\over 4\pi d_L^2(z)K_{{\rm opt}}},{L_{{\rm rad}}\over 4\pi d_L^2(z)K_{{\rm rad}}},z\right)\nonumber\\ &&\times\,\left({1\over 4\pi d_L^2}\right)^{2} {1\over K_{{\rm opt}}K_{{\rm rad}}V^{\prime }(z)}, \end{eqnarray} \tag{ A5 }$$

so that

$$\begin{eqnarray} &&G_T(R,z) = (K_{{\rm opt}}K_{{\rm rad}}V^{\prime })^{-1}\int ^\infty _0 n_{\rm obs}(f_{{\rm opt}},{\hat{R}}f_{{\rm opt}},z) \, f_{{\rm opt}} \, df_{{\rm opt}},\nonumber\\ \end{eqnarray} \tag{ A6 }$$

$$\begin{eqnarray} &&G_{\rm obs}(R,z) = (K_{{\rm opt}}K_{{\rm rad}}V^{\prime })^{-1}\int ^\infty _{f_{{\rm lim}}} n_{\rm obs}(f_{{\rm opt}},{\hat{R}}f_{{\rm opt}},z) \, f_{{\rm opt}} \, df_{{\rm opt}},\nonumber\\ \end{eqnarray} \tag{ A7 }$$

where V' = dV(z)/dz, and f_lim = f_{m, opt} for ${\hat{R}}>R_{\rm obs}$ and $f_{{\rm lim}}=f_{m,{\rm opt}}(R_{\rm obs}/{\hat{R}})$ for ${\hat{R}}<R_{\rm obs}$. Note that ${\hat{R}}$ depends on redshift to the extent that the optical and radio K-corrections are different. However, since the radio and optical spectra can be well approximated by power laws with almost the same spectral index the two K-corrections are almost equal. In what follows we ignore this small difference and set ${\hat{R}}=R$. The sources with R > R_obs can be called optically limited because their optical flux is the main determining factor for their inclusion in the sample. Similarly sources with R < R_obs can be called radio limited.

Clearly the intrinsic and observed distributions of R are different. In reality the situation is slightly more complicated because the lower limit of the integration does not extend to zero. The samples of available quasars are also truncated by minimum luminosities, say L_{m, opt} and L_{m, rad}, which introduce a second critical value for R, namely R_int ≡ L_{m, rad}/L_{m, opt}. The above equations are valid for redshifts $z>z_{\rm min,{\rm opt}}$ or z_{min, rad} defined as

$$\begin{equation} L_{m,i}=4\pi d_L^2(z_{\rm min,i})K_i(z_{\rm min,i})f_{m,i}. \end{equation} \tag{ A8 }$$

For z < z_min (defined as the lower of $z_{\rm min,{\rm opt}}$ and z_{min, rad}) there is no truncation due to flux limits and G_obs = G_T.

It is convenient to define

$$\begin{equation} \Phi (R,z;x)=\int ^\infty _i\Psi (L_{{\rm opt}},R \,\, L_{{\rm opt}},z) \, L_{{\rm opt}} \, dL_{{\rm opt}}, \end{equation} \tag{ A9 }$$

so that the true distribution can be written as

$$\begin{eqnarray} G_T(R,z) & = & \Phi (R,z;L_{m,{\rm opt}}) \,\,\,\,\,{\rm for}\,\,\,\,\, R>R_{\rm int}, \end{eqnarray} \tag{ A10 }$$

$$\begin{eqnarray} G_T(R,z) & = & \Phi \left(R,z; {{R_{\rm int} \, L_{m,{\rm opt}}} \over {R}} \right) \,\,\,\,\,{\rm for}\,\,\,\,\,R<R_{\rm int}.\nonumber\\ \end{eqnarray} \tag{ A11 }$$

The observed distribution depends on the relative values of R_int and R_obs.

For a sample where R_obs > R_int, which means that $z_{\rm min,rad}>z_{\rm min,{\rm opt}}$, such a sample may be classified as mainly optically selected. (For a purely optically selected sample R_obs = 0.) In this case, the observed distribution will be the same as the true distribution (G_obs = G_T) for all redshifts z < z_{min, rad}. Obviously, the reverse is true in the opposite case with radio exchanged for optical. At higher redshifts

$$\begin{eqnarray} {G_{\rm obs}(R,z)\over G_T} & = & {{ \Phi (R,z ;L_{\rm min,{\rm opt}}(z))} \over {\Phi (R, z; L_{m,{\rm opt}}) }}\,\,\,\,\,{\rm for}\,\,\,\,\, R>R_{\rm int},\nonumber\\ \end{eqnarray} \tag{ A12 }$$

$$\begin{eqnarray} {G_{\rm obs}(R,z)\over G_T} & = & {\Phi \left(R,z ;{R_{\rm obs}\,L_{\rm min,{\rm opt}}(z)\over R} \right) \over \Phi \left(R, z; {R_{\rm int} \, L_{m,{\rm opt}} \over {R}} \right) } \,\,\,\,\,{\rm for}\,\,\,\,\,R<R_{\rm int}.\nonumber\\ \end{eqnarray} \tag{ A13 }$$

A simple example. Let us assume that the radio and optical luminosities are uncorrelated and do not evolve so that we have Ψ(Lopt, Lrad, z) = ψopt(Lopt) ψrad(Lrad) ρ(z), where ρ(z) describes the density evolution. Furthermore if we assume simple power-law LFs ψopt(Lopt) = AoptL−αoptopt and ψrad(Lrad) = AradL−αradrad, it is easy to show that

$$\begin{eqnarray} &&G_{T}(R,z)\propto (R/R_{\rm int})^{1-\alpha _{\rm rad}} \,\,\,\,{\rm for}\,\,\,\, R>R_{\rm int},\nonumber\\\\[-9pt] && G_{T}(R,z)\propto (R/R_{\rm int})^{\alpha _{\rm opt}-1} \,\,\,\,{\rm for}\,\,\,\, R<R_{\rm int}.\nonumber \end{eqnarray} \tag{ A14 }$$

Similarly, it is easy to show that for redshifts z > z_min

$$\begin{eqnarray*} G_{\rm obs}(R,z)&=&G_T{\left(L_{\rm min,{\rm opt}}(z) \over L_{m,{\rm opt}}\right)^\beta }\nonumber\\ && \times\, \left\lbrace \begin{array}{@{}ll}1 &{\rm for}\ R>R_{\rm obs}>R_{\rm int}, \cr (R_{\rm obs}/R)^{\beta } &{\rm for}\ R_{\rm obs}<R<R_{\rm int},\cr (R_{\rm obs}/R_{\rm int})^\beta &{\rm for} \ R_{\rm obs}<R_{\rm int}<R, \end{array}\right. \end{eqnarray*}$$

where β = α_opt + α_rad − 2.

The observed and true distributions have different shapes in the range between the intrinsically and observationally limiting values of R. The shapes become identical when these two values are equal. Larger differences will be the case if the radio and optical luminosities are correlated nonlinearly and undergo different kinds of luminosity evolution. As a result the fraction of RL or RQ sources (arbitrarily chosen at some value or R) will vary with redshift and/or luminosities. For example, if the radio and optical luminosities were correlated linearly making the radio loudness independent of both luminosities then we will be dealing with a separable LF ψ(R, L_opt, z) = G(R)ψ_opt(L_opt, z). In such a case, because of flux limits at any redshift the observed range extends to R_min(L_opt, z) which decreases with increasing L_opt but increases with increasing z. This will cause the fraction of RL sources to decrease with luminosity but increase with redshift. All such trends can be determined by proper accounting of the correlations and evolution as described in this paper.

Luminosity correlations and flux limits. As pointed out by the referee, it has been suggested by Antonucci (2011) that the observed correlation between the luminosities may be induced by observational selection effects. Other works that have examined the issue include Khembavi et al. (1986), Feigelson & Berg (1983), and Chanan (1983). It is not clear whether the correlation seen in this work between L_rad and L_opt is inherent in the quasar population or is introduced by the selection effects of the surveys. In some sense this question does not matter for the analysis presented here, because the rotation to L_crr is a technique required to achieve independent variables (L_opt and L_crr) in the context of the data present so that we can recover the inherent redshift evolutions. This is independent of the underlying luminosity–luminosity correlation of the population. However, it is of general interest whether the correlation seen between radio and optical luminosity is inherent.

In order to begin an investigation of this, we simulate via MC techniques two cases of quasar populations, distributing the objects in radio luminosity, optical luminosity, and redshift, and applying the flux limits of the White et al. data set to achieve for both MC populations "observed" data sets of a similar size to the real White et al. data used in this analysis. For the MC populations, we use a redshift distribution with ρ(z)∝z⁴ to redshift 1.2 and then constant density with redshift for 1.2 ⩽ z ⩽ 3.2. In the first, uncorrelated, case, we draw the optical and radio luminosities of the population from separate distributions with power-law slopes in luminosity of −2. In the second, correlated, case, we draw the optical luminosities of the population from a distribution with a power-law slope in luminosity of −2, then the radio luminosities of the population are assigned according to L_rad∝L_opt^1.3 and then randomized about that value by a factor of 10^{2.5 × A}, where A is a normally distributed random number with mean of zero and standard deviation of one. For this simple analysis we have not included the effects of luminosity evolution, because it would require an orders of magnitude larger simulation.

Our results show that that indeed the radio and optical luminosities of the "observed" data sets in both cases are correlated. An analysis identical to that presented in Section 4.1 reveals that the observed data for the correlated MC population have a value for the correlation index α of 1.3 ± 0.3 while those of the uncorrelated MC population have a value of 0.85 ± 0.15. However, these "observed" sets can be compared to the real White et. al data set, and it is seen that the observed data set for the correlated MC population more closely resembles the real data set than that of the uncorrelated MC population does, for instance comparing the redshift distributions of the observed objects, shown in Figure 15. From this analysis, we conclude that the correlation seen in the White et al. data set between L_rad and L_opt could be inherent in the population. Exact determinations of what fraction of the observed correlations is inherent in the population and what fraction is due to selection effects depend on the values of the many parameters describing the LFs and evolutions. This is an important issue in many areas of astrophysics and requires considerable work which is beyond the scope of this paper. In a forthcoming work we will address this general question.

Zoom In Zoom Out Reset image size