© ESO, 2011

1. Introduction

At high redshift, most of the baryons are located in the intergalactic medium (IGM; e.g. Petitjean et al. 1993) where they are highly ionized by the UV-background produced by galaxies and QSOs (Gunn & Peterson 1965), at least from z ~ 6 (Fan et al. 2006; Becker et al. 2007). The large absorption cross section of the H i Lyman-α transition implies that the small fraction of neutral hydrogen in the IGM produces the so-called Lyman-α forest, which is composed of numerous absorption lines detected in the spectra of high-redshift quasars (see Lynds 1971; Rauch 1998,for a review).

Analytical models (Bi et al. 1992) and numerical N-body simulations (Cen et al. 1994; Petitjean et al. 1995; Zhang et al. 1995; Hernquist et al. 1996; Theuns et al. 1998; Riediger et al. 1998) have been very successful at reproducing the properties of the Lyman-α forest as measured from high spectral resolution, high SNR data obtained with the Ultraviolet and Visual Echelle Spectrograph (UVES) on the Very Large Telescope (VLT, e.g. Bergeron et al. 2004; Kim et al. 2007) and HIRES on the Keck telescope (e.g. Hu et al. 1995). The overall picture shows that lower column-density H i absorption lines trace the filaments of the “cosmic web”, and higher column-density absorption lines trace the surroundings of galaxies. Detailed studies of absorption line properties and of their clustering properties along one or several adjacent lines of sight give additional constraints on the ionization history, correlation length, matter power spectrum etc. (see e.g. Petitjean et al. 1998; Croft et al. 1998; McDonald et al. 2005; Theuns & Srianand 2006). The next generation of quasar surveys from BOSS (SDSS-III, Schlegel et al. 2007; Eisenstein et al. 2011) to Big-BOSS (Schlegel et al. 2009) should provide the first detection of baryonic acoustic oscillations in the IGM at z ~ 2−3 (Slosar et al. 2009; White et al. 2010).

An important quantity to measure in a quasar spectrum is the mean amount of absorption in the Lyman-α forest, D_A, defined as: D_A  =  1 − ⟨F⟩ (Oke & Korycansky 1982), where F is the quasar normalized flux, F = F_obs/F_cont, F_obs is the observed flux and F_cont is the estimated unabsorbed continuum flux. The absorption can be defined by the mean effective optical depth as well, τ_eff  =  −ln⟨F⟩. These quantities are sensitive to the physical properties of the IGM and have been used to constrain Ω_b (Rauch 1998; Tytler et al. 2004), the ionization history (Rauch et al. 1997; Kirkman et al. 2005; Bolton et al. 2005; Bolton & Haehnelt 2007; Prochaska et al. 2009), and in particular the He ii reionization (Bernardi et al. 2003; Theuns et al. 2002). The latter could possibly induce a dip in the evolution with redshift of τ_eff at z ~ 3.2 (Schaye et al. 2000).

Bernardi et al. (2003) first discovered such a dip in the evolution of the effective optical depth in the Lyman-α forest using SDSS spectra. The existence of a feature was later confirmed from high-resolution studies (Faucher-Giguère et al. 2008; Dall’Aglio et al. 2008) at a more modest statistical significance but at a coincident redshift.

Intermediate resolution data have been used as well (McDonald et al. 2005; Dall’Aglio et al. 2009) but the feature was not detected. However, the methods used may not be entirely appropriate and we come back to this point in the present paper. Note that Faucher-Giguère et al. (2008) cautioned that (i) the dip interpretation is only valid if one insists on fitting a single power-law to the background evolution, for which there is no clear physical motivation and (ii) this feature is not necessarily caused by He ii reionization, and other interpretations are also possible.

The definition of the unabsorbed quasar continuum over the Lyman-α forest is a critical point (e.g. Tytler et al. 2004; Kim et al. 2007). For low-resolution spectra, most analyses first define a continuum redwards of the QSO Lyman-α emission, where there are only few absorption lines, and extrapolate the shape of the continuum in the Lyman-α forest region (see Sect. 2.1 for more details). It is mostly commonly assumed that the QSO continuum in regions without an emission line is a power-law that can be easily extrapolated. However, this assumption usually neglects weak emission lines both in the red and, more importantly, within the Lyman-α forest region. Because of this, Suzuki et al. (2005, S05) have applied a principal component analysis (PCA) to HST spectra of quasars at z ≤ 1. Quasar continua are described with a limited set of eigenvectors and a controlled sample is used to define a projection matrix that allows one to recover the continuum in the Lyman-α forest from the shape of the continuum in the red part.

The shape of the quasar continuum can evolve from z ≤ 1 to z ~ 3. In this case, a PCA at z ≤ 1 would not give a fair representation of quasar continuum at z ~ 2−3. In this paper, after describing the procedures in Sect. 2, we take advantage of the large database provided by SDSS-DR7 to define a large enough sample of quasars at z ~ 3 to which we can apply the same procedure as in S05 (Sect. 3). New eigenvectors and a projection matrix are generated and then used to predict the continuum of all SDSS-DR7 spectra. We apply the method to determine the evolution with redshift of the mean flux in the Lyman-α forest (Sect. 4) and discuss the significance of the bump at z ~ 3.2−3.4 before drawing our conclusions in Sect. 5.

2. Procedures

2.1. Different methods to estimate the QSO continuum in low-resolution spectra

The methods used to estimate the continuum in low-resolution spectra can be broadly classified as below:

1. A direct estimate of the continuum in the Lyman-α forest region:

   • Using a spline interpolation: a cubicspline is interpolated on adaptative intervals betweenobserved data points in the forest to construct a lo-cal continuum. A correction is thenapplied to take into account that these data pointscan be affected by some absorption. Dall’Aglioet al. (2008, 2009)applied a systematic correction that accounts for resolutioneffects and line blending and is estimated from idealizedMonte-Carlo simulated spectra. This approach, however,neglects the possibility of continuous absorption from thesmooth IGM (rather than discrete absorbers) as well as cor-relations from large-scale structure. At highredshift, continuous absorption can be important and maycause the true continuum to be underestimated even afterapplying the Monte-Carlo method to high-resolution data.Faucher-Giguère et al. (2008)developed an alternative method to correct the continuumplacement using cosmological simulations and showed thatthis effect is indeed important at the >10% level at z  =  4.

   • Taking into account the difference between the continuum and absorption wavelength dependencies (e.g. Bernardi et al. 2003; Prochaska et al. 2009): the continuum is a property of the quasar and depends only to first order on the restframe wavelength (λ_r = λ/(1 + z_em)); while the absorption depends on the redshift (z_abs = λ/1215.6701 − 1) only. Thus, if one separates the dependencies in the flux, F(λ_r,z) = C(λ_r)exp(−τ(z)), both quantities can be recovered in principle.

2. Using the red part of the QSO spectrum to predict the blue part:

   • A power-law is adjusted to the red part of the spectrum inregions free of emission and absorption lines(see Sect. 2.2). The power-lawis then simply extrapolated over the Lyman-α forest wavelength range. This procedure does not account for weak emission lines in the Lyman-α forest region.

   • A PCA applied to a reference sample describes the continuum of a quasar spectrum as a linear combination of eigenvectors. A projection matrix is generated and used to translate the weights of the eigenvectors describing the red side of the spectrum into the weights of the eigenvectors for the whole spectrum.

We focus on the two last methods in which information provided by absorption-free regions redwards of the Lyman-α emission line are used to predict the continuum in the forest. Both methods require first to estimate the quasar emission redshift. Our procedure to estimate this redshift and to provide a power-law continuum is described in the next subsection. We then describe in Sect. 2.3 the PCA method.

2.2. Determination of the position of the emission lines and the power-law continuum

We derive a redshift of the quasar using C iv and C iii] emission lines to be able to estimate the position of these lines and to avoid them when fitting the power-law. This is therefore not an attempt to derive the exact systemic quasar redshift, which is known to be shifted compared to the C iv-C iii] redshift (e.g. Vanden Berk et al. 2001; Hennawi & Prochaska 2007).

Here we assume that the continuum redwards of the Lyman-α emission line can be described as the sum of a power-law component and a Gaussian function for each of C iv and C iii] emission lines. The redshift and the power-law component are estimated as follows (see Fig. 1 for a typical example at z_em = 2.862):

• 1.

  After convolving the whole spectrum with a Gaussian filter offixed FWHM = 250   km   s^-1, the position where the flux is maximum is associated to the QSO Lyman-α emission line. This gives a first rough estimate of the redshift, z₁.

• 2.

  The average of the positions of the maximum flux within a window of 100 Å in the restframe of the quasar (at z = z₁) around C iv and C iii] emission lines provides a second redshift estimate, z₂. This redshift should be more accurate than z₁ because the peak of the Lyman-α emission is a poor estimate of the redshift owing to the Lyman-α forest and the blending with N vλ1240.

• 3.

  A power-law component of the continuum is fitted using windows devoid of emission lines between 1430−1500 Å, 1600−1830 Å and 2000−2500 Å in the restframe (see grey windows in Fig. 1). This component is subtracted from the spectrum.

• 4.

  Finally, C iv and C iii] emission lines are simultaneously fitted with Gaussian functions. The width and amplitude are independent parameters, but the two Gaussian functions are bound to have the same redshift, z₃.

The red line in Fig. 1 is the sum of the power-law and the two emission lines. The extrapolation of the power-law component bluewards of the Lyman-α emission provides a first estimate of the quasar continuum (red line in Fig. 1).

2.3. Principal component analysis

We summarize here the main steps of the method as described in Francis et al. (1992) and Suzuki et al. (2005).

2.3.1. Reconstructed continuum

A representative sample of quasar spectra at the redshift of interest must be gathered for which it is possible to define a true continuum, q(λ), i.e. unspoiled by intervening absorption. Suzuki et al. (2005) used HST spectra at z   ≤ 1 because the IGM is sparse at these redshifts and thus the continuum can be easily interpolated above absorption lines. The sample of SDSS-DR7 quasars we used has a mean emission redshift of z ~ 2.9 and is defined in Sect. 3. We derived the true continuum, q(λ), by eye and used these fitted quasar continua in the following.

A covariance matrix V is first calculated for the N QSOs in the sample as (1)where is the mean quasar continuum.

The principal components are found by decomposing the covariance matrix V into the product of the orthonormal matrix P, which is composed of eigenvectors, and the diagonal matrix Λ containing the eigenvalues: (2)We call the eigenvectors (i.e., the columns of the matrix P) the principal components, ξ_j. The principal components are ordered according to the amount of variance in the training set they can accommodate, such that the first principal component is the eigenvector that has the highest eigenvalue.

The distribution of the weights, c_j, of the jth principal component in Eq. (4) is found from the distribution of the c_ij for all i = 1..N QSOs of the sample: (3)Note that the upper limit of the integration is higher here than that of S05. This is discussed further in Sect. 3. A mock continuum can now be constructed using Eq. (4) over the rest frame wavelength range of the spectra in the sample. (4)

2.3.2. Predicted continuum

The goal is to quantify the relationship between the red and blue sides of the spectra in the sample. The first m principal components ξ_j(λ) and their weights, c_ij, are derived as described above, using the whole rest wavelength range, 1020 to 2000 Å. Another set of m principal components, ζ_j(λ) and their weights, d_ij, are defined using only the red rest wavelength range, 1216 to 2000 Å. Finally, we solve linear equations to find a projection matrix relating c_ij and d_ij. Weights can be written in the N  ×  m matrix form C  =  c_ij and similarly for D. We then use singular value decomposition techniques (Press et al. 1992) to derive the m  ×  m projection matrix X  =  x_ij translating weights found with the red side only into the weights for the whole spectrum: (5)Once matrix X is known, we can estimate the continuum over the Lyman-α forest for any quasar spectrum from the red part of the spectrum. We proceed in three steps. The weights for the red spectrum are found, (6)The weights from the red side b_j are translated to weights for the whole spectrum, using (7)Then the continuum for the whole spectrum is built as (8)

3. New PCA continuum at z  ~  3

The eigenvectors and coefficients as derived by S05 at low redshift were used to generate mock spectra to test different analyses at high redshift (Dall’Aglio et al. 2008; Kirkman et al. 2005). To our knowledge, one attempt has also been made to derive a similar decomposition at high redshift using SDSS spectra by McDonald et al. (2005). These authors comment that this continuum determination is robust enough to infer the mean flux evolution, but unstable as far as the Lyman-α power spectrum is concerned. However, the decomposition is not performed on a well controlled training sample as in S05 and in the present work (see below), therefore, they do not provide a projection matrix. The analysis performed by Yip et al. (2004) is closer to our purpose. They applied a PCA to the 16 707 Sloan Digital Sky Survey DR1 quasar spectra (0.08 < z < 5.41) and reported that the spectral classification depends on redshift and luminosity. No compact set of eigenspectra succeeds in describing the variations observed over the whole redshift range. Besides, it seems that there is a differential evolution with redshift of the coefficients. Since these authors are interested in the quasar continuum only, they do not try to recover the exact continuum over the Lyman-α forest and consider only the observed flux (that is, continuum plus absorption).

Owing to the large difference of redshift between the quasars used in S05 and those involved in any Lyman-α forest study using SDSS spectra, one may wonder if there is any evolution in the continuum of quasars or any change in the correlation between the shapes of the continuum over the Lyman-α forest and redwards of the Lyman-α emission line compared to what is found by S05. Deriving new components at redshift 3 should answer this question.

The other motivation for this work is to provide principal components and distributions of coefficients over a larger wavelength coverage than in previous studies: the S05 matrix allows one to generate continua from Lyman-β to C iv emission lines, while in the present work we will extend the wavelength coverage beyond the C iii] emission line (until 2000 Å in the restframe). This should in principle facilitate the extrapolation in the blue.

3.1. Deriving new principal components from a sub-sample of z  ~  3 SDSS-DR7 quasar spectra

The difficulty of this analysis is to estimate the continuum in the Lyman-α forest, where absorption can be neither neglected nor easily removed because of the low resolution of the SDSS spectra. In particular it would be very difficult to define the true continuum automatically.

We first selected spectra with a signal-to-noise ratio per pixel greater than 14 redwards of the Lyman-α emission (the SNR is computed around 1280 Å in the restframe). We require the redshift of the quasars to be greater than 2.82 and lower than 3.00. The lower limit is chosen as such for two reasons: the Lyman-α forest has to be complete (z > 2.7) and we noticed that there are some problems with the flux calibration at the very blue end of the SDSS spectra, which is why we would like to avoid this part of the spectrum. To illustrate these problems and estimate the exact wavelength at which to start the study, we selected spectra where a damped Lyman-α system (DLA) is observed with aborption redshift greater than 3.7 and with a column density N(H i)  ≥  10^20.5   cm^-2 (the list is available in Noterdaeme et al. 2009). In these spectra, and because of the presence of the DLA, the flux is expected to be equal to zero for λ_obs  ≤  4280 Å. When stacking the selected lines of sight (Fig. 2), we note instead that the flux increases for wavelength lower than 4000 Å. The difference with zero is as high as 0.05 for a normalized spectrum, meaning that this part of the spectrum should probably not be used for the analysis. Because the study of the forest is usually limited to beyond the O vi emission line (λ_rest > 1050 Å), the minimum emission redshift will be 2.82. The upper limit is a compromise between the number of spectra needed for the analysis and our ability to estimate confidently a continuum by eye: it has been set to z = 3. BAL quasars, lines of sight containing a DLA or any spectrum for which fitting a continuum is too risky because of missing pixels or reduction issues are removed from the analysis. The SDSS spectra are observed with two cameras, one for the blue and one for the red part of the spectrum and we were concerned about the presence of a possible discontinuity or break at the merging point. We avoided any spectra that were possibly affected by this. After we applied all those constraints, 78 spectra remained in the training set (Table 1).

Once the training set was defined, spectra were smoothed and the continuum was “hand-fitted”. Redwards of the Lyman-α emission line, we followed the different emission lines and ignored isolated absorption lines. In the Lyman-α forest, the continuum cannot be uniquely defined. We assumed that the continuum has a smooth shape, that it roughly follows the peaks of the spectrum, and that the blending of lines at the SDSS resolution is large. Points are placed around the peaks of the flux, and a spline interpolation is used to connect them. After a first try, we minimized the number of points used and checked that the continuum is indeed located above the blends of lines, but at about the level of “flat” regions. A typical example is given in Fig. 3. Note that with this procedure we de facto take into account that the chosen points can be affected by some absorption. To check that this is indeed the case, the spectra were then rebinned with 0.5 Å restframe pixels and the mean flux evolution was computed and compared to Faucher-Giguère et al. (2008) measurements from high and medium resolution and high signal-to-noise spectra. Both estimates agree which gives us confidence in our hand-fitted continua (Fig. 4).

The procedure described in Sect. 2.3 was then performed and the correlation matrix was computed and is displayed in Fig. 5. In agreement with S05, a moderate correlation (0.3−0.6) is found between the shape of the continuum in the forest and the region between Lyman-α and C iv emission lines. Thanks to the larger restframe wavelength coverage of this study, a moderate anti-correlation (from −0.6 to −0.4) between the shape of the continuum in the forest and in the region between C iv and C iii] emission lines is found. Suzuki et al. (2005) noticed that a PCA continuum has the good shape in the forest but that the amplitude of the power-law component is unstable. The anti-correlation in that extrapart of the continuum may improve the stability of the prediction of the continuum amplitude in the Lyman-α forest. This is discussed in more details in Sect. 3.2.1.

The mean continuum is shown in Fig. 6 together with the mean continuum of S05 (z  ≤  1) and the composite spectrum derived by Vanden Berk et al. (2001). In the wavelength range of interest here, quasars contributing to the Vanden Berk et al. (2001) composite cover a redshift range from 2.13 to 4.789 for the Lyman-α region and from 1.5 to 4.789 for C iv. Our mean continuum excellently agrees redwards of the Lyman-α emission line with the Vanden Berk et al. (2001) composite. The discrepancy in the blue is simply caused by the absorption in the Lyman-α forest that Vanden Berk et al. (2001) did not try to remove. When comparing to S05, one can see that the amplitude of C iv, Lyman-α and Lyman-β emission lines relative to the continuum are less important in the SDSS spectra. While the determination of the C iv emission is relatively straightforward, the presence of absorption at the position of the Lyman-β emission line and in the blue wing of the Lyman-α emission line makes the continuum difficult to estimate. Thus, the main and robust difference between the mean continua in SDSS and HST spectra is the variation of the C iv equivalent width. This evolution of the QSO emission line equivalent widths has been noted for a long time (Baldwin 1977) and has also been reported by Zheng et al. (1997). The latter authors used 101 HST spectra to compute a low-z composite spectrum (90% of the quasars had a redshift lower than 1.5) and compared it to the Francis et al. (1991) composite (z ~ 3). This evolution is probably related to the quasar luminosity. Note that no evolution is found by Fan (2009) from z ~ 2 to z  >  6.

The first ten eigenvectors are displayed in Fig. 7 when derived from the full wavelength coverage and in Fig. 8 when derived from the region redwards of the Lyman-α emission together with the components provided by S05. The distributions of the coefficients, c_i,j and d_i,j, computed on our sample using Eq. (3), are also shown. The first component looks very similar in the two decompositions and is dominated by the amplitude of the Lyman-α and C iv emission lines. The most significant difference with S05 lies in the shape of the distribution of the associated coefficients: in their study, the distribution is Gaussian whereas in ours this distribution is log-normal. This trend (concerning the eigenvectors and the distributions of coefficients) agrees with what has been found by Francis et al. (1992). The discrepancy between the shapes of the coefficient distributions could mean that the S05 sample is more homogeneous than ours in terms of the emission line amplitudes.

The second component is dominated by the continuum slope and it seems that there is a difference between what is found here and in S05, which seems to be somewhat compensated by the difference in the third component. Other components show small differences but they are less pronounced and the coefficient distributions are very similar.

To estimate more quantitatively the similarity of the low and high redshift sets of eigenspectra, we follow Yip et al. (2004) and compute the sum of the projection operators of each set of eigenvectors |ξ_j > : (9)and then the trace of the products of the projection operators: (10)where D will be the common dimension of both sets. The two sets are disjoint if the trace is zero. If the bases are completely alike, D should be equal to their dimension, therefore D = 10 in our case. To compute this number, we cut our eigenspectra at 1600 Å (rest) and we find D = 7.6. This means that the two decompositions are similar but not exactly the same, confirming the slight evolution of the decomposition with redshift.

We provide the first 10 eigenvectors of the PCA in an electronic form. The distributions of associated coefficients are very close to Gaussian functions, except for the first coefficient, the distribution of which is fitted with a log-normal distribution. Their characteristics are listed in Table 2.

3.2. Quality of the predicted continuum

The decomposition of the quasar emission investigated at z ≤ 1 by S05 with HST spectra and at z ~ 3 in this paper with SDSS spectra yield two similar but different bases of eigenvectors, as shown in the previous section. One would like to know if the larger wavelength coverage of our eigenvectors provides any advantage and how far the new determination is required to reproduce the correct quasar continuum at z ~ 3. In other words, is the prediction of quasar continuum at z ~ 3 better if one uses the new PCA eigenvectors derived in this paper? To answer this question, we applied three tests to the predicted continua.

3.2.1. Error on the predicted PCA continuum in the Lyman-α forest

In order to estimate the difference between the true and predicted continuum in the Lyman-α forest, we derived a set of eigenvectors and a projection matrix using 77 spectra of the training set (out of 78) and estimated the continuum over the Lyman-α forest for the remaining quasar using these parameters and following the method described in Sect. 2.3. This procedure was repeated on each of the 78 spectra in the training set.

Following S05, we estimate for each spectrum the absolute fractional flux error |δF|, defined as (11)where p and q are the predicted and real continua respectively.

This was computed over the restframe wavelength ranges, 1050−1170 Å in the forest, and 1280−2000 Å (or 1280−1600 Å) in the red. The cumulative distributions of the absolute fractional flux errors are plotted in Fig. 9 (in the red) and Fig. 10 (in the forest) using two different wavelength coverages (black line for 2000 Å and grey line for 1600 Å).

In the red, the median error is 5.8% and 5% using the 2000 Å and the 1600 Å decompositions (Fig. 9, black dashed and grey dotted vertical lines respectively) and the 90th percentile error with the 2000 Å decomposition is larger than the error using the 1600 Å set of eigenspectra. This indicates simply that the range up to 1600 Å is easier to fit.

The important result is that the opposite trend is observed in the forest (Fig. 10): median values (around 5%) are very close for the two wavelength coverages (black and grey dashed vertical lines) whereas the 90th percentiles are different (black and grey dot-dot-dash vertical lines) with 9.5% and 12.4% errors for the 2000 and 1600 Å decompositions. This means that using the full coverage reduces the number of outliers with more than 12% error in the Lyman-α forest. Errors in S05 are similar to those we find here but, by using the 2000 Å decomposition, outliers are less frequent than in the previous study. To illustrate what those numbers mean on the continuum level, examples of spectra with their predicted continua are displayed in Fig. 11.

3.2.2. Distribution of spectral indices

We can also compare the characteristics of mock continua generated from the set of computed eigenvectors and the distribution of weights c_ij given in Table 2 with those of real spectra from SDSS-DR7. For this, we fitted a power-law to mock continua and to SDSS-DR7 spectra in the same way (see Sect. 4.2 for more details). The normalized (sum equal 1) distributions of the derived power-law index are displayed in Fig. 12. The distribution from mock spectra is more peaked than the one from SDSS spectra but is centred around the same value. This behaviour is expected because the PCA gives us a mean description of the whole quasar population.

3.2.3. Prediction of the continuum: fitting coefficients versus using the projection matrix

Our goal is to estimate the quasar continuum over the Lyman-α forest. For this, we first estimate the weights of the red part of the spectrum in the basis obtained from PCA of the red part of the spectra. We then multiply these weights by the projection matrix (see Eq. (5)) to compute the weights to be used in the basis obtained from the overall spectrum. This gives us the spectrum reconstructed in the whole wavelength range (method 1). One may be tempted to directly use the coefficients obtained from the red part of the spectrum as representative of the whole spectrum by just replacing the eigenvectors obtained from the red part with those obtained from the full wavelength coverage (method 2). To test how useful the projection is, both methods were used to predict the continuum and the distributions of the absolute fractional flux error (Eq. (11)) were computed and are displayed in Figs. 13 and 14.

When using method 2, not surprinsigly, the red part of the spectrum is very well fitted with a median error of less than 2.5% (Fig. 13, dashed grey line) and method 1 (projection) leads to larger errors (dashed black line in Fig. 13, median error ~6%). In the forest, the trend is opposite with a median error of less than 5% when method 1 is applied (Fig. 14, dashed black line) and more than 7% when method 2 is applied (dashed grey line). In addition, with method 2, 10% of the spectra have an error larger than 15% error (dot-dot-dash grey line) for only one percent with method 1. Clearly the projection matrix should be used.

4. Evolution of the mean flux

An important application of the quasar continuum estimate over the Lyman-α forest wavelength range is the determination of the redshift evolution of the mean flux in the IGM. Numerous authors have performed this measurement using high and/or intermediate spectral resolution data (e.g. Songaila 2004; Bernardi et al. 2003; Dall’Aglio et al. 2008, 2009; Faucher-Giguère et al. 2008). The evolution is smooth except for a possible bump at z ~ 3.2, which could be related to the He ii reionization (Schaye et al. 2000), although this is not the only possible explanation (Faucher-Giguère et al. 2008). In this section we reinvestigate this question applying the method developed in Sect. 3.

4.1. Comparison with B03

We first would like to check if we can recover the Bernardi et al. (2003) results. The sample used in the Bernardi et al. (2003) study is a subsample of SDSS-DR7 that contains all spectra observed up to the end of 2001 (corresponding to a modified Julian day MJD = 52   274). Some selection was applied to remove the most prominent BALs and DLAs. These objects are not clearly defined in Bernardi et al. (2003) so that we had to apply our own selection. We avoided all BALs and DLAs as defined by Noterdaeme et al. (2009). Our final sample has 837 QSOs when B03 had 1041. The comparison to B03 is shown in Fig. 15 (left hand-side panel) for power-law (triangles) and PCA (squares) estimates of the continuum. We also show in the figure the Faucher-Giguère et al. (2008) results. As expected (see next section), the power-law estimate is lower at z < 3 than other estimates. The evolutions found by us and B03 agree and a departure from a smooth evolution is visible at z ~ 3.2.

We randomly drew from SDSS-DR7 a large number (500) of samples identical in size and redshift distribution to the B03 sample. For each sample, we derived the mean flux observed at each redshift and calculated the mean over the 500 samples. The result is shown in Fig. 15 (right hand-side panel). The feature at z ~ 3.2 is still visible and could be a “bump” or a break in the evolution. Note that, as emphasized by Faucher-Giguère et al. (2008), a bump is visible only if one insists on fitting a single power-law.

4.2. Redshift evolution of the mean flux in the IGM using SDSS-DR7

Continua of SDSS-DR7 spectra were fitted using the two different methods we described above. We derived a power-law and a PCA (using our and S05 sets of eigenspectra) continua. We restricted our study to spectra with a signal-to-noise ratio higher than 8 around 1280 Å in the restframe to avoid instabilities in the fit of the power-law. We restricted the analysis to quasars with redshifts higher than z > 2.45 to avoid the blue end of the spectra.

Lines of sight containing damped Lyman-α systems (DLAs) were removed following the lists provided by Noterdaeme et al. (2009) and Prochaska et al. (2005). Broad absorption line quasars (BALs) flagged by Shen et al. (2010) were avoided as well. After this selection, we are left with 2576 quasars. Following Bernardi et al. (2003), we computed the mean flux from the Lyman-α forest between 1080 and 1160 Å in the restframe to avoid the O vi-Lyman-β and Lyman-α emission lines.

The mean flux in the Lyman-α forest was then computed using three different continua: two PCA continua obtained from S05 principal components and the principal components derived in this work, and a power-law continuum. Figure 16 shows the redshift evolution of this quantity in redshift bins of size Δz = 0.1. Error bars were computed from a bootstrap resampling. There is apparently no difference in the results when using the two PCA decompositions derived at low (z ~ 1) and high (z ~ 3) redshifts. The mean flux derived from the power-law continuum is systematically lower at z < 3. This is expected because the power-law tends to overestimate the continuum in the forest (see Sect. 4.3).

Evidently there is a change in the evolution of the mean flux at z ~ 3 with a steepening of the evolution at large redshift.

To estimate what kind of feature we are able to recover with our procedure, we construct in the following DR7-like samples of simulated intermediate-resolution quasar spectra probing an IGM with a mean flux evolution as seen in the high-resolution data and test if we can recover this evolution by applying our procedures.

4.3. Detectability of the bump with DR7

In order to test the detectability of a feature equivalent to the bump seen at z ~ 3.2 by Faucher-Giguère et al. (2008), we performed simulations of this effect and its measurement with mock spectra.

4.3.1. Mock spectra

We aim at simulating the Lyman-α forest that will give a mean flux redshift evolution similar to what is seen in high-resolution data. This evolution was fitted by Faucher-Giguère et al. (2008) as (12)where A = 0.00153, B = 4.060, C = −0.0969, D = 4.267 and E = 0.0769.

We assume that the Lyman-α forest is made up of absorption lines with a column density distribution (13)where α = 4.9  ×  10⁷ and β = 1.46 over the column density range 10¹³−10¹⁷ cm^-2 and a Doppler parameter distribution given by (14)with K = 6.82 and b_σ = 24.09   km   s^-1 (Kim et al. 2001). Therefore the evolution in F(z) is supposed to be caused by an evolution in the number of clouds per unit redshift. This is probably oversimplistic because, if any, the feature at z ~ 3.2 is claimed to be possibly caused by an ionization process but this is probably fine for what we aim to estimate.

We first constructed the relation that gives the mean flux versus the number of clouds present in a redshift bin of Δz = 0.1 at a given redshift. This relation is parametrized by (15)where n_abs is the number of clouds drawn at random from the above population of clouds and a_z and ϵ_z are determined by the simulation and depend on redshift.

Combining Eqs. (12) and (15), we then derived at each redshift the actual number of clouds that is needed to reproduce the relation given by Faucher-Giguère et al. (2008). We finally fitted the redshift evolution of the number of clouds by a similar function: (16)Our best fit gives the values n₀ = 8.281  ±  0.080 and γ = 3.076  ±  0.006, C = −121.3  ±  4.5, D = 4.267  ±  0.003 and E = −0.081  ±  0.003.

Once the number of clouds per unit redshift was correctly calibrated, we generated 50 000 mock spectra with a uniform emission redshift distribution in the range 2.3−4.5. For a given emission redshift, the number of absorption lines was computed from the line number density and was modulated to introduce Poisson noise. For each absorption line, the column density N_HI and the Doppler parameter were randomly chosen following Eqs. (13) and (14). The spectrum was then degraded at the SDSS resolution (R ~ 1800) and a PCA continuum was added using 10 principal components and choosing the weights at random within the calculated distributions. The wavelength scale was binned as for SDSS spectra and noise was added following the SDSS g-magnitude distribution.

To check the validity of our procedures, the mean flux evolution was computed from spectra with no noise and no continuum added (see Fig. 17). The mean flux evolution recovered by our procedure (grey diamonds) excellently agrees with the theoritical input (black dashed line) assumed to follow the evolution as derived by Faucher-Giguère et al. (2008).

4.3.2. Should we detect any feature with DR7?

We computed 100 mock samples with the same number of quasars as the SDSS-DR7 and the same distributions of emission redshift and signal-to-noise ratio. For each sample, we computed the mean flux evolution fitting quasar continuum with a power-law as performed on real data. At each redshift, we computed the scatter in the mean flux derived from these samples. The recovered evolution is shown as black points in Fig. 17, to be compared with the dashed line showing the input assumed for the simulation. As already mentioned, obviously the power-law fit of the QSO continuum underestimates the mean flux. However, the bump at z ~ 3.2, introduced in the input, is recovered although slightly smoothed out by the procedure. The errors derived from the simulations are shown as a grey area. They have to be compared with errors expected from the data. The latter were estimated using the errors obtained in SDSS-DR7 from B03-like samples (as in Fig. 15) but scaled by the square root of the ratio of the number of quasars in the B03 and SDSS-DR7 samples. These errors are shown in Fig. 17 by vertical error bars. They are larger than the errors from the simulations, as expected. Mock spectra are indeed idealized and additional sources of uncertainty are present in the calibration of the data and the consequences on the continuum fit of the somewhat odd shape of some quasars.

4.3.3. Feature at z  ~  3.2

We summarize in Fig. 18 the mean flux evolution measured from the SDSS-DR7 data for a continuum estimated with a PCA (blue squares) or a power-law (red triangles) corrected for the systematic bias as seen in Fig. 17. Indeed, the power-law continuum systematically overestimates the amount of absorption in the Lyman-α forest. The mean flux evolution derived with a power-law continuum is corrected by the expected difference seen in mocks between the measured mean flux and the input of the simulation.

The results of both methods agree excellently. Overplotted as black circles in Fig. 18 are the results by Faucher-Giguère et al. (2008), which agree very well with the PCA estimate except, maybe, in the bin around z ~ 3.2. Note that the discrepancy is less than 2σ. In any case, obviously the smooth redshift evolution of the mean flux becomes steeper around redshift z ~ 3.

The slight difference between high and intermediate resolution data at z ~ 3.2 may be explained by a different selection of the quasars. Indeed, Worseck & Prochaska (2011) have argued that SDSS preferentially selects 3 < z_em < 3.5 quasars with intervening H i Lyman limit systems. However, this should have little influence on the mean flux in the overall forest and cannot explain the discrepancy by itself. This could also be because our procedure smoothes out a sharp feature. However, this would be surprising given the width of the bins and the results of our simulations (see Fig. 17).

5. Conclusion

The first goal of this paper was to provide a new PCA decomposition of the continuum of z ~ 3 quasars. This should be useful for studies of the Lyman-α forest and to generate mock quasar continua that are to be implemented in simulations constructed to search for systematic effects in future analyses or surveys. We took the opportunity to enlarge to 1020−2000 Å the wavelength range over which the spectra are decomposed to compare them with the previous wavelength range 1020−1600 Å used by S05. The mean spectrum at z ~ 3 has a similar shape as the mean spectrum derived at z ~ 1 by S05, except that the strength of the Lyman-α and C iv emission lines relative to the continuum is less at high redshift. We concentrated on the estimate of the continuum in the Lyman-α forest and provided all outputs of this analysis that are required to generate mock continua.

We used this decomposition to revisit the evolution with redshift of the mean flux in the Lyman-α forest and compared two methods to estimate the quasar flux in the Lyman-α forest: the extrapolation of a power-law fitted to the red part of the spectrum and an estimate of the flux from PCA coefficients.

We find that

(i): the power-law method systematically underestimates themean flux by an amount decreasing with redshift. Whencorrecting for this bias, as estimated with simulations, we findthat the method gives similar results to thePCA method;

(ii): the PCA method yields results very similar to what is measured by Bernardi et al. (2003) and from high spectral resolution data (Faucher-Giguère et al. 2008);

(iii): from our simulations a bump at z ~ 3.2, if present, should be marginally detected with the data set of SDSS-DR7;

(iv): finally, from our analysis, we find that there is a definite break in the evolution of the mean flux at z ~ 3 in the sense of a steeper decrease of the mean flux at high redshift. We caution that this could be a consequence of a more prominent bump, which can be slightly smoothed out by our procedures.

The increase of the statistics but most importantly of the quality of the data that will soon be provided by the BOSS survey should definitely settle this point.

Acknowledgments

We thank an anonymous referee for important and very useful comments. This project was supported by the Agence Nationale de la Recherche under contract ANR-08-BLAN-0222.

References

All Tables

All Figures