Table Of ContentApJ,in press
PreprinttypesetusingLATEXstyleemulateapjv.6/22/04
CONSTRAINTS ON REIONIZATION AND SOURCE PROPERTIESFROM THE ABSORPTION SPECTRA OF
Z >6.2 QUASARS
Andrei Mesinger1,2 & Zolta´n Haiman2
ApJ, in press
ABSTRACT
We make use of hydrodynamical simulations of the intergalactic medium (IGM) to create model
quasar absorptionspectra. We compare these model spectra with the observedKeck spectra of three
z > 6.2 quasars with full Gunn-Peterson troughs: SDSS J1148+5251 (z = 6.42), SDSS J1030+0524
7 (z = 6.28), and SDSS J1623+3112 (z = 6.22). We fit the probability density distributions (PDFs)
0 of the observed Lyman α optical depths (τα) with those generated from the simulation, by adopting
0 a template for the quasar’s intrinsic spectral shape, and exploring a range of values for the size of
2 the quasar’ssurrounding HII region, R , the volume-weightedmean neutral hydrogenfraction in the
S
n ambient IGM, x¯HI, and the quasar’sionizing photon emissivity, N˙Q. In order to avoidaveragingover
a possibly large sightline–to–sightline fluctuations in IGM properties, we analyze each observed quasar
J independently. We find the following results for J1148+5251, J1030+0524, and J1623+3112: The
0 best–fit sizes RS are 40, 41, and 29 (comoving) Mpc, respectively. For the later two quasars, the
3 valueissignificantlylargerthantheradiuscorrespondingtothewavelengthatwhichthequasar’sflux
vanishes. These constraints are tight, with only 10% uncertainties, comparable to those caused by
2 redshift–determinationerrors. Thebest–fitvalue∼sofN˙ are2.1,1.3,and0.9 1057 s−1,respectively,
Q
v with a factor of 2 uncertainty in each case. These values are a factor of ×5 lower than expected
8 fromthetemplat∼espectrumofTelferetal. (2002). Finally,thebest–fitvalues∼ofx¯ are0.16,1.0,and
HI
5
1.0,respectively. The uncertainty in the case of J1148+5251is large,and x¯ is not well constrained.
2 HI
However,for both J1030+0524and J1623+3112 we find a significantlower limit of x¯ >0.033. Our
0 HI
method is different from previous analyses of the GP absorption spectra of these quas∼ars, and our
1
results strengthen the evidence that the rapid end–stage of reionization is occurring near z 6.
6
∼
0 Subject headings: cosmology: theory–earlyUniverse–galaxies: formation–high-redshift–evolution
/ – quasars: spectrum
h
p
-
o 1. INTRODUCTION ization. TheLyαopticaldepthsassociatedwithsuchGP
r troughs, and more importantly their accelerated evolu-
The epoch of reionization, when the radiation from
st early generations of astrophysicalobjects ionized the in- tionatz >6,suggeststhatthereionizationepochisend-
a tergalacticmedium(IGM),offersawealthofinformation ing at z ∼6, with a lower limit on the volume-weighted
v: aboutcosmologicalstructureformationandphysicalpro- IGM neu∼tral fraction, x¯HI(z 6)>10−3, along at least
∼
i cessesintheearlyuniverse. Onlyrecentlyhavewebegun thelineofsight(LOS)totwoquasa∼rs(White et al.2003;
X
to gather clues concerning this epoch. The Thomson Fan et al. 2006a). The sharp decline in flux at the edge
r scattering optical depth, τ = 0.09 0.03 (Page et al. of the GP trough of quasar J1030+0524 can be used to
a e
2006; Spergel et al. 2006b), measured±from the polariza- infer the stronger constraint of x¯HI(z 6) > 0.2, along
∼
tion anisotropies in the cosmic microwave background thatLOS(Mesinger & Haiman2004). Theo∼bservedsize
(CMB) by the Wilkenson Microwave Anisotropy Probe of the quasars’ surrounding HII regions can be used to
(WMAP), suggests that cosmic reionization began at infer a similar, independent lower limit (Wyithe & Loeb
redshift z > 10. However, the current uncertainty of 2004; Mesinger & Haiman 2004), although this requires
this value, a∼nd the fact that the optical depth measures additionalassumptionsconcerningthequasaranditsen-
only an integrated electron column density, means that vironment, and the inferred neutral fraction scales di-
many possible competing reionization scenarios can still rectly with these assumptions (Bolton & Haehnelt 2006;
be accommodated. Maselli et al. 2006). The rapid evolution in the sizes of
In order to chart the progress of reionization, several the HII regions for different quasars between z 5.7
∼
attempts have been made to measure the neutral frac- and z 6.4 suggests that the average neutral frac-
∼
tion of the IGM at high redshifts (z > 6). The detec- tion increased by a factor of 10 during this interval
∼
tion of dark spectral regions, so-calle∼d Gunn-Peterson (Fan et al.2006a),althoughtherearealternativewaysto
(GP) troughs, in the Lyman line absorption spectra of accountfor the steep size–evolution(Bolton & Haehnelt
high-redshiftquasarsdiscoveredinthe SloanDigitalSky 2006). The distribution of spectral sizes of dark gaps
Survey (SDSS) offers a probe of the late stages of reion- should help distinguish these scenarios in the future
(Gallerani et al. 2006; Fan et al. 2006a). Upper limits
1 YaleCenterforAstronomyandAstrophysics,YaleUniversity, ontheneutralfractionatz 6havealsobeensuggested
260WhitneyAvenue, NewHaven, CT06520 by other means. The lack o∼f significant evolution in the
2 Department of Astronomy, Columbia University, 550 West
luminosity function (LF) of Lyman α emitters between
120thStreet,NewYork,NY10027
z 5.7 and z 6.5 was used to infer x¯ (z 6) < 0.3
HI
≈ ≈ ∼
∼
2
(Malhotra & Rhoads 2004; Haiman & Cen 2005), if the the sightline-to-sightline fluctuations in these quantities
galaxies are assumed to lie in isolation, or x¯ (z 6)< (Mesinger & Haiman 2004).
HI
∼
0.5(Furlanetto et al.2006),ifananalyticprescriptiono∼f The purpose of this paper is to derive these param-
the clustering ofHII bubbles (Furlanetto et al.2004b)is eters, i.e. the hydrogen neutral fraction, the size of
taken into account. However, a recent, more accurate surrounding HII region, and the ionizing emissivity, by
measurementofthe LFatz 6.5revealedevidence ofa separately modeling the spectra of each of the z > 6.2
decrease,weakeningtheseup≈perlimits(Kashikawa et al. quasars. Specifically,weanalyzeJ1148+5251(z =6.42),
2006). Finally, the lack of noticeable absorption due J1030+0524 (z = 6.28) and J1623+3112 (z = 6.22). At
to the Lyα damping wing in the afterglow spectrum of present, a 4th quasar is known with a full GP trough.
the gamma-ray burst 050904 at z = 6.3 can constrain This source, however, is a broad absorption line (BAL)
x¯ (z 6)<0.2 along that LOS (Totani et al. 2006). quasarJ1048+4637(z =6.2),withcomplexspectralfea-
HI
The∼prep∼onderance of constraints mentioned above tures (Fan, X., private communication) that preclude a
suggeststhatthe combinationofallexistingdatais con- simple interpretation in terms of IGM properties. We
sistent with x¯ (z 6) 0.1. While this is plausible, have therefore omitted this quasar from our study.
HI
the constraints liste∼d abo∼ve allrely onvarious model as- Our analysis differs from our previous work
sumptions. For example, the conversion of an observed (Mesinger & Haiman 2004) where we modeled the
average flux decrement to a neutral fraction (Fan et al. gross transmission features in Lyα and Lyβ close to
2002; White et al. 2003; Fan et al. 2006a) is difficult, the edge of the HII region surrounding J1030+0524. In
and requires precise ab–initio knowledge of the IGM the present study, we perform a much more detailed
density distribution at high-redshifts (Oh & Furlanetto analysis, by fitting the quasar’s intrinsic emission, and
2005; Fan et al. 2006a). All of the other constraints modeling the distribution of optical depths, pixel–by–
mentioned above rely on the IGM density distribution, pixel, blueward of the Lyα line center, following the
as well, albeit less sensitively. Techniques involving the procedure proposed in Mesinger et al. (2004). Unlike in
sizes of HII regions (Wyithe & Loeb 2004; Maselli et al. Mesinger & Haiman (2004), here we do not make use
2006; Fan et al. 2006a) are subject to degeneracies be- of the Lyβ absorption spectrum, except to set a lower
tweenquasarlifetime,ionizingflux,andtheIGMneutral limit in the allowed size of the HII region surrounding
fraction. Further complicating matters for the HII bub- J1030+0524, the only one among the three quasars
ble size techniques is the fact that the onset of the GP that has a noticeable offset between the Lyα and Lyβ
troughneednotcorrespondtotheedgeoftheHIIregion. GP troughs. Our reason for ignoring the Lyβ region is
The GP through just corresponds to the region where that the majority of the pixels we analyze (see Figure 2
the flux drops below the detection threshold, and thus below) have flux detected in the Lyα absorption region,
merely provides a lower limit on the size of the HII re- which can directly be converted to a Lyα optical depth.
gion (Mesinger & Haiman 2004; Maselli et al. 2006; see This avoids the need for a statistical modeling of the
also Bolton & Haehnelt 2006 for a recent, comprehen- foregroundLyα absorption, necessary to convert Lyβ to
sive review of these issues). Furthermore, determining a Lyα optical depth.
thenumberdensityandevolutionofLyαemittinggalax- This paper is organized as follows. In 2, we describe
§
ies is complicated and depends on knowledge of galaxy our analysis technique, including generating the simu-
clustering(Haiman & Cen2005;Furlanetto et al.2004a, lated absorption ( 2.1), fitting for the quasars’intrinsic
§
2006),absorberclustering,aswellasthelarge-scaleden- emission spectrum ( 2.2), and the statistical method
§
sity (e.g. Barkana & Loeb 2004b,a) and velocity fields used to compare the observed and simulated spectra
(e.g. Cen et al. 2005), which are all unconstrained at ( 2.3). In 3, we present our results for each of the
§ §
the relevant high redshifts. quasars we analyze. In 4, we discuss various aspects
§
More generally, the drawback of many of these meth- of our results, such as the origin of the constraints we
ods is that they rely on the combination of data from obtain, and the robustness of our results in the face of
multiple LOSs to obtain statistical significance. How- uncertaintiesinthequasartemplatespectrum. In 5,we
§
ever, assuming constancy in the IGM properties along summarizeour key findings andpresentour conclusions.
disparateLOSsis dangerousathigh-redshifts,especially
approachingtheendofreionization,whenlargesightline- 2. ANALYSIS
to-sightline fluctuations are expected in the density field In this section, we describe the three components of
andthe ionizingbackground(see,e.g., recentreviewsby our analysis: (i) a hydrodynamical simulation to model
Djorgovskiet al. 2006and Lidz et al. 2006). Indeed, the theabsorption( 2.1);(ii)modelingthequasars’intrinsic
current sample of high-redshift quasars already reveals emission spectru§m ( 2.2); and (iii) the statistical com-
large fluctuations in the IGM properties along different parisonofobserveda§ndsimulatedspectra( 2.3). Unless
LOSs (Fan et al. 2006a). This would suggest that the stated otherwise, we adopt the backgroun§d cosmologi-
most rewarding technique would be able to analyze and cal parameters (Ω , Ω , Ω , n, σ , H ) = (0.76, 0.24,
Λ M b 8 0
obtain constraints from each source independently from 0.0407, 1, 0.76, 72 km s−1 Mpc−1), consistent with the
other sources. This is statistically less powerful (e.g. three–yearresults by the WMAP satellite (Spergel et al.
Mesinger et al. 2004), but the rewards can be rich, po- 2006a). Weuseredshift,z,andtheobservedwavelength,
tentially providing additional constraints on the proper- λ , interchangeably throughout this paper as a mea-
obs
tiesofindividualquasarsandtheirenvironments,suchas sureofthedistancealongthe lineofsight. Unlessstated
thequasar’semissivityofionizingphotons,sizeofitssur- otherwise, we quote all quantities in comoving units.
rounding HII region, sharpness of the HII regionbound-
ary, the quasars’X-ray emission, as well as a measure of 2.1. Simulating Absorption in the IGM
3
In order to simulate the absorption spectra, we make Since high-redshift quasars are surrounded by their
useofahydrodynamicalsimulationboxinaΛCDMuni- own highly ionized HII regions, it is useful to express
verse at redshift z = 6. The details of the simulation τ as a sum of contributionsfrom inside (τ ) andoutside
R
are described in Cen et al. (2003), and we only briefly (τ ) the HII region, τ = τ +τ . Defined as such, the
D R D
discuss the relevant parameters here. The box is 11 h−1 resonant optical depth, τ , would be given by Eq. (1),
R
Mpc on each side, with each pixel being about 25.5 h−1 with the lower limit of integration changed to the red-
kpc. This scale resolves the Jeans length in the smooth, shiftcorrespondingto the edge ofthe HII region,z(R ).
S
partially-ionized3 IGM by more then a factor of 10. The The damping wing optical depth τ , would be given by
D
box also has a slightly different set of cosmological pa- Eq. (1), with the upper limit of integration changed to
rameters: (Ω ,Ω ,Ω ,n,σ ,H )=(0.71,0.29,0.047,1, z(R ) 5.
Λ M b 8 0 S
0.854 , 70 km s−1 Mpc−1). Since cosmic hydrogen is highly ionized at low red-
Density and velocity information was extracted from shifts, one can replace the lower limit of integration in
the simulation box along randomly chosen lines of sight. Eq. (1) with z , denoting the redshift below which HI
end
The step size along each LOS was taken to be 5.1 absorption becomes insignificant along the LOS to the
h−1 kpc, which resolves the Lyman α Doppler width source. In our analysis, we chose the value of z sep-
end
in the partially-ionized IGM by more than a factor of arately for each quasar to correspond approximately to
40. The exact value of the step size was chosen some- the blue edge of its GP trough. This roughly translates
what arbitrarily, and does not influence the results as toz 0.2–0.3forourthreequasars. We findthatthe
end
∼
long as it adequately resolves the Doppler width. At exactchoiceofz isnotimportant,becausemostofthe
end
each step, the density and velocity values were aver- contributiontoτ comesfromneutralhydrogencloseto
D
aged for the neighboring pixels, and weighted by the z(R ). Inside the HII region, we assume an IGM tem-
S
distance to the center of the pixels. We extended each perature of T = 2 104 K. Outside the HII region, for
×
LOS by the common practice of randomly choosing a thecaseofamostlyionizeduniverse,weassumeanIGM
LOS through the box, and stacking several segments temperature of T = 104 K, while the temperature in a
together (Cen 1994; Mesinger et al. 2004; Kohler et al. neutraluniversewastakentobeT =2.73 151(1+z)2K,
× 151
2005; Bolton & Haehnelt 2006). The r.m.s. mass fluc- valid for z < 150 (Peebles 1993). Our results are fairly
tuation on the scale the box is σ(z = 6.3) 0.15. insensitive to the exact temperature used. The Doppler
∼
Thisisasmallenoughvaluethatneglectingwavelengths width of the Lyman-α absorption cross section scales as
longer than the box size is justified in our analysis. ν T1/2, but the totalintegratedarea under the cross
D
Note however, that this is not the case for modeling sect∝ion is independent of temperature.
highly non-linear processes, such as reionization (e.g. Assumingionizationequilibrium,wecomputetheneu-
Barkana & Loeb2004b). WhenusingsuchaLOSinsim- tralhydrogenfractionateachpointalongthe LOSwith:
ulatingagivenquasarspectrum,densityvaluesfromthe
simulation were scaled as (1+z)3. nx (Γ +Γ )=n2(1 x )2α , (2)
HI Q BG HI B
The total optical depth due to Lyman α absorption, −
τ, between an observerat z = 0 and a source at z = zQ, where ΓQ and ΓBG are the number of ionizations per
at an observedwavelengthof λ =λ (1+z ), is given second per hydrogen atom due to the quasar’s and the
obs 0 Q
by: background’s ionizing fluxes, respectively, and αB is the
case-B hydrogen recombination coefficient. The LHS of
Eq. (2) accounts for the number of hydrogenionizations
zQ dt λ per cm3, while the RHS accounts for the number of hy-
obs
τ(λobs)=Z0 dzcdz n(z)xHI(z)σ"(1+z)(1− v(cz))# donroegoebntareincos:mbinations per cm3. Solving Eq. (2) for xHI
(1)
where n(z), v(z), xHI(z) are the total hydrogen number b √b2 4
density,thecomponentofthepeculiarvelocityalongthe xHI = − − − , (3)
2
LOS, and the local hydrogen neutral fraction, respec-
tively, at redshift z. The Lyα absorption cross section, where b 2 (Γ + Γ )/(nα ). We further
Q BG B
≡ − −
to first order in v(z)/c, is given by σ[λ /(1+z)/(1 parametrize the quasar’s ionization rate with:
obs
−
v(z)/c)].
Γ f Γe (4)
3ThetemperatureofaneutralIGMbeforeanyionizationwould Q ≡ Γ Q
beseveral orders of magnitude lower than inthe partially-ionized where Γe = (4πr2)−1 ∞1.55 1031 (ν/ν )−1.8 [(1+
case,andtheJeanslengthintheIGMwouldbeunder-resolvedby Q νH × H
af4acTtohreofth∼r6e.e–year WMAP results yielded a smaller value of zb)e/in(1g+thzeQi)o]−n0iz.8a(tσio/nhνfr)edRqνuse−n1c,ywoifthhyνdHro=ge3n.2a9n×d1r0b15eiHngz
σ8=0.76±0.05(Spergeletal.2006b). Asmallerσ8 wouldresult the luminosity distance between the source and redshift
in a narrower distribution of optical depths than those obtained z. Γe results from redshifting a power–law spectrum
from the box we use. Although this does not exactly mimic the Q
effect of increased damping wing absorption from a more neutral with a slope of ν−1.8, normalizedsuch that the emission
universe,whichshiftstheopticaldepthdistributiontohigherval-
ues, it is likely that if confirmed, a smaller σ8 would weaken our 5 Notethattheterms“resonant”and“dampingwing”areonly
lowerlimitsontheIGMneutral fractionbelow. Estimatesforthe
appropriate in describing the contributions of τR and τD inside
sisenqsuitainvtitityatoifveoluyrdreiffisuclutsltowniothuoructhoruicnenoinfgcoasmsuoiltoegiocfalmpoadrealms.etWeres the HII region (i.e. at λobs >∼ λα[1+z(RS)]). Outside of the
note however that analysis combining the three–year WMAP re- HIIregion(λobs<∼λα[1+z(RS)]),τD actuallyintegratesoverthe
sultswithdatafromtheLyman–αforestsuggestsahighervalueof resonantpartoftheabsorptioncrosssectionandτRintegratesover
σ8=0.86±0.03(Lewis2006),whichisconsistentwithourchoice. thedampingwing(seeFig. 1).
4
rate of ionizing photons per second is 1.3 1057, match- 1. R , the radius of the quasar’s HII region,
S
×
ing the emissivity one obtains (Haiman & Cen 2002) for
SDSS J1030+0524 using a standard template spectrum 2. xIGM, the neutralhydrogenfractionin the IGM at
HI
(Elvis et al. 1994). Note that the normalization inferred mean density,
from the more recent template of Telfer et al. (2002) is
3. f , the quasar’s ionizing luminosity, in units of
a factor of of 5 higher. This roughly translates to an Γ
∼ 1.3 1057 s−1.
emission rate of ionizing photons of
×
To obtainanintuition forthe effect of varyingthese free
N˙ =f 1.3 1057 s−1 . (5)
Q Γ parameters, in Figure 1 we plot τ (solid curve) and
× × R
τ (dashed curve) for a typical simulated LOS, for pa-
We ignore radiative transfer along the LOS (and also D
theionizationofhelium). Instead,outsideofthequasar’s rameter choices (RS, xIHGIM, fΓ) = (40 Mpc, 1, 1). As
HII region, we compute x assuming the quasar con- stated above, the total Lyman α optical depth is the
HI
tributes nothing to the ionizing flux, i.e. ΓQ(R>RS)= sum of these two contributions, τ = τR + τD. Roughly
0,while inside the HII region,weassume the gasis opti- speaking,changing RS movesthe dashed (τD) curve left
callythininthe ionizingcontinuum. Thisassumptionof and right,while changingxIGM moves this curve up and
HI
a sharp HII region boundary is an approximation, valid down in Figure 1. Changing fΓ approximately corre-
in the regime where either the quasar’s spectrum is soft, spondstomovingthesolid(τR)curveupanddown. The
or where a strong soft ionizing background flux domi- dotted line demarcates roughly the maximum total op-
nates over that of the quasar already at radii smaller tical depth, τ 5.5, at which a signal can be detected
∼
that the mean free path of the typical ionizing quasar withaconfidencegreaterthan1-σ,inthespectralregion
photon (as is the case of the proximity effect at lower just bluewardof the Lyα emissionline (see discussion in
redshift; Bajtlik et al. 1988) If the quasar has a hard 2.3). Finally, after the optical depth and the associ-
spectrum, and the background intensity is low, then the §ated transmission (e−τ) is calculated according to the
edge of the HII region will be “blurred”. A similar blur- procedure outlined above, the transmission is averaged
ring would arise if additional ionizing sources inside the over 3.3 ˚A wavelength bins to match the Keck ESI
≈
HII region, with an extended spatial distribution, would spectra with which the simulated spectra are compared
contributesignificantflux (Wyithe & Loeb2006),a pos- (note that the actual smoothing length ∆λobs differs by
sibility we do not address here. However,already strong a few percent for each source, since they are at different
upper limits can be placed on the extent of such blur- redshifts).
ring from the offset between the wavelengths of the Lyα ThefreeparameterRS wasvariedinlinearincrements
and Lyβ GP troughs. Since Lyβ absorption is less sen- of∆RS =2Mpc,withthelowerlimitsetbytherededge
sitive than Lyα, the amount of offset between the GP of the Lyα GP through in the observed spectra. The
troughs gives a measure of the spatial (i.e. wavelength) otherparameterswerevariedintheranges: 3.2 10−4
evolution of τ close to the edge of the HII region, with xIHGIM ≤1, in multiples of 5, and 0.1≤fΓ ≤10,×in linea≤r
a large offset indicating slow τ evolution, and a small increments of ∆fΓ = 0.3 between 0.1 fΓ 2 and
≤ ≤
offset indication rapid τ evolution. Mesinger & Haiman ∆fΓ =2 between 2<fΓ 10.
≤
(2004) studied this forthe case ofJ1030+0524,but even
stronger constraints could be placed on the spectra of 2.2. Modeling the Intrinsic Emission Spectrum
J1148+5251 and J1623+3112, since they do not show As observational input in our analysis, we make
any offsets between the GP troughs (Kramer et al., in use of Keck ESI spectra of the three highest redshift
preparation). Furhtermore, neglecting possible obscura- quasars discovered to date: J1148+5251 (z = 6.42),
Q
tioneffectsarisingfromopticallythickgasanddustclose J1030+0524 (z = 6.28), and J1623+3112 (z = 6.22).
Q Q
to the quasarcan bias our determination of the quasar’s These quasars are the only quasars discovered to date
luminosity to lower values. This means that our deter- with z > 6.2, as well as the only quasars (aside from
Q
mination of the quasar’s luminosity should be viewed as the BAL quasar J1048+4637) exhibiting complete GP
theobservedluminosityalongthe LOS,andequaltothe troughs, with no detectable flux over a wide wavelength
intrinsic luminosity in the low opacity limit. range (Fan et al. 2006b). Readers interested in the de-
Finally, we define the neutral hydrogenfraction in the tailsofthecorrespondingobservationsareencouragedto
IGM, xIHGIM, with Eq. (3) using the mean density of the consult Fan et al. (2006b) and White et al. (2003).
universeattheedgeoftheHIIregion,n=n0[1+z(RS)]3, Inordertocomparethequasars’observedspectrawith
and ΓQ = 0. Although ΓBG is the more physically rel- our simulated spectra, we must know the intrinsic emis-
evant quantity and we use it in our analysis, henceforth sionspectrumfromeachofthe observedsources. Obser-
we willuse xIHGIM as the proxyfor ΓBG in orderto match vations of lower redshift quasars (e.g. Telfer et al. 2002)
convention. As defined above, xIGM represents the neu- suggestthattheirintrinsiccontinuumemissioniswellfit
HI
tral fraction at the mean density. The more common by a broken power–law, with a mean spectral slope of
representation of the neutral fraction in the literature is α 0.7(f να)redwardoftheLyαline,andaslope
ν
∼− ∝
inthe formofthevolume(ormass)weightedmean,x¯ . ofα 1.8bluewardof the Lyα line. The high-redshift
HI
∼−
This definition, of course, depends on the assumed den- SDSS quasars do not show any evidence for deviating
sity distribution. Hence, for comparison purposes, we from this trend (e.g. White et al. 2003). The precise
express our results both in terms of xIGM and x¯ , with location of the wavelength of the break in the power–
HI HI
the latterobtainedby averagingoveroursimulationbox law is difficult to determine for any given source, due
densities and assuming ionization equilibrium. to a strong broad Lyα emission line superimposed on
To summarize, our analysis has three free parameters: the continuum spectrum. However,for our purposes, we
5
the fits to the Lyα line of J1148+5251and J1623+3112.
However,the Lyα line ofJ1030+0524requirestwogaus-
sians for an accurate fit, as it contains a broad and nar-
row component, similar to the shape of a Voigt profile
(Ho et al. 1997; Greene & Ho 2004). The redshifted line
centersofthegaussianswereheldfixedinthefittingpro-
cedure, set to (1+z )λ , where λ =1215.67˚A for Lyα
Q 0 0
andλ =1240.81˚AforNV.Theheightandwidthofthe
0
gaussianswerethenfoundbyfitting the fluxinpixelson
the red side of the line center. Care was taken, in addi-
tion, to excise the spectral region immediately redward
of the Lyα line center, as this region could be suscepti-
ble to non-negligibleabsorptionby dense infalling gas in
the foreground, close to the quasar (e.g. Mesinger et al.
2004; Barkana & Loeb 2004a). Indeed, the quasar spec-
tra do exhibit sharp drops in the observed flux as one
approaches the line center from the red side (see Fig. 2;
see also Barkana & Loeb 2003). This region was excised
by eye from the fitting procedure. The amplitude of the
continuumwasfoundseparately,byfittingthefluxinre-
gionsofthespectrumfurthertothered,thatarea-priori
knownnottosuffermetallineabsorption,corresponding
to roughly 1300 – 1350 ˚A in the rest frame. We explic-
itly avoid spectral sub–regions with obvious emission or
absorption lines. More specifically, the spectral regions
Fig. 1.— Model from a hydrodynamical simulation for the in the observed frame used in the fit for the continuum
opticaldepthcontributionsfromwithin(τR)andfromoutside(τD)
the local ionized region for a typical lineof sight towards a zQ = were (demarcated with horizontal lines in Figure 2):
6.42quasar. Thisfigure was created withparameter choices (RS,
xanIHGdIMth,efΓs)ol=id(c4u0rvMepcco,rr1e,s1p)o.nTdshetodaτRsh.eTdhceurtvoetacloLrryemspaonndαsotpotiτcDal, • J9510104˚A8+,59275010:–10090100˚A0–;9125˚A, 9128–9190˚A, 9220–
depth is the sum of these two contributions, τ = τR + τD. The
dotted linedemarcates roughlythemaximumtotal optical depth, J1030+0524: 8880–8990˚A, 9020–9100˚A, 9840–
tτh∼an51.5-,σa,tinwthhicehsapescigtrnaallrceagniobnejduesttebctluedewwairtdhoafctohnefiLdeynαceemgriesasitoenr • 9930˚A, 9960–10020˚A;
line(seediscussionin§2.3). Inouranalysis,thetransmission(e−τ)
isaveragedover∼3–4˚Awavelengthbins(1/10oftheresolution J1623+3112: 8830–8895˚A, 8980–9040˚A, 9790–
inthefigureabove),thusdecreasingthefluctuationintheeffective • 9863˚A, 9880–9907˚A;
optical depth. For reference, the redshifted Lyα wavelength is at
9020˚A,offtheplot.
InFigure2,weshowthedatafortheobservedspectra,
F(λ ), with 1-σ error bars as well as the continuum +
obs
areinterestedinthe opticaldepthatwavelengthswithin Lyα + N V fit to the intrinsic emission, F0(λobs), gen-
just a few tens of rest-frame Angstroms blueward of the erated by our procedure (solid curve). The Lyα optical
Lyαlinecenter,wheretheLyαemissionlineisstilldom- depth can then be inferred in each pixel separately as:
inant over the continuum. Therefore, we do not model τobs(λobs)= ln[F(λobs)/F0(λobs)].
−
thebreak,anduseasinglepower–lawcontinuumcompo-
nent instead in our emission template. Throughout our 2.3. Comparing Simulated and Observed Spectra
analysis, we keep the logarithmic slope α = 0.7 fixed After the observed optical depths, τ (λ ), are ob-
− obs obs
(aswewilldemonstratebelow,ourresultsareinsensitive tained as discussed in 2.2, and the simulated optical
to this choice). depths,τ (λ ),areo§btainedforeveryLOSandpoint
sim obs
We deriveourconstraintsbelow entirelyfromthe blue in our parameter space as discussed in 2.1, we wish
side of the Lyα line in our analysis, where the effects of to statistically compare the simulated an§d observed op-
theIGMcanbefeltmoststrongly(Mesinger et al.2004). tical depths for each quasar. More precisely, for each
To this end, we make use of “clean” spectral regions on quasarandeverypoint in our three–dimensionalparam-
the redsideofthe Lyαline,wherethe resonantattenua- eter space, (R , xIGM, f ), we want to answer the ques-
S HI Γ
utiaotnioins niseglelisgsibthlean(τR10∼%0()e,−aτnDd>th0e.9d)aemvpeninignwthinegcaatsteeno-f utieosn:anwdhathteislisthteofliτkelih(oλod t)hvaatlutheeslwisetreofdτroabwsn(λforbosm) vtahle-
sim obs
afullyneutralIGMandasmallRS. WefittheNVemis- same underlying probability distribution? Note that we
sionlineat(1+zQ)1240.81˚Awithasinglegaussian,and have a fixed list of about 30 τobs(λobs) values for each
thecontinuumwithasinglepower–law. Additionally,we quasar,butwegenerateadifferentlistof3000τ (λ )
sim obs
fit the Lyα line with a single gaussian for J1148+5251 values, combining 100 LOSs, for each choice of (R ,
S
((wwiitthh aa fifitttteeddwwiiddtthhooff∆∆λλoobbss==8111˚A3˚A) )anadndwiJt1h62a3s+u3m11o2f xIHGWIMe, fuΓs)e. the Kolmogorov-Smirnov (K-S) test (e.g.
two independent gaussians for J1030+0524 (with fitted Press et al.1992)to comparethe cumulativeprobability
widths of ∆λ = 45˚A and ∆λ = 113˚A). We have distributions (CPDFs) – i.e., histograms – of τ (λ )
obs obs obs obs
found that adding a second gaussian does not improve and τ (λ ). The usefulness of the K-S test is that
sim obs
6
4 6 10
3 8
4
2 6
1 2 4
2
0 0 0
8900 8950 9000 8750 8800 8700 8720 8740
Fig. 2.—Observedspectra(shownaspointswith1-σerrorbars)andthecontinuum+Lyα+NVfitoftheintrinsicemissionspectrum
generated byourprocedure(solid curve). Thespectralregions(ontheredsideoftheLyαlinecenter)that weusedinthefitaremarked
with horizontal lines. The inlaid panels zoom in on the spectral region (on the blue side of the Lyα line center) that were then used to
compare the absorption optical depth τ with hydrodynamical simulations. The vertical dashed lines bracket the wavelength ranges used
toconstructτ histograms(see§2.3). ThespectrashownarethoseofJ1148+5251, J1030+0524, J1623+3112 (left to right).
it provides a figure of merit which can be directly in- formation from the wavelength dependence is fully pre-
terpreted as the likelihood that two distributions were served; however,only a single observedpixel is available
drawn from the same underlying distribution, without to compare with a corresponding pixel from each LOS,
theneedforthead-hocbinningofdata,andtheresulting rendering the use of the K-S test unreliable. The con-
loss of information, required by such methods as the χ2 verse is true in the other extreme, where the wavelength
test. ThedrawbackofusingtheK-Stestisthatcomput- binsizeapproachestheentirespectralregionweanalyze:
ing absolute confidence levels around our best fit values the shape of the CPDF is accurately captured, but in-
would require prohibitively expensive Monte-Carlo sim- formationfromthewavelength–dependenceisdestroyed.
ulations. We caution the reader that, as a result, the Here we (somewhat arbitrarily) divide the spectral re-
relative likelihoods we quote below do not representtra- gion of analysis into wavelength bins of width 25–30
ditional “error–bars”. ˚A, so that each bin contains between 6–8 pixe∼ls. The
Our statistical analysis is quite similar to that em- two blue–most bins were chosen bracket the edge of the
ployed by Rollinde et al. (2005) to model density struc- GP trough in the observed spectrum. Furthermore, we
turesurroundingmoderateredshift(z 2)quasars. The remove the region 20–30 ˚A blueward of the line center
∼
main difference in technique is that while our present from our analysis, corresponding to the expected mean
approach uses a parameterized model, Rollinde et al. radius of the large–scale overdense region surrounding
(2005) fit for the intrinsic spectrum directly from the such quasars (Barkana & Loeb 2004a). This step was
low–τ regions of the spectrum. Using simulated spec- necessary,as our simulation box size is too smallto con-
tra, they show that by applying the K-S test on optical tain a statistical sample (or even a single occurrence) of
depth CPDFs, one can statistically discriminate among suchlargeoverdensities,requiredfortheiraccuratemod-
different models of the distance–evolution of the optical eling. The wavelength bins we use in our analysis are
depth. In Mesinger et al. (2004), we use a similar statis- shownbythe verticaldashedlinesinthe inlaidpanelsof
ticaltesttoshowthatthemoderatedegeneracybetween Figure2. We verifythatourresultsarefairlyinsensitive
R and xIGM can be broken with a single spectrum, if to our choice of bin size, i.e. the relative K-S likelihoods
S HI
R ismoderatelyconstrained. Thisprioranalysisdidnot between models change little with decreasing bin size.
S
take into account the additional constraint arising from In order to incorporate an uncertainty into our emis-
the evolution of the CPDF with distance, as our current sion template, we add Gaussian-distributed, uncorre-
analysis does. Nevertheless, we caution the reader that lated flux variations to each pixel in the template. We
thesetestsonlyconfirmtheabilitytostatisticallyrecover do this by regarding our fitted value of the template,
input parameters using our statistical analysis. They do F (λ ), asthe meanvalue ofthe true flux distribution,
0 obs
not test the validity of our model. andthe“true”emissionatλ isdrawnfromaGaussian
obs
The optical depth CPDF should have a strong depen- distributionaroundthatmean,withastandarddeviation
dence on wavelength (since the gas is more highly ion- σ (λ ). Hence, we replace every τ (λ ) value by
flux obs sim obs
ized near the quasar). In principle, the information con- a set of 1000 values, defined by
tainedinthisdependencecouldbeextractedandusedto
further constrain model parameters. For example, one τsampled(λ ) ln[e−τsim(λobs) A] , (6)
could compare the CPDFs of τ (λ ) and τ (λ ) sim obs ≡− ×
obs obs sim obs
in several wavelength bins for each quasar. Unfortu- where A represents our emission uncertainty, and is
nately, in practice, one also needs to have many points drawn from a Gaussian distribution with a mean of 1,
in each individual bin, to be able to sample the under- and a standard deviation of σ (λ )/F (λ ). We
flux obs 0 obs
lying distribution in that bin. In one extreme, as the set σ (λ )/F (λ ) = 0.3, which corresponds ap-
flux obs 0 obs
wavelength bin size approaches the pixel width, the in- proximately to the typical quasar-to-quasar scatter in
7
the continuum level immediately redward of Lyman α bottom–right with increasing f . Both of these trends
Γ
(Vanden Berk et al. 2001). However, we have verified areevidentinFigure3. Forexample,thecoloredsquares
that our results are not sensitive to this choice. For ex- (correspondingtoK-Slikelihoodsthatarewithinafactor
ample, in the case of J1623+3112, we find that in the of27ofthepeaklikelihood)inthemiddlerowcovermost
range 0.1 < σ (λ )/F (λ ) < 0.5, the peak K-S of the displayed parameter space in the left–most panel,
flux obs 0 obs
probability∼remains unchanged, and∼the confidence lim- though there are no pixels with high likelihoods. In-
its quoted below change at most by a single pixel in our creasingthequasar’sionizingflux(i.e. goingthroughthe
parameter space. panels from left to right), shifts the colored squares in-
Finally,incomparingspectra,oneneedstodefinewhat creasinglytowardsthe bottom–rightcornerofthe panel.
oneconsiderstobe“zeroflux”,assimulationshavenear-
infinite precision and observations do not. To this end, 3.1. J1148+5251
wedefine amaximumallowedvalueoftheopticaldepth, Likelihood estimates for J1148+5251 (z = 6.42) can
Q
MAX[τ] 5.5. This approximatelycorrespondstoa1-σ beseeninthetoprowofthreepanelsinFig. 3. Theleft,
≡
detectionin ourobservedspectra atwavelengthsaround center,andrightpanels correspondto valuesoff = 0.7,
Γ
the edge of the GP trough, where the majority of pixels 1.6, and 1.9, respectively. The peak likelihood of 3.0%
with fluxes less than this value are located (see Fig. 2). occursat(R ,xIGM,f )=(40Mpc,0.2,1.6). Thisisthe
S HI Γ
In our comparison, all values of τ > 5.5 are treated the smallest peak likelihood we obtain (c.f. peak likelihoods
same, regardless of their actual value. of 30% for the two quasars below), and most likely
To summarize, for each quasar we divide the spectral mea∼ns that our model isn’t as accurate in describing the
region roughly blueward of (1+zQ)λα–25˚A and red- physical spectra of this particular quasar. The volume–
∼
ward of ∼ [1+z(RS)]λα–25˚A into 3–5 bins of approx- averagedneutralfractionat the peak is x¯HI =0.16. The
imately equal wavelength width, containing 6–8 pixels rangeofparametervalueswhoselikelihoodestimatesare
each. In each bin we compare the CPDFs of τobs(λobs) withinafactorof3ofthepeaklikelihoodisencompassed
andτsample(λ )(thelatterconcatenatedover100simu- by:
sim obs
latedLOSs). Foreverypointinour3Dparameterspace,
40 Mpc <R < 42 Mpc,
these two distributions are compared via the K-S test, • S
∼ ∼
which produces a probability to be interpreted as the
xIGM < 1.0 or x¯ < 1.0,
likelihood that the two distributions were drawn from • HI HI
∼ ∼
the same underlying distribution. Finally, the K-Sprob-
0.7 <f < 1.9.
abilities for each wavelength bin are multiplied to give a • Γ
∼ ∼
likelihood estimate for that point in the 3D parameter
Note that here as well as below, the end points of the
space.
quoted ranges correspond to the minimum or maximum
value of the parameterfor whichthere exists a combina-
3. RESULTS tionoftheothertwoparametersgivingaK-Sprobability
Our results for the likelihoods of various parameter– within a factor of three from the peak probability. Note
combinations are shown in Figure 3 for J1148+5251 that these ranges are similar in spirit to “marginalized
(top), J1030+0524 (middle), and J1623+3112 (bottom). errors”,inthattheotherparametersareallowedtovary,
Thecrossineachcasemarksthepointinthe3Dparam- rather than fixed (but they are not actual marginalized
eter space where the likelihood peaks. In order from errors;see caveat mentioned in the footnote above).
lighter to darker, the yellow, green, and blue squares
correspond to points in our parameter space with like- 3.2. J1030+0524
lihoods that are a factor of 27, 9, and 3 below the peak Likelihood estimates for J1030+0524 (z = 6.28) are
Q
likelihood, respectively.6 Parameter combinations with shown in the middle row of panels in Figure 3. The left,
likelihoodssmallerthan1/27thofthe peak areshownas center,andrightpanels correspondto valuesoff = 0.4,
Γ
empty squares. 1.0, and 1.3, respectively. The peak likelihood of 34%
Below, we will discuss each quasar in turn. However, occurs at (R , xIGM, f ) = (41 Mpc, 1, 1.0). The range
S HI Γ
it is encouraging to note first that the isocontours in ofparametervalueswhoselikelihoodestimatesarewithin
Figure 3 behave as expected, following the moderate de- a factor of 3 of the peak likelihood is encompassed by:
generacies in the parameters. Namely, since a stronger
damping wing (τ ) contribution can be obtained either 37 Mpc <R < 45 Mpc,
D S
•
by increasing the neutral fraction or by decreasing the ∼ ∼
sizeoftheHII region,theisocontoursshouldbe oriented xIGM < 1.0 or x¯ < 1.0,
• HI HI
frombottom-left to upper-rightineachpanel. Similarly, ∼ ∼
sinceasmallerionizingflux(i.e. greaterτR)canroughly • 0.7 <fΓ < 1.9.
be compensated for in the total optical depth, τ +τ , ∼ ∼
R D
by increasing the size of the HII region or decreasing As mentioned earlier, J1030+0524 has a unique (at
the neutral fraction (both decrease τ ), the isocontours least among our sample of three quasars) feature.
D
should shift monotonically from the upper–left to the Namely,theLyβGPthroughisnoticeablyoffsetfromthe
LyαGPthrough. Themerepresenceofflux(irrespective
6 In order to obtain absolute confidence limits, we would need of its properties) in the Lyβ region of the spectrum can
to performa set of computationally expensive Monte-Carlosimu- be used to set a minimum size of the surrounding HII
lations. For simplicity, we avoid such computations, and instead region,R >40 Mpc. The regionexcluded by this prior
contendourselvesbyquotingresultsintermsoftherelativeoffset S
fromthepeaklikelihood. isshadedov∼erinthemiddlerowinFigure3. Takingthis
8
Fig. 3.— Likelihoods for J1148+5251 (top), J1030+0524 (middle), and J1623+3112 (bottom) in the 3–dimensional parameter space of
(RS, xIHGIM, fΓ). Each panel shows contours ina2D sliceof the 3D likelihood surface. Thecross marks the parameter combination with
the peak likelihood value (see § 2.3 for how the likelihoods were derived). In order from lighter to darker, the yellow, green, and blue
squares correspond to points inour parameter space withlikelihoods that are afactor of 27, 9, and 3below the peak value, respectively.
Parametercombinationswithlikelihoodsbelow1/27thofthepeakareshownasemptysquares. J1148+5251: Thepeaklikelihoodof3.0%
occurs at (RS, xIHGIM, fΓ) = (40 Mpc, 0.2, 1.6). The left, center, and right panels correspond to values of fΓ= 0.7, 1.6, 1.9, respectively.
J1030+0524: Thepeaklikelihoodvalueof34%occursat(RS,xIHGIM,fΓ)=(41Mpc,1,1.0). Theleft,center,andrightpanelscorrespond
tovaluesoffΓ=0.4,1.0,1.3,respectively. TheshadedregionshowsvaluesruledoutbythepresenceoffluxintheLymanβ region,which
setstheadditionalconstraintRS >∼40Mpc. J1623+3112: Thepeaklikelihoodvalueof39%isat(RS,xIHGIM,fΓ)=(29Mpc,1,0.7). The
left,center,andrightpanelscorrespondtovaluesoffΓ=0.4,0.7,1.0,respectively.
lowerlimitonR intoaccount,therangeofvalueswhose 3.3. J1623+3112
S
likelihood estimates are within a factor of 3 of the peak
Likelihood estimates for J1623+3112 (z = 6.22) are
Q
likelihood shrinks to:
shown in the bottom row of panels in Figure 3. The
left,center,andrightpanelscorrespondtovaluesoff =
41 Mpc <R < 45 Mpc, Γ
• S 0.4, 0.7, and 1.0, respectively. The peak likelihood of
∼ ∼
39% occurs at (R , xIGM, f ) = (29 Mpc, 1, 0.7). The
0.04 <xIGM or 0.033 <x¯ , S HI Γ
• HI HI rangeofparametervalueswhoselikelihoodestimatesare
∼ ∼
withinafactorof3ofthepeaklikelihoodisencompassed
0.7 <fΓ < 1.3. by:
•
∼ ∼
9
solidblue,long-dashedred,short-dashedgreen,anddot-
tedpurple curvescorrespondtoparametercombinations
of (R , xIGM, f ) = (29 Mpc, 1.0, 0.7), (29 Mpc, 0.04,
S HI Γ
0.7), (33 Mpc, 1.0, 0.7), and (29 Mpc, 1.0, 0.1), respec-
tively. The crosses are located at the set of observed
optical depth values. The solid blue curve corresponds
to the parameter combination with the peak likelihood
value, while the parameter combinations represented by
the other curves have likelihood values below 1/27th of
the peak. The sharp jump corresponding at the largest
τ binandtheoverlayingofsevencrossesinthetoppanel
areduetothefinitedetectionthreshold,wherebyallval-
ues of τ > 5.5 are indistinguishable from one another,
and hence fall in the same bin in this figure. The fig-
ure clearly demonstrates that our lower limit on xIGM
HI
arises from the low–τ tail predicted in models with low
neutralfractionsorthosecorrespondingtoregionsofpa-
rameter space which could mimic a low neutral fraction
(seethediscussiononisocontoursin 3). Whilethesolid
§
bluecurvesareconsistentwiththedistributionofpoints,
these other curves all have tails that are “too long”: no
points corresponding to these tails are seen in the data.
Note also that while some models might perform well in
onewavelengthbin,onlythesolidbluecurveisconsistent
with the observed data in all wavelength bins.
Fig. 4.— Distributions of simulated optical depths for the
three wavelength bins used in our analysis of J1623+3112. The 4.1. Sanity Checks
solidblue,long-dashedred,short-dashedgreen,anddottedpurple
curvescorrespondtoparametercombinationsof(RS,xIHGIM,fΓ)= As our analysis technique is as yet untested, it is use-
(29Mpc,1.0,0.7),(29Mpc,0.04,0.7),(33Mpc,1.0,0.7),and(29 ful to take a step back, and verify that our results make
Mpc, 1.0, 0.1), respectively. The crosses are located at the set of sense. Asalreadymentioned,thelikelihoodiso–contours
observedopticaldepthvalues. Thesolidbluecurvecorrespondsto
behave as expected. Here we discuss several other en-
the parameter combination with the peak likelihood value, while
the parameter combinations represented by the other curves have couragingchecksonthevalidityofouranalysis. Wecau-
likelihoodvaluesbelow1/27thofthepeak. Thesharpjumpatthe tion the reader, however, that these are not meant to
largestτbinandtheoverlayingofsevencrossesinthetoppanelare be a proof of the accuracy of our analysis, and instead
duetothefinitedetectionthreshold,wherebyallvaluesofτ >5.5
should be regarded merely as consistency checks.
areindistinguishablefromoneanother,andhencefallinthesame
bin in this figure. The figure demonstrates that our lower limit First, the procedure here recovers fairly well the re-
on xIGM comes from the low–τ tail predicted in models with low sult of our previous analysis in the case of J1030+0524
HI
neutral fractions or those corresponding to regions of parameter (Mesinger & Haiman 2004). This need not be the case,
spacewhichcouldmimicalowneutralfraction(seethediscussion
since the nature of the two analysis procedures are dif-
on isocontours in § 3). No points corresponding to these tails are
seeninthedata. Notealsothatwhilesomemodelsmightperform ferent; namely, our previous work focused only on gross,
wellinone wavelength bin, only the solidbluecurve isconsistent approximateLyαandLyβ transmissionfeaturesnearthe
withtheobserveddatainallwavelengthbins. edge of the HII region. Nevertheless, in both cases, the
peak likelihood occurs at the same value of xIGM and
HI
25 Mpc <RS < 29 Mpc, RS, with only the value of fΓ being somewhat larger
• ∼ ∼ here (fΓ=1.0 instead of fΓ=0.4 in Mesinger & Haiman
• 0.04 <∼xIHGIM or 0.033 <∼x¯HI, 2cl0u0s4io).nWofewaatvtreilbeuntgeththbisinsshfiuftrtihnefrΓinpwriamrdarfirlyomtotthheeHinII-
0.4 <f < 1.3. region edge, where the effects of R and xIGM are felt
• Γ S HI
∼ ∼ less strongly and f wields the most statistical weight.
Γ
In particular,the likelihoods in the two blue-most wave-
4. DISCUSSION length bins are comparable for f = 1.0 and f = 0.4
Γ Γ
In this section, we discuss a few important aspects of (given xIGM =1.0 and R =41Mpc); however, the like-
HI S
the above results. First, it is encouraging that the sim- lihoods in the two red-most wavelength bins are 7–8
∼
ple methodology we outlined can deliver statistically in- times greater for f =1.0 than f =0.4.
Γ Γ
teresting constraints, even from the few dozen available Second, our results on the neutral fraction agree in
spectral pixels. Clearly, it will be interesting to apply relative terms with estimates obtained from transmis-
a similar analysis to a larger sample of sources in the sion in GP troughs (Fan et al. 2006a). J1030+0524and
future. J1623+3112 both have very dark GP troughs in Lyα
Perhaps our most interesting result is the lower limit and Lyβ, and thus are only able to provide lower limits
we obtain, for two of the quasars, on the neutral frac- onthe opticaldepth andcorrespondingneutralfraction;
tion. We would like to thereforeunderstand where these however, J1148+5251 has transmission gaps in the GP
constraints actually come from. To this end, in Figure 4 troughs,thus favoringasmalleropticaldepthandcorre-
we show the simulated optical depth PDFs for the three spondingneutralfraction. Thisisqualitativelysimilarto
wavelengthbinsusedinouranalysisofJ1623+3112. The our results: we were unable to constrain xIGM from the
HI
10
spectra of J1148+5251 obtaining only a shallow peak in the Lyα line width. This does not, of course, preclude
the likelihoodatxIGM =0.2;however,bothJ1030+0524 the possibility that the intrinsic Lyα lines of the z 6
HI ≈
and J1623+3112 had peak likelihoods at xIGM = 1, and quasars possess some deviations from a gaussian shape.
HI
strongly prefer xIGM >0.04. For example, if there is an asymmetry, in the sense of
HI a steeper drop of flux on the blue side of the line than
Third, the maximu∼m likelihoods of f (i.e. the rela-
Γ on the red side, our technique would have mistakenly
tive strength of the quasars’ ionizing luminosity) match
attributed this steeper drop to absorption. This implies
therelativestrengthsofthequasars’continuaredwardof
that we might have overestimated the neutral fractions.
Lyα. The values of f corresponding to the peak likeli-
Γ At present, there is, however, no evidence of such an
hoods are 1.6, 1.0, and 0.7 for J1148+5251 J1030+0524
asymmetry in the emission lines of low redshift quasars
and J1623+3112 respectively. Looking at Figure 2,
(where both sides of the lines can be seen unabsorbed),
one can note that this is also the order of the relative
nor is there evidence, from the red side of the lines, that
strengthsofthequasars’continuaredwardofLyα,asone
thehigh–zquasars’lineshaveshapesdifferentfromthose
wouldexpectiftheUVpower-lawindicesandamountof
of their lower–redshiftcounterparts.
obscuration didn’t vary greatly among the quasars.
Finally, we also note that our spectral fits, describing
Finally, the quasar with the smallest HII region,
the intrinsic emission line on the unabsorbed red side
J1623+3112 with a peak likelihood at R = 29 Mpc,
S of the line, leave typical flux residuals of only 1%.
also has the highest favored value of xIHGIM = 1.0 and Similarresidualsarepresumablypresentontheblu∼eside
the lowest favored value of f . This makes sense, since
Γ of the line, as well. However, the absorption exp( τ )
a fainter quasar, embedded in a more neutral medium − D
causedbytheIGMdampingwingonthebluesideofthe
shouldhaveasmallerHIIregion,givenasimilarlifetime.
line (see Fig. 1), has an amplitude of (0.1 1)x ,
HI
exceeding the residuals for x > 0.01.∼Furthe−rm−ore, in
HI
4.2. Uncertainties in the Intrinsic Emission Spectra orderfor the residualsto mimic∼the IGMdamping wing,
Althoughwehaveaddedgaussiannoisetoourinferred theywouldhavetobecoherent(monotonicallyincreasing
template spectrum, as described above, we wish to get with wavelength).
a sense of how correlated errors, modifying the broader
shapeofthetemplatespectrum,wouldaffectourresults.
4.3. Biases in Local Density and Velocity Fields
The general concern is that mis–estimates of the flux,
correlatedovermany pixels, couldpotentially mimic the As mentioned above, we remove the region 20–30 ˚A
effect of shifts in our parameter space. To get a sense of bluewardofthelinecenterfromouranalysis,correspond-
the magnitude of this effect, we chose to vary the width ing to the expected mean radius of the large–scale over-
ofthefittedLyαgaussianforquasarJ1623+3112andsee denseregionsurroundingsuchquasars(Barkana & Loeb
how this modifies our results. Naively,one wouldexpect 2004a). This is a necessary step as present-day simu-
a narrower Lyα emission line and an accompanying de- lations cannot statistically model such rare, large–scale
crease in effective τ to cause a shift in the parameter overdensities. However, if the density bias extended
obs
space to regions favoring smaller optical depths (larger into our region of analysis, our inferred xIHGIM lower lim-
R , smaller xIGM, larger f ), with the opposite trends its would be conservative and our inferred fΓ values
S HI Γ would most likely be underestimated. If the source is
for a wider Lyα emission line.
hosted by a large–scale overdensity, then density would
We find that such degeneracies are very weak. For a
decrease with decreasing λ (increasing distance from
10% narrower line, we find that the peak likelihood re- obs
the source). If such a decreasing density extended into
mains at the same parameter–combinations. The range
ourregionofanalysis,itwouldtranslateintoashallower
ofparametervalueswhoselikelihoodestimatesarewithin
rise in τ n2/Γ. In matching the total optical depth
afactorof3ofthepeaklikelihoodchangesbyatmostone R
∝
pixelinour3Dparametergrid: 25Mpc<R <29Mpc; τLyα = τR +τD, our analysis, which assumes constant
S
mean density, would translate this shallower rise in τ
0.008 < xIGM < 1.0; 0.4 < f < 1.6. ∼This i∼s a rather R
HI Γ to a shallower rise in τ , i.e. a lower value of xIGM than
small c∼hange, a∼nd we emp∼hasiz∼e that a 10% narrower appropriate if the densDity bias was taken into aHcIcount7.
Lyα line would actually be an unacceptably poor fit to Note that our results on xIGM are driven by the wave-
the(redsideofthe)observedspectrum. AwiderLyαline HI
length dependence of τ; a constant offset in the density
ismorephysicallymotivated, sincestrongdampingwing
field can be absorbed by a shift in Γ (i.e. f ), which
Q Γ
absorptionmightbe able tosomewhatattenuateflux far
dominates over Γ throughout most of our region of
BG
into the red side of the line. Hence, we vary the width
interest and is predominantly set by the amount of ab-
of the Lyα line, well in excess of the absorbing potential
sorptioninbinsfurtherfromtheHIIregionedgewhereit
of a strong damping wing. Again, we find little change
has the most statistical weight. Furthermore, it is possi-
in our results. Even for a 50% wider line, the param-
ble that an infalling IGM could introduce a similar bias,
eter combination at peak likelihood remains the same,
buttheeffectislikelytobesmallforabrightquasarand
while the range of parameter values whose likelihood es-
far away from the line center on the blue side. More-
timates are within a factor of 3 of the peak likelihood
over,the effect from the infall would have the same sign
changes by at most one pixel to: 25 Mpc < R < 31
S as the density bias above, making the flux fall off on
Mpc; 0.2<xIGM <1.0;0.4<f <1.3. Ase∼xpected∼,an
HI Γ
intrinsical∼ly wider∼Lyα line s∼hows∼a slight preference for 7 Note that a shallower rise in τD could also be achieved if we
a more neutral universe. were overestimating RS, instead of underestimating xIHGIM. How-
ever,our inferredRS values arealreadyatthelowestlimitsetby
Overall, the above exercise suggests that our results
the onset of GP troughs for two of our quasars and thus can not
are robust and conservative,at least to mis–estimates of bedecreasedfurther.