NEW CONSTRAINTS ON THE QUASAR BROAD EMISSION LINE REGION

Models that predict hydrogen and helium emission line strengths were generated for a range of hydrogen number density (n_H) and hydrogen ionizing flux (Φ_H) values. The plane formed by these parameters represents the range of cloud densities and distances from the ionizing continuum that are expected to exist within the BELR. The output emission of the photoionization simulations is computed as an equivalent width in angstroms relative to the incident continuum at 1216 Å (W₁₂₁₆). The equivalent width is a useful description of the line emission, as it measures how efficiently the line is produced from the incident continuum.

For each point in the grid, Cloudy's standard AGN-ionizing continuum (similar to the Mathews & Ferland 1987 continuum) was used. Despite a range in observed ionizing continua, the H and He emission lines are relatively insensitive to the shape of the ionizing continuum, unlike the metal emission lines (Korista et al. 1997). A constant number density within each cloud was assumed and a constant column density of N_H = 10²³ cm⁻² was also assumed. Although a range of column densities is expected within the BELR, the hydrogen and helium emission lines are not very sensitive to the cloud column density in the range 10²² < N_H(cm⁻²) < 10²⁴ (Korista et al. 1997; A. J. Ruff et al. 2012, in preparation). Note that solar chemical abundances were used for all photoionization simulations.

The range of values for each of the physical parameters in the LOC model is motivated by observations and physical considerations. Broad, semi-forbidden C iii] λ1909 is observed, implying the presence of ∼10⁹ cm⁻³ gas (Davidson & Netzer 1979), as it is collisionally suppressed at higher densities. However, broad forbidden lines are not observed, indicating that the electron density in the BELR is higher than the critical density of the forbidden lines. This gives a lower limit on the number density within the BELR and only clouds with n_H> 10⁷ cm⁻³ were considered. While none of the broad emission lines can be used to determine an upper limit on the number density, analyses of the Balmer lines, whose ratios are sensitive to high-density gas, and strong Fe ii emission indicate that high-density (∼10¹² cm⁻³) gas is also present within the BELR (e.g., Collin-Souffrin et al. 1982; Rees et al. 1989). At these high densities, collisional excitation becomes important. The collisional excitation depends on density, temperature, and ionization state in a fundamentally different manner to the recombination Case B (Ferland et al. 2009). Therefore, detailed photoionization simulations are required.

Reverberation mapping measurements show that broad emission lines originate at a range of radii from the central ionizing source (Peterson 1994). At large distances from the continuum source, the ionizing flux is low and the temperature decreases. Below temperatures of ∼1500 K dust grains condense, making the line emissivity low (e.g., Netzer & Laor 1993). This temperature corresponds to a hydrogen ionizing flux of Φ_H ∼ 10¹⁸ cm⁻² s⁻¹ (Nenkova et al. 2008, after substituting a fiducial bolometric luminosity). The H i and He i emissivities are strongly dependent on the incident ionizing flux, unlike most of the prominent metal lines (Korista & Goad 2004). This is because the optical depth of the H i and He i excited states increases with ionizing flux (Ferland et al. 1979, 1992), until the cloud is fully ionized. Close to the ionizing source, the hydrogen and helium atoms will be fully ionized and will not produce line emission. As the H i and He i line emissivity is low at large Φ_H values, the upper limit on Φ_H is considered to be unimportant.

Considering the limits discussed above, a large range in BELR parameter space was simulated, spanning seven decades in both n_H and Φ_H: 10⁷ ⩽ n_H (cm⁻³) ⩽10¹⁴ and 10¹⁷ ⩽ Φ_H (cm⁻² s⁻¹) ⩽ 10²⁴, each stepped in 0.25 decade intervals. Each model was calculated separately by Cloudy, giving a total of 841 photoionization calculations. Contour plots of each line, showing the W₁₂₁₆ as a function of n_H and Φ_H will be presented in A. J. Ruff et al. (2012, in preparation). A lower limit on the ionization parameter,

$$\begin{equation} U\equiv \frac{\Phi _{{\rm H}} }{c\,n_{{\rm H}} }, \end{equation} \tag{ 1 }$$

was also considered. This is motivated by observations of the calcium infrared triplet, Ca ii xyz λλλ8498, 8542, 8662. If clouds with very low values of the ionization parameter (log U ⩽ −6) existed in the BELR, stronger Ca ii xyz (Joly 1989; Ferland & Persson 1989) and Na i λ5895 emission would be observed (see Korista et al. 1997, for a range of W₁₂₁₆ distributions in the density–flux plane). In addition, the Fe ii emission would be totally dominated by resonance transitions clustered near 2400 Å and 2600 Å (Verner et al. 1999), which is dissimilar to the pseudocontinuum observed. A physical interpretation of the lower limit on the ionization parameter, U_min is that very dense clouds are not found far from the ionizing continuum. Although this limit may not exactly represent the true physical cloud distribution, it is a convenient representation in the LOC model.

These considerations provide only weak constraints on the upper n_H limit, and lower Φ_H and U_min limits in the BELR. However, they are very important because of the large emissivity of H i and He i close to these limits. The H i and He i line ratios are particularly sensitive to the lower Φ_H and upper n_H values (Φ_min and n_max, respectively). Thus, rather than using fiducial values, a range of upper and lower limits on Φ_H and n_H was considered.

Following the standard LOC prescription, the flux from each point in the simulated parameter space was summed to calculate the emission from each species of interest. As shown by Baldwin et al. (1995), the line luminosity is given by

$$\begin{equation} L_{\rm line} \propto \int ^{r_{\rm max}}_{r_{\rm min}}\!\!\int ^{n_{\rm max}}_{n_{\rm min}} W_{1216} (r,n_{{\rm H}})\,f(r)\,g(n_{{\rm H}})\,dn_{{\rm H}} \,dr, \end{equation} \tag{ 2 }$$

where W₁₂₁₆(r, n_H) is the equivalent width of a particular line (relative to the continuum at 1216 Å) from a single cloud at radius r, with number density n_H. The cloud covering fractions as functions of radius and number density are given by f(r) and g(n_H), respectively. For a more complete description of the LOC parameters, see Bottorff et al. (2002). In the simplest case of LOC integration, the covering fractions are given by f(r)∝r⁻¹ and g(n_H)∝n⁻¹_H (Baldwin et al. 1995; see also Matsuoka et al. 2007). The line luminosity can then be simplified to

$$\begin{equation} L_{\rm line} \propto \int ^{\log \Phi _{\rm max}}_{\log \Phi _{\rm min}}\!\!\int ^{\log n_{\rm max}}_{\log n_{\rm min}} W_{1216} (\Phi _{{\rm H}} ,n_{{\rm H}})\, d\log n_{{\rm H}} \,\,d\log \Phi _{{\rm H}} , \end{equation} \tag{ 3 }$$

where the integration limits are free parameters of the model and Φ_H∝r⁻². This is simply a sum over each point in parameter space between the integration limits, since the simulated photoionized clouds in the logarithmic grid are weighted evenly per decade. We note that this weighting scheme maximizes the distribution of cloud characteristics in gas density and incident ionizing flux—a key characteristic of the LOC model. Power-law indices in the density and distance weighting functions (f(r), g(n_H)) that deviate substantially from −1 will take on characteristics of a model with a single gas density and ionizing continuum source distance.⁶

This sum (Equation (3)) then allows us to compute the predicted total equivalent width (W₁₂₁₆) for each emission line, assuming full geometric coverage. The model-integrated covering fraction is then determined from the ratio between the observed W₁₂₁₆(Lyα) and its predicted value, above. While in the present work we are mainly interested in comparing hydrogen Balmer and Paschen emission line flux ratios (which are independent of the integrated covering fraction), Table 5 lists the value of the integrated covering fraction that results in a predicted equivalent width of Lyα of 100 Å, a typical value (see also Section 4.4).

Broad emission line flux ratios were predicted for several emission line pairs using Equation (3). Rather than using fiducial values of the integration limits to predict ratios, a range of upper and lower bounds on each of the free parameters: Φ_min, n_max, and U_min was considered. From the discussion in Section 2.2, the Φ_min, n_max, and U_min integration limits were varied over the ranges

$$\begin{equation} \begin{array}{llc}\log \Phi _{\rm min} &=& [17.00, 20.00] \\ \log n_{\rm max} &=& [10.25, 14.00] \\ \log U_{\rm min} &=& [-6.5, -4.5] \end{array}. \end{equation} \tag{ 4 }$$

To limit the number of assumptions in the calculation, these integration limits span a larger range than suggested by the observational limits discussed in Section 2.2. The Φ_min and n_max integration limits were stepped in 0.25 decade intervals and the U_min limit in 0.5 decade intervals. The Φ_max and n_min limits were also varied, however the effect on the total emission was small, because the H i and He i lines have low emissivity in this region of parameter space. The values of Φ_max and n_min were fixed at 10²³ cm⁻² s⁻¹ and 10⁸ cm⁻³, respectively. The change in predicted ratio caused by varying these limits was less than 1%.

To visualize how the ratios change as the Φ_min, n_max, and U_min integration limits are varied, several flux ratios are plotted against each of the integration limits. Figure 1 plots several predicted ratios: Pa α/Pa β, Pa γ/Pa β, He i 1.083μ/Pa β, Hα/Hβ, and He i 1.083μ/He i λ5876 as a function of Φ_min, n_max, and U_min integration limits. The three panels show ratios plotted as a function of each of the integration limits over the range given in Equation (4), while the other integration limits are held constant at fiducial values (given in the legend).

Zoom In Zoom Out Reset image size

Four quasar spectra from the Glikman et al. (2006) sample have high enough S/Ns to measure the strongest Paschen emission line fluxes. For each emission line in the Sloan Digital Sky Survey (SDSS) and Glikman et al. spectra, the local continuum was measured over a wavelength range ±200 to 500 Å, centered on the line's expected wavelength. Linear fits to the continuum were performed after removing the line itself and any 2σ outliers. To measure line strengths, both Gaussian fits and simple integrations of the flux per wavelength bin were performed for each line. While Gaussians are an imperfect match to the line shapes, they minimize parameterization while allowing for de-blending of the He i and Pa γ lines. Integrating the flux per wavelength bin over ±2 FWHM reproduced the Gaussian line fluxes within the errors. Summing over ±3 FWHM admits smaller outlying features and so tends to result in an overestimate of the line strength. Measured line fluxes are shown in Table 1 together with their 1σ measurement uncertainties. The SDSS line strengths were scaled by f_SDSS to allow direct comparison. An additional 10% systematic error was conservatively added in quadrature to measurement errors in all figures and calculations.

For each of the four objects in the sample, SDSS spectra and Galaxy Evolution Explorer (GALEX) ultraviolet flux measurements were obtained from archival data. Color excess values, E(B − V) (Schlegel et al. 1998) are given in Table 2, and these were used to correct all the data for Galactic extinction following Fitzpatrick (1999). Note that although the continuum fluxes inferred in this section have large uncertainties, the inferred physical conditions do not depend on these values. The measured continuum flux values are only used in Section 4.4 to check that the model produces enough energy to be consistent with previous observations.

For each object, the continuum level of the SDSS spectrum was scaled to match the NIR continuum across the region of overlap. The scaling factor for each SDSS spectrum, denoted f_SDSS, was typically measured to better than 2% accuracy. The continuum flux at 5100 Å (S₅₁₀₀) was measured from the SDSS spectrum, which can be scaled to match the NIR data using f_SDSS. An extrapolation of the best linear fit to the optical continuum was made to give an estimate of the rest-frame 1216 Å continuum flux (S₁₂₁₆(SDSS)). S₁₂₁₆(SDSS) is likely to overestimate the true continuum flux at 1216 Å due to the turnover in the spectrum. An alternative estimate of the 1216 Å continuum flux was obtained from the GALEX far-ultraviolet (FUV) measurement. Due to the shape, breadth, and variable position of Lyα in the GALEX bandpass for our sample, we cannot use the equivalent width of Lyα alone to estimate the 1216 Å continuum flux. Instead we convolve an appropriately redshifted composite quasar spectrum (Vanden Berk et al. 2001) with the GALEX FUV bandpass to estimate the continuum level. The SDSS and GALEX estimates of the 1216 Å continuum flux level are in good agreement, considering the significant uncertainty in extrapolating the continuum. Figure 2 shows the NIR spectra, scaled SDSS spectra, and GALEX photometry.

Zoom In Zoom Out Reset image size

A set of "best-fit" integration limits was found for each object using the measured Pa α/Pa β, Pa γ/Pa β, and Hα/Pa β ratios. Ratios involving optical emission lines were not used, as they suffer differential extinction from dust and the amount of dust is unknown. There are only three independent infrared hydrogen ratios, so the model must be limited to fitting to two parameters. Since the ratios are fairly constant with the lower limit on the ionization parameter, log U_min = −5.5 is fixed, and Φ_min and n_max are free parameters.

The measured emission line ratios of four objects: J005812+160201, J010226−003904, J015950+002340, and J032213+005513 were used to find the most likely set of Φ_min and n_max limits by minimizing χ²:

$$\begin{equation} \chi ^2 = \sum _{i=1}^N \frac{(R_{{\rm obs},i} -R_{{\rm model},i})^2}{ \sigma _{{\rm stat},i}^2 + \sigma _{{\rm sys},i}^2}, \end{equation} \tag{ 5 }$$

where R_obs and R_model are the observed and predicted ratios, and σ_stat and σ_sys are the statistical and systematic uncertainties on each measured ratio. The sum is over each ratio used to find the best model: Pa α/Pa β, Pa γ/Pa β, and Hα/Pa β. R_model is the predicted ratio for each combination of Φ_min and n_max values calculated using Equation (3), over the ranges given in Equation (4). Table 3 shows the minimum χ² value for each object. Figure 3 shows the χ² contours for each object in the Φ_min–n_max plane. This plot shows that the χ² value is minimized where the density upper integration limit is high and the lower limit on the ionizing flux is low.

Zoom In Zoom Out Reset image size

Figure 4 shows the best-fit simulated ratios as a function of the two free parameters: n_max and Φ_min. The best-fit n_max and Φ_min limits are given in the legend for each object. The n_min and Φ_max integration limits were fixed at fiducial values given in Section 2.4. The measured infrared ratios are overplotted as a band, indicating a 1σ uncertainty. Note that the He i 1.083μ/Pa β ratio was not used as a constraint, as this ratio predicts Φ_min and n_max values that are systematically offset from those given by the hydrogen line ratios, as shown in Figure 4.

Zoom In Zoom Out Reset image size

The limits on n_max and Φ_min that were calculated by analyzing the NIR emission lines were then used to predict Hα/Hβ and Hγ/Hβ ratios. The predicted ratios are shown with the measured ratios in Table 3. Two of the four objects are in the Dong et al. (2008) sample, which selected AGNs with blue optical continua in order to minimize the effect of dust extinction. If blue continua do indicate less dust along the line of sight, the simulated ratios should be in good agreement with the measured ratios for these objects. This is true for the two objects listed, however, a much larger sample with a larger number of ratios studied is required to test this hypothesis. J015950+002340 also has a blue continuum, however simulated values cannot be compared to measurements, due to the poor fit to the Hα line.

For all objects, the predicted Balmer decrements are slightly larger than the measured values, which is the opposite effect to dust. However, the predicted values are within measurement uncertainties, which may indicate that there is little internal reddening in these objects. The predicted Hγ/Hβ ratios are consistent with the measured data. The results for each object are summarized briefly in Section 4.5.

Although UV emission lines have not been measured for these objects, a selection of predicted UV emission line ratios is presented in this section. Typical measured strengths of stronger UV lines were compared to the predicted values to check the consistency of our spectral model derived from the measured hydrogen spectrum. We note that the UV metal emission line strengths depend more on the shape of the incident ionizing continuum (and to varying degrees the gas metallicity) than do the Balmer and Paschen hydrogen lines under study (e.g., see Korista et al. 1998 for a discussion). The predicted UV emission line strengths reported here are those that result from the best fit in the physical parameters (n_max, Φ_min) determined from the three hydrogen emission line ratios (as noted in Table 3) in each of the four AGNs.

Predicted Lyα/Hβ ratios for each object are listed in Table 4. Most of the Lyα/Hβ predicted ratios are lower than the simple, classical photoionization prediction (Osterbrock & Ferland 2006). This is because of the inclusion of high-density gas. Lyα thermalizes at high densities, whereas the Balmer lines continue to emit efficiently, causing the predicted Lyα/Hβ ratio to decrease with the inclusion of high n_H gas. Nevertheless, these ratios remain larger than many of those measured in AGNs, and this remains an unsolved problem. However, uncorrected extinction due to grains extrinsic to the Milky Way may explain much of the remaining discrepancies (see Netzer et al. 1995).

A selection of predicted UV metal lines (relative to Lyα) is listed in Table 4. These lines were chosen as they are relatively strong, and they span a large range in ionization parameter. The O vi, C iv, and Mg ii emission lines are also the most important coolants within their ionization zones (excepting emission by Fe ii). As these lines are major coolants, the line strengths are not very sensitive to the gas metallicity, due to the thermostatic effects of such lines. The predicted C iv values are somewhat higher than the ratio measured from the Vanden Berk et al. (2001) composite spectrum, however, the predicted C iv/Lyα ratio is strongly dependent on the hardness of the ionizing continuum. The observed value shows significant variation; for example, the C iv/Lyα ratio is 0.62 in the Zheng et al. (1997) composite, and is closer to 1 in Seyfert 1 spectra (Korista & Goad 2000; see also Osmer et al. 1994). Although additional free parameters, such as metallicity and the shape of the ionizing continuum, have not been considered here, the predicted strengths of the UV metal lines relative to Lyα are generally consistent with typical measured values.

A simple check, albeit with a substantial uncertainty, can establish whether the best-fit models produce enough energy to reproduce the observed total integrated line flux. The largest uncertainty in comparing measured equivalent width values to predicted line strengths is knowing which continuum flux to compare the line flux to. The NIR continuum is a poor indicator of the incident ionizing continuum (i.e., photons with wavelengths shorter than 912 Å), and is contaminated with hot thermal dust emission as well as that of the host galaxy. To test whether the models produce enough energy, the predicted and measured Pa β equivalent widths were compared, relative to two different points in the continuum: 1216 Å and 5100 Å. Note that the predicted equivalent width depends on the integrated covering factor. In this section, the integrated covering factor was calculated such that the simulated equivalent width of Lyα is 100 Å. This is a typical measured value for W₁₂₁₆(Lyα), and was chosen to ensure that the model produces enough energy to be consistent with the observed value. For each object, this value of f_c was used to calculate all equivalent width values presented in this section.

The equivalent width relative to the continuum flux at 1216 Å is often chosen because it is easily observable and close to the ionizing continuum. Contaminants such as light from the host galaxy and dusty torus are negligible in this region. The measured value was estimated using both the GALEX and SDSS UV continuum flux values (S₁₂₁₆) listed in Table 2. The observed rest-frame equivalent width of Pa β relative to the incident continuum at 1216 Å (in Å) is given by

$$\begin{equation} W_{1216}({\rm {\rm Pa}\,\beta }) = \frac{F({\rm {\rm Pa}\,\beta })}{S_{1216}}\times \frac{1}{(1+z)}, \end{equation} \tag{ 6 }$$

where S₁₂₁₆ is the measured incident continuum flux at 1216 Å and F(Pa β) is the total integrated emission line flux, given in Table 1. The dominant measurement uncertainties are quasar variability and scaling the UV measurement to the NIR. There is an additional uncertainty from dust attenuation. The measured S₁₂₁₆ values show a large variance between the GALEX and SDSS estimates, therefore this value is estimated to be accurate within a factor of two.

The measured Pa β equivalent width value relative to the incident continuum at 5100 Å (W₅₁₀₀(Pa β)) was also calculated. At this wavelength the continuum flux is directly measurable from the scaled SDSS spectrum (unlike the SDSS continuum at 1216 Å, which is an extrapolated estimate). However, the 5100 Å continuum is much further from the ionizing continuum and host galaxy contamination is not negligible at this wavelength, making it a poorer proxy for the ionizing continuum. There may also be some thermal emission from the BELR clouds and reflected light from the dusty torus contributing to the measured continuum flux at this wavelength. The calculation of W₅₁₀₀(Pa β) is analogous to the calculation of W₁₂₁₆(Pa β) given in Equation (6). Given the large uncertainties in scaling the continuum measurements to the NIR, quasar variability, and host galaxy contamination in the S₅₁₀₀ measurement, the measured values listed in Table 5 should only be considered accurate to within a factor of two.

Using the best-fit n_max and Φ_min integration limits determined in Section 4.1, the simulated W₁₂₁₆(Pa β) and W₅₁₀₀(Pa β) values were calculated for each object. Both measures of equivalent width can be output by Cloudy directly. (Note that these estimates assume the standard Cloudy continuum described in Section 2.1.) To compare the simulated equivalent width values to the measured values, predicted W₁₂₁₆(Pa β) and W₅₁₀₀(Pa β) values must be scaled by the integrated cloud covering fraction (f_c). For each object, the integrated covering fraction required to produce W₁₂₁₆(Lyα) = 100 Å (a typical value) was computed for the best-fit parameters. This value of f_c was then used to scale the predicted equivalent width values.

The measured and predicted W₁₂₁₆(Pa β) and W₅₁₀₀(Pa β) values, and the value of f_c are listed in Table 5. As shown in this table, the computed integrated covering fraction required to produce W₁₂₁₆(Lyα) = 100 Å for each object is not unrealistically large. The measured and predicted W₁₂₁₆(Pa β) values are also in good agreement (with the exception of J005812+160201, which has a shallower UV continuum slope than the Cloudy AGN continuum). Therefore, Table 5 shows that the best-fit models produce enough energy to be consistent with both the measured Pa β and previously observed Lyα equivalent width values.

A brief description of the results for each object is given below.

J005812+160201. The model is a good fit to the data, although the predicted Hγ/Hβ ratio is 2σ larger than the measured ratio. The measured Hγ/Hβ ratio is low compared with the other objects. The predicted total equivalent width barely produces enough energy to be consistent with the observations. The low predicted equivalent width of Pa β is at least partially due the assumed Cloudy continuum, which is a poor approximation in this case, as the measured UV spectrum is much shallower than the Cloudy continuum. The measured W₁₂₁₆(Pa β) is high compared to the other objects, again because the continuum slope is very shallow.

J010226−003904. The model is a good fit to the data. This object has a strange broad absorption line (BAL)-like feature in He i 1.083μ, probably due to the atmosphere. However, this feature could be real, similar to that reported by Leighly et al. (2011). The predicted Φ_min value is quite high, which results from the large measured Pa γ/Pa β ratio.

J015950+002340. The model is a fair fit to the data, as the reduced χ² value is 3.1. This source, also known as Mrk 1014, is an ultraluminous infrared galaxy (ULIRG; see Armus et al. 2004, and references therein). The Hα emission line is highly asymmetric and is not well fit by a Gaussian, making comparisons to simulations difficult. The predicted high density could be a consequence of the underestimation of the Hα flux. Note that the measured Hα/Pa β and Hα/Hβ ratios are low compared with the other objects. Also note that the optical and NIR He i line flux (relative to the hydrogen line flux) is high compared with the other objects.

J032213+005513. The model produced a good fit to the emission line ratios, although the predicted Hα/Hβ ratio is 2σ above the measured value.

In three of the four cases, a good fit to the hydrogen emission lines is obtained by the LOC model. In the fourth case the fit is fair. In addition, the shape of the χ² contour distributions in the density–flux plane is similar for all objects. High Φ_min limits are generally highly disfavored. This suggests that log Φ_min < 18, and that low-incident-ionizing fluxes are required to reproduce observed hydrogen emission line flux ratios. The determined limits on Φ_min are also consistent with the dust sublimation radius (Nenkova et al. 2008). This is in agreement with several recent studies that suggest the Balmer lines and Mg ii in large part come from the outer BELR and near or at the dust sublimation radius (Zhu et al. 2009; Czerny & Hryniewicz 2011; Mor & Netzer 2012).

Figure 3 shows that lower n_max values in combination with low Φ_min values are also highly disfavored (i.e., χ² values are high in the lower left corner of the plot), meaning that BELRs with a low high-density cutoff and low low-incident-ionizing flux cutoff are unlikely solutions. The minimum of the χ² contours tend to lie on a diagonal in the n_max–Φ_min plane: if n_max is larger, then Φ_min is smaller. If the BELR includes regions where the incident ionizing flux is very low, then gas with a very high number density is required to reproduce the observed emission line ratios. This solution is quite unlikely though, as dust grains condense at low Φ_H, destroying the line emission. The validity of the inferred limits is uncertain. A larger sample, with higher signal to noise is required to test these results.

The predicted Hα/Hβ and Hγ/Hβ ratios are generally in very good agreement with the measured ratios, as shown in Table 3. However, the measured He i 1.083μ/Pa β ratio does not match the simple LOC model described here. We also find that the predicted He i 5876/Hβ line ratio is significantly higher than the measured values. The simulated values are typically a factor of ∼2 greater than the measured values. This offset could be caused by the simplicity of the covering fractions assumed in this model. The shape of the W₁₂₁₆(He i 1.083μ) contours is quite different to the Balmer and Paschen lines, making ratios involving He i 1.083μ more sensitive to the power laws assumed in Equation (2). Although the slope of the power laws is not expected to deviate significantly from −1, refinement of the models using different power-law distributions in f(r) and g(n_H) will be the subject of future studies. The predicted overestimate of both the He i λ5876/Hβ and 1.083μ/Pa β line ratios could also indicate that Cloudy overpredicts the He i emission spectrum relative to H.

The amount of attenuation from internal dust (i.e., dust within the host galaxy) is unknown. The Balmer decrement is the most widely measured ratio is AGN spectroscopy (e.g., Osterbrock 1977; Groves et al. 2012), and has long been used as an indicator of the amount of extinction (e.g., Berman 1936). Dong et al. (2008) measured the Balmer decrement for a sample of 446 Seyfert 1s and QSOs with blue continua, which were selected to minimize the effect of dust extinction. They found only a small intrinsic scatter in the measured ratio, which suggests that the intrinsic Balmer decrement could be relatively consistent between objects, and the observed variance is due to different amounts of dust along the line of sight within the host galaxy. However, this result relies on the slope of the continuum being a dust indicator (see Davis et al. 2007 for a discussion).

The continua of our subset of quasars are reasonably blue, with the exception of J005812+160201. If dust attenuation caused reddening of the continuum in this case, a larger observed Hα/Hβ would be expected. Unfortunately, it is not possible to test this assertion with such a limited sample size. The presence of uncorrected extinction due to dust along the line of sight would also increase the observed He i 1.083/He i 5876 ratio above that of the predicted one (see Table 3).

With better quality spectra, an estimate of the amount of dust attenuation could be made. The Lyα/Hβ problem could be at least partially explained by internal extinction (Shuder & MacAlpine 1979) and constraining the amount of dust attenuation is important to resolve this matter.

The value of the Pa α/Pa β line ratio measured from the Glikman et al. (2006) composite spectrum was 0.64 ± 0.01. However, this value is much lower than the values measured from individual objects, which are consistently above 1.0 (see Table 3). Other observed Pa α/Pa β emission line ratios in quasars are consistent with this result (e.g., Soifer et al. 2004, who measured a Pa α/Pa β ratio of 1.05 ± 0.2 in a z = 3.91 quasar). Ratios of Pa α/Pa β that are above 1.0 have also been measured in lower-luminosity AGNs (Landt et al. 2008).

While it is possible that systematic effects cause objects with larger Pa α/Pa β ratios to be preferentially selected in our sample, a bias caused by the composite methodology is a more likely explanation. Composite spectra are useful for describing generic properties of a sample and for identifying emission lines that are too weak to be observed in individual spectra. However, care must be taken in ascribing physical significance to individual features of composite spectra. It is not clear whether composite spectra are useful for investigating detailed physical properties of the line-emitting gas. Emission line ratios in the composite spectrum can also be affected by poor S/Ns of individual spectra, and wavelength ranges where there are few objects contributing to the final spectrum. See Sulentic et al. (2002) for a discussion of how emission line ratios are affected, and an alternative method for producing composite spectra.

The NIR spectroscopy was not simultaneous with the optical and UV measurements. As quasars are intrinsically variable, this will introduce a systematic uncertainty. Even after continuum matching, the line fluxes will be altered; this effect is known as the intrinsic Baldwin effect (e.g., Pogge & Peterson 1992). There is often a considerable difference between GALEX and SDSS S₁₂₁₆ values. The objects in Table 2 SDSS are all observed during 2000. NIR observations are from 2003 September to 2004 September and GALEX observations are in 2003 (and one in 2005).

Fitting the emission lines with Gaussians introduces another systematic uncertainty. Glikman et al. (2006) point out that Gaussians do not fit the broad wings well; however, Gaussians generally fit the bulk of the broad emission line flux, while requiring only a small number of free parameters. The estimated contribution from the narrow component is less than 5%. Dong et al. (2008) deconvolved the emission lines into both a broad and narrow component, and in each case the equivalent width of the narrow component is only a few percent of the broad component. Fe ii contamination was not considered. Fitting the Fe ii contamination would introduce additional errors, given the low signal to noise of the data.

Fixing the lower limit on the ionization parameter, U_min is not ideal, however, the ratios considered are insensitive to U_min and with only a few measured ratios, the model must be limited in complexity. The U_min limit would be better determined with the Ca ii infrared triplet (Ca ii xyz) or other U_min-sensitive ratios (see A. J. Ruff et al. 2012, in preparation for a discussion). Single point solutions in the density–flux plane were not considered for reasons discussed in Section 2.2.

Since the Balmer and Paschen lines are optically thick, the geometry will affect the predicted emission. If the BELR geometry is such that emission from the "back" of the clouds is not observed, the predicted emission line ratios will change (Ferland et al. 1992; O'Brien et al. 1994; A. J. Ruff et al. 2012, in preparation).

Any broad emission line flux measurement is from a snapshot in time. The ratios will vary with time as the continuum luminosity varies, and this will be dependent on the distribution of the gas—both in geometry and physical state. Comparing measured flux ratios can be problematic if the source is highly variable and lines are formed at very different radii. The recombination lines are emitted relatively efficiently over a broad range of physical conditions, and their continuum reprocessing efficiencies are largely dependent on the incident flux, and therefore the distance from the ionizing source, due to optical depths in the excited states (Korista et al. 1997; A. J. Ruff et al. 2012, in preparation). This manifests in different measured reverberation mapping signals and line profile variability (e.g., Sergeev et al. 2001; Kollatschny & Zetzl 2010). For example, not only do the Balmer lines show different measured time lags, but the line responsivity is observed to decrease from line center to line wing (see Korista & Goad 2004, for a discussion). This could be corrected for if the distribution of broad emission line gas could be determined using high-quality spectral and temporal data in conjunction with photoionization models (Horne et al. 2003).