How Much Information Does the Sodium Doublet Encode? Retrieval Analysis   of Non-LTE Sodium Lines at Low and High Spectral Resolutions

Chloe Fisher; Kevin Heng

arXiv:1906.07035·astro-ph.EP·August 14, 2019

How Much Information Does the Sodium Doublet Encode? Retrieval Analysis of Non-LTE Sodium Lines at Low and High Spectral Resolutions

Chloe Fisher, Kevin Heng

PDF

TL;DR

This study investigates how much atmospheric information the sodium doublet lines encode in exoplanet transmission spectra, revealing limitations in current retrieval methods at different spectral resolutions and the challenges in distinguishing LTE from NLTE conditions.

Contribution

The paper provides a comprehensive analysis of sodium line retrievals in exoplanet atmospheres, highlighting the limitations of current low- and high-resolution spectra in constraining atmospheric parameters and LTE/NLTE states.

Findings

01

High-resolution spectra alone cannot precisely determine atmospheric pressure levels.

02

Low-resolution retrievals face degeneracies due to Rayleigh scattering and normalization issues.

03

Current data cannot reliably distinguish LTE from NLTE sodium line interpretations.

Abstract

Motivated by both ground- and space-based detections of the sodium doublet in the transmission spectra of exoplanetary atmospheres, we revisit the theory and interpretation of sodium lines in non-local thermodynamic equilibrium (NLTE), where collisions are not efficient enough to maintain a Boltzmann distribution for the excited and ground states of the sodium atom. We consider non-Boltzmann distributions that account for the ineffectiveness of collisions. We analyze the sodium doublet in transmission spectra measured at low (HAT-P-1b, HAT-P-12b, HD 189733b, WASP-6b, WASP-17b and WASP-39b) and high (WASP-49b) spectral resolutions. Nested-sampling retrievals performed on low-resolution optical/visible transmission spectra are unable to break the normalization degeneracy if the spectral continuum is associated with Rayleigh scattering by small cloud particles. Using mock retrievals, we…

Tables3

Table 1. Table 1: Relevant physical properties of Na D lines

\begin{matrix} Quantity & Value \\ λ_{D_{1}} ​ (Å) & 5895.92424 \\ ν_{D_{1}} ​ (Hz) & 5.085 \times 10^{14} \\ λ_{D_{2}} ​ (Å) & 5889.95095 \\ ν_{D_{2}} ​ (Hz) & 5.090 \times 10^{14} \\ A_{D_{1}} ​ (s^{- 1}) & 6.14 \times 10^{7} \\ A_{D_{2}} ​ (s^{- 1}) & 6.16 \times 10^{7} \\ E_{D_{1}} ​ (erg) & 3.369 \times 10^{- 12} \\ E_{D_{1}} ​ (eV) & 2.1 \\ E_{D_{2}} ​ (erg) & 3.373 \times 10^{- 12} \\ E_{D_{2}} ​ (eV) & 2.1 \\ f_{D_{1}} & 0.320 \\ f_{D_{2}} & 0.641 \end{matrix}

Table 2. Table 2: Retrieved parameters and their prior ranges

Quantity	Symbol	Units	High-res Range	Low-res Range	Prior Type
Temperature	$T$	K	0–10,000	0–6000	uniform
Sodium volume mixing ratio	$X_{Na}$	–	$10^{- 13}$ –1	$10^{- 13}$ –1	log-uniform
Broadening parameter	$f_{broad}$	–	$10^{- 1}$ – $10^{4}$	$10^{- 5}$ – $10^{10}$	log-uniform
Gray cloud opacity	$κ_{cloud}$	cm² g^-1	$10^{- 4}$ – $10^{5}$	$10^{- 10}$ – $10^{5}$	log-uniform
Rayleigh scattering scaling parameter (clouds)	$f_{cloud}$	–	–	1– $10^{7}$	log-uniform
Reference pressure	$P_{0}$	bar	1–100	1–100	log-uniform

Table 3. Table 3: Summary of retrieval outcomes (both LTE and NLTE models)

Object	$T_{eq}$ (K)	$T$ (K)	$\log X_{Na}$	$\log f_{broad}$	$\log κ_{cloud}$ (cm² g^-1)	$\log f_{cloud}$	$C$	$\log P_{0}$ (bar)	$R_{0}$ ( $R_{J}$ )
HAT-P-1b $†$	1320	351 ${}_{- 172}^{+ 200}$	-8.11 ${}_{- 3.15}^{+ 4.2}$	0.26 ${}_{- 3.31}^{+ 4.67}$	-6.32 ${}_{- 2.34}^{+ 2.75}$	0.95 ${}_{- 0.73}^{+ 2.67}$	1.18 ${}_{- 0.16}^{+ 1.67}$	0.5 ${}_{- 0.35}^{+ 0.75}$	1.26
NLTE		339 ${}_{- 152}^{+ 185}$	-8.54 ${}_{- 2.85}^{+ 4.03}$	0.2 ${}_{- 3.35}^{+ 4.89}$	-6.4 ${}_{- 2.35}^{+ 2.96}$	1.0 ${}_{- 0.76}^{+ 2.35}$	1.22 ${}_{- 0.2}^{+ 1.81}$	0.5 ${}_{- 0.35}^{+ 0.75}$	1.26
HAT-P-12b $†$	960	984 ${}_{- 302}^{+ 308}$	-9.61 ${}_{- 2.24}^{+ 3.59}$	-0.56 ${}_{- 2.9}^{+ 4.63}$	-6.57 ${}_{- 2.15}^{+ 2.1}$	0.85 ${}_{- 0.55}^{+ 1.13}$	1.18 ${}_{- 0.17}^{+ 1.01}$	0.57 ${}_{- 0.38}^{+ 0.64}$	0.9
NLTE		1018 ${}_{- 283}^{+ 278}$	-9.67 ${}_{- 2.24}^{+ 3.57}$	-0.49 ${}_{- 2.84}^{+ 4.43}$	-6.74 ${}_{- 2.09}^{+ 2.13}$	0.72 ${}_{- 0.44}^{+ 0.86}$	1.24 ${}_{- 0.21}^{+ 0.86}$	0.57 ${}_{- 0.38}^{+ 0.64}$	0.9
HD189733b $†$	1200	2797 ${}_{- 211}^{+ 208}$	-5.44 ${}_{- 1.79}^{+ 1.77}$	2.24 ${}_{- 1.86}^{+ 1.94}$	-5.37 ${}_{- 3.11}^{+ 3.65}$	4.78 ${}_{- 0.75}^{+ 0.84}$	4.78 ${}_{- 0.73}^{+ 0.99}$	1.04 ${}_{- 0.67}^{+ 0.65}$	1.12
NLTE		4782 ${}_{- 472}^{+ 467}$	-6.29 ${}_{- 2.41}^{+ 2.36}$	1.2 ${}_{- 2.35}^{+ 2.51}$	-5.97 ${}_{- 2.67}^{+ 2.9}$	2.92 ${}_{- 0.58}^{+ 0.63}$	5.92 ${}_{- 1.03}^{+ 1.39}$	1.04 ${}_{- 0.67}^{+ 0.65}$	1.12
WASP-6b $†$	1150	880 ${}_{- 312}^{+ 362}$	-9.01 ${}_{- 2.53}^{+ 3.84}$	0.07 ${}_{- 3.26}^{+ 5.21}$	-5.7 ${}_{- 2.72}^{+ 2.97}$	2.55 ${}_{- 1.48}^{+ 2.22}$	1.63 ${}_{- 0.6}^{+ 3.64}$	0.89 ${}_{- 0.58}^{+ 0.68}$	1.18
NLTE		1004 ${}_{- 287}^{+ 272}$	-8.94 ${}_{- 2.48}^{+ 3.77}$	0.04 ${}_{- 3.21}^{+ 5.32}$	-5.89 ${}_{- 2.58}^{+ 2.8}$	1.95 ${}_{- 1.11}^{+ 2.31}$	1.91 ${}_{- 0.82}^{+ 4.6}$	0.89 ${}_{- 0.58}^{+ 0.68}$	1.18
WASP-17b $†$	1740	1489 ${}_{- 399}^{+ 336}$	-9.32 ${}_{- 2.37}^{+ 3.74}$	0.38 ${}_{- 3.51}^{+ 4.45}$	-6.91 ${}_{- 2.01}^{+ 2.39}$	0.75 ${}_{- 0.47}^{+ 0.91}$	1.2 ${}_{- 0.18}^{+ 0.71}$	0.58 ${}_{- 0.36}^{+ 0.54}$	1.73
NLTE		1568 ${}_{- 329}^{+ 307}$	-9.25 ${}_{- 2.45}^{+ 3.77}$	-0.36 ${}_{- 3.04}^{+ 4.67}$	-6.78 ${}_{- 2.1}^{+ 2.24}$	0.58 ${}_{- 0.38}^{+ 0.57}$	1.2 ${}_{- 0.18}^{+ 0.56}$	0.58 ${}_{- 0.36}^{+ 0.54}$	1.73
WASP-39b $†$	1120	835 ${}_{- 132}^{+ 131}$	-5.56 ${}_{- 2.29}^{+ 2.3}$	2.48 ${}_{- 2.4}^{+ 2.14}$	-6.04 ${}_{- 2.71}^{+ 2.84}$	2.22 ${}_{- 0.89}^{+ 1.05}$	1.39 ${}_{- 0.08}^{+ 0.1}$	1.05 ${}_{- 0.68}^{+ 0.65}$	1.18
NLTE		1057 ${}_{- 125}^{+ 205}$	-5.52 ${}_{- 2.98}^{+ 2.96}$	1.6 ${}_{- 3.06}^{+ 2.68}$	-5.94 ${}_{- 2.7}^{+ 2.31}$	1.48 ${}_{- 0.64}^{+ 0.7}$	1.42 ${}_{- 0.1}^{+ 0.12}$	1.05 ${}_{- 0.68}^{+ 0.65}$	1.18
WASP-49b^♡	1400	$7209_{- 1892}^{+ 1763}$	$- {5.64}_{- 2.90}^{+ 2.55}$	$- {5.70}_{- 2.62}^{+ 3.61}$	$- {0.18}_{- 2.77}^{+ 3.14}$	$-$	${1.52}_{- 0.37}^{+ 0.84}$	${1.02}_{- 0.69}^{+ 0.63}$	$1.198$
NLTE		$8415_{- 1526}^{+ 1020}$	$- {4.44}_{- 2.51}^{+ 2.02}$	$- {6.56}_{- 2.30}^{+ 2.91}$	$- {0.90}_{- 2.06}^{+ 2.82}$	$-$	${1.36}_{- 0.24}^{+ 0.42}$	${1.04}_{- 0.68}^{+ 0.64}$	$1.198$

Equations86

C = \frac{Δ R _{fit, CF}}{Δ R _{fit}},

C = \frac{Δ R _{fit, CF}}{Δ R _{fit}},

\frac{\partial n _{2}}{\partial t} = (n_{e} + n_{total}) n_{1} C_{12} + \frac{n _{1} N _{γ} g _{2} A _{21}}{g _{1}} - (n_{e} + n_{total}) n_{2} C_{21} - n_{2} A_{21} (1 + N_{γ}),

\frac{\partial n _{2}}{\partial t} = (n_{e} + n_{total}) n_{1} C_{12} + \frac{n _{1} N _{γ} g _{2} A _{21}}{g _{1}} - (n_{e} + n_{total}) n_{2} C_{21} - n_{2} A_{21} (1 + N_{γ}),

C_{21} = 8.629 \times 1 0^{- 8} \mbox c m^{3} \mbox s^{- 1} T_{4}^{- 1/2} \frac{Ω _{21}}{g _{2}},

C_{21} = 8.629 \times 1 0^{- 8} \mbox c m^{3} \mbox s^{- 1} T_{4}^{- 1/2} \frac{Ω _{21}}{g _{2}},

N_{γ} = (e^{h ν / k_{B} T} - 1)^{- 1},

N_{γ} = (e^{h ν / k_{B} T} - 1)^{- 1},

\frac{n _{2}}{n _{1}} = \frac{C _{12}}{C _{21}} (1 + f_{e}) (1 + f_{e} + \frac{n _{crit}}{n _{total}})^{- 1},

\frac{n _{2}}{n _{1}} = \frac{C _{12}}{C _{21}} (1 + f_{e}) (1 + f_{e} + \frac{n _{crit}}{n _{total}})^{- 1},

\frac{n _{2}}{n _{1}} = \frac{g _{2}}{g _{1}} e^{- E_{12} / k_{B} T} (1 + \frac{n _{crit}}{n _{total}})^{- 1} .

\frac{n _{2}}{n _{1}} = \frac{g _{2}}{g _{1}} e^{- E_{12} / k_{B} T} (1 + \frac{n _{crit}}{n _{total}})^{- 1} .

\frac{C _{12}}{C _{21}} = \frac{g _{2}}{g _{1}} e^{- E_{12} / k_{B} T} .

\frac{C _{12}}{C _{21}} = \frac{g _{2}}{g _{1}} e^{- E_{12} / k_{B} T} .

n_{crit} \equiv \frac{A _{21}}{C _{21}} .

n_{crit} \equiv \frac{A _{21}}{C _{21}} .

σ_{Na} = S Φ.

σ_{Na} = S Φ.

S = \frac{n _{2} A _{21} h ν}{ϵc n _{total}} .

S = \frac{n _{2} A _{21} h ν}{ϵc n _{total}} .

ϵ = \frac{4 π B _{ν}}{c},

ϵ = \frac{4 π B _{ν}}{c},

B_{ν} = \frac{2 h ν ^{3}}{c ^{2}} (e^{h ν / k_{B} T} - 1)^{- 1}

B_{ν} = \frac{2 h ν ^{3}}{c ^{2}} (e^{h ν / k_{B} T} - 1)^{- 1}

A_{21} = \frac{8 π ^{2} e ^{2} ν ^{2}}{m _{e} c ^{3}} \frac{g _{1}}{g _{2}} f_{12},

A_{21} = \frac{8 π ^{2} e ^{2} ν ^{2}}{m _{e} c ^{3}} \frac{g _{1}}{g _{2}} f_{12},

S = \frac{π e ^{2} f _{12}}{m _{e} c} \frac{n _{1}}{n _{total}} (1 - e^{- h ν / k_{B} T}) (1 + \frac{n _{crit}}{n _{total}})^{- 1} .

S = \frac{π e ^{2} f _{12}}{m _{e} c} \frac{n _{1}}{n _{total}} (1 - e^{- h ν / k_{B} T}) (1 + \frac{n _{crit}}{n _{total}})^{- 1} .

S = \frac{π e ^{2} f _{12}}{m _{e} c} (1 + \frac{n _{crit}}{n _{total}})^{- 1},

S = \frac{π e ^{2} f _{12}}{m _{e} c} (1 + \frac{n _{crit}}{n _{total}})^{- 1},

Φ_{L} = \frac{Γ _{L} / π}{( ν - ν _{0} ) ^{2} + Γ _{L}^{2}},

Φ_{L} = \frac{Γ _{L} / π}{( ν - ν _{0} ) ^{2} + Γ _{L}^{2}},

Φ_{D} = Γ_{D}^{- 1} \frac{ln 2}{π} e^{- (ν - ν_{0})^{2} / σ_{th}^{2}} .

Φ_{D} = Γ_{D}^{- 1} \frac{ln 2}{π} e^{- (ν - ν_{0})^{2} / σ_{th}^{2}} .

Φ = \frac{ln 2}{π} \frac{H _{V}}{Γ _{D}},

Φ = \frac{ln 2}{π} \frac{H _{V}}{Γ _{D}},

H_{V} \equiv \frac{a _{0}}{π} \int_{- \infty}^{+ \infty} \frac{e ^{- y^{2}}}{( x - y ) ^{2} + a _{0}^{2}} d y,

H_{V} \equiv \frac{a _{0}}{π} \int_{- \infty}^{+ \infty} \frac{e ^{- y^{2}}}{( x - y ) ^{2} + a _{0}^{2}} d y,

H_{V} = e^{a_{0}^{2} - x^{2}} erfc (a_{0}) cos (2 a_{0} x) + \frac{2}{π} \int_{0}^{x} e^{y^{2} - x^{2}} sin [2 a_{0} (x - y)] d y,

H_{V} = e^{a_{0}^{2} - x^{2}} erfc (a_{0}) cos (2 a_{0} x) + \frac{2}{π} \int_{0}^{x} e^{y^{2} - x^{2}} sin [2 a_{0} (x - y)] d y,

H_{V} = ℜ [w (z)],

H_{V} = ℜ [w (z)],

z = x + i a_{0} .

z = x + i a_{0} .

R = R_{0} + H [γ + ln (\frac{P _{0} κ}{g} \frac{2 π R _{0}}{H})],

R = R_{0} + H [γ + ln (\frac{P _{0} κ}{g} \frac{2 π R _{0}}{H})],

κ = \frac{X _{Na} π e ^{2} f _{12}}{m m _{e} c} Φ + κ_{cloud},

κ = \frac{X _{Na} π e ^{2} f _{12}}{m m _{e} c} Φ + κ_{cloud},

m = 2.4 X_{H_{2}} m_{amu} + X_{Na} m_{Na},

m = 2.4 X_{H_{2}} m_{amu} + X_{Na} m_{Na},

1.1 X_{H_{2}} + X_{Na} = 1,

1.1 X_{H_{2}} + X_{Na} = 1,

κ = κ_{LTE} (1 + \frac{n _{crit} k _{B} T}{P})^{- 1},

κ = κ_{LTE} (1 + \frac{n _{crit} k _{B} T}{P})^{- 1},

κ_{LTE} \equiv \frac{X _{Na} π e ^{2} f _{12}}{m m _{e} c} Φ + κ_{cloud} .

κ_{LTE} \equiv \frac{X _{Na} π e ^{2} f _{12}}{m m _{e} c} Φ + κ_{cloud} .

R = R_{LTE} - H ln (1 + \frac{n _{crit} k _{B} T}{P}),

R = R_{LTE} - H ln (1 + \frac{n _{crit} k _{B} T}{P}),

R_{LTE} \equiv R_{0} + H [γ + ln (\frac{P _{0} κ _{LTE}}{g} \frac{2 π R _{0}}{H})] .

R_{LTE} \equiv R_{0} + H [γ + ln (\frac{P _{0} κ _{LTE}}{g} \frac{2 π R _{0}}{H})] .

Peer Reviews

No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.

Videos

No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.

Full text

How Much Information Does the Sodium Doublet Encode?

Retrieval Analysis of Non-LTE Sodium Lines at Low and High Spectral Resolutions

Chloe Fisher11affiliation: University of Bern, Center for Space and Habitability, Gesellschaftsstrasse 6, CH-3012, Bern, Switzerland. Emails: [email protected], [email protected] 22affiliation: University of Bern International 2021 Ph.D Fellowship

Kevin Heng11affiliation: University of Bern, Center for Space and Habitability, Gesellschaftsstrasse 6, CH-3012, Bern, Switzerland. Emails: [email protected], [email protected]

Abstract

Motivated by both ground- and space-based detections of the sodium doublet in the transmission spectra of exoplanetary atmospheres, we revisit the theory and interpretation of sodium lines in non-local thermodynamic equilibrium (NLTE), where collisions are not efficient enough to maintain a Boltzmann distribution for the excited and ground states of the sodium atom. We consider non-Boltzmann distributions that account for the ineffectiveness of collisions. We analyze the sodium doublet in transmission spectra measured at low (HAT-P-1b, HAT-P-12b, HD 189733b, WASP-6b, WASP-17b and WASP-39b) and high (WASP-49b) spectral resolutions. Nested-sampling retrievals performed on low-resolution optical/visible transmission spectra are unable to break the normalization degeneracy if the spectral continuum is associated with Rayleigh scattering by small cloud particles. Using mock retrievals, we demonstrate that un-normalized ground-based, high-resolution spectra centered on the sodium doublet alone are unable to precisely inform us about the pressure levels probed by the transit chord and hence to identify the region (i.e., thermosphere, exosphere) of the atmosphere being probed. Retrievals performed on the HARPS transmission spectrum of WASP-49b support this conclusion. Generally, we are unable to distinguish between LTE versus NLTE interpretations of the sodium doublet based on the computed Bayesian evidence with the implication that LTE interpretations tend to under-estimate the temperature probed by the transit chord. With the current low-resolution data, the sodium line shapes are consistent with Voigt profiles without the need for sub-Lorentzian wings. The retrieved sodium abundances are consistent with being sub-solar to solar.

planets and satellites: atmospheres

1 Introduction

The use of the resonant lines of the sodium atom to probe the atmospheres of exoplanets has a rich history. The detectability of the resonant doublets of sodium and potassium in transmission spectra was first predicted by Seager & Sasselov (2000). In a contemporaneous study, Sudarsky et al. (2000) reached the same conclusion. Later, sodium was detected for the first time in an exoplanet by Charbonneau et al. (2002) for the hot Jupiter HD 209458b. Since the predictions and discovery, the detection of sodium has become routine from space at low spectral resolution (e.g., Sing et al. 2016) and from the ground both at low (e.g., Nikolov et al. 2018) and high spectral resolution (e.g., Redfield et al. 2008; Snellen et al. 2008; Jensen et al. 2012; Wyttenbach et al. 2015, 2017; Khalafinejad et al. 2017; Casasayas-Barris et al. 2018). (But see Gibson et al. 2019 for a retraction of the detection of potassium.) These detections motivate a series of theoretical and phenomenological questions concerning the sodium lines, which we address in the current study.

1.1 Are sodium lines in non-local thermodynamic equilibrium?

Local thermodynamic equilibrium (LTE) involves a set of four assumptions: the radiation field is described by a Planck or blackbody function, the velocity field follows a Maxwellian distribution, the distribution of neutrals versus electrons and ions is a solution of the Saha equation, and the energy levels of the ground versus excited states of the atom follow a Boltzmann distribution. One of the few studies to examine non-LTE (NLTE) effects associated with the sodium lines in exoplanetary atmospheres is Fortney et al. (2003), who considered non-Saha distributions of the electron density caused by photoionization.

Here, we take a different approach: we study the NLTE effect associated with non-Boltzmann distributions of the ground and excited states of the sodium atom. Physically, this departure from LTE arises from the diminished role of collisions. To maintain a Boltzmann distribution, collisions between the sodium atoms and the atoms or molecules of the bulk gas they are embedded in (e.g., atomic or molecular hydrogen) need to be efficient. There are two aspects to this NLTE effect. The first is that the sodium doublet, also known as the Fraunhofer D1 and D2 lines, is a pair of resonant lines, where the difference in transit radii between line center and wing spans about 20 scale heights (e.g., Heng 2016), corresponding to about 9 orders of magnitude in pressure. It implies that the importance of collisions varies between the peak and wings of each sodium line—they are more important in the line wings, which probe deeper into the atmosphere, and less influential at the line peak.

Figure 1 shows that the ratio of number densities of the first excited state of the sodium atom to its ground state may depart from the LTE value at the order-of-magnitude level, depending on the temperature and pressure. The departure from LTE is non-negligible for pressures of 1 $\mu$ bar and 1 mbar. The LTE and NLTE ratios of number densities agree only at a pressure of 1 bar, as shown in the top panel of Figure 1. Since transmission spectra are expected to probe pressures much less than 1 bar, this NLTE effect is worthy of investigation. The departure from NLTE produces a correction factor, which reduces to unity in the limit of LTE, to the integrated strength of the sodium lines (as we will elucidate in §2). When the sodium lines are in NLTE, the correction factor exceeds unity and produces weaker opacities. Figure 2 shows that the LTE model tends to over-estimate the overall strength of the sodium doublet in the transmission spectrum, regardless of whether clouds are present.

When applied to low-resolution Hubble Space Telescope (HST) transmission spectra (Sing et al., 2016), we find this NLTE effect produces minor corrections to the retrieved atmospheric properties, except for HD 189733b. However, the interpretation of a high-resolution HARPS transmission spectrum of the hot Jupiter WASP-49b (Wyttenbach et al., 2017) is biased towards lower temperatures if LTE is assumed, as we will demonstrate in §3.

1.2 Are sub-Lorentzian line wings of sodium needed to interpret spectra?

A classical approach to modeling the spectral lines of atoms or molecules is to assume that the shapes of these lines are described by the Voigt profile, which is a mathematical convolution of the Lorentz and Doppler profiles (e.g., Chapter 6.5 of Draine 2011, Chapter 5.2 of Heng 2017). For the sodium line observed at $\sim 1000$ K, the damping parameter is low enough ( $\sim 10^{-3}$ ) that the Voigt profile consists of a Doppler core and Lorentzian wings, if pressure broadening is ignored (Heng et al., 2015). However, the studies of Burrows et al. (2000) and Allard et al. (2012) suggest that a more realistic description of the sodium line shape consists of sub-Lorentzian wings, even though these studies do not agree on the quantitative details of how to compute these line wings.

From a phenomenological viewpoint, we wish to ask if sub-Lorentzian wings are a necessary ingredient for fitting transmission spectra of sodium lines, both at low and high resolutions? In the current study, we adopt an agnostic and data-driven approach, which is to describe the line wings of sodium by a dimensionless broadening parameter ( $f_{\rm broad}$ ). A classical Voigt profile with no pressure broadening present has $f_{\rm broad}=1$ . If pressure broadening is present, we have $f_{\rm broad}>1$ . ** In high-resolution transmission spectra, a resolved sodium line with $f_{\rm broad}>1$ may also encode information associated with the Rossiter-McLaughlin effect and large-scale atmospheric winds (e.g., Louden & Wheatley 2015)**. Any value of $f_{\rm broad}<1$ may be interpreted as the presence of sub-Lorentzian wings. Figure 2 shows that as the broadening parameter increases, the sodium lines become stronger in the wings (while the line-center strengths remain invariant). One expects a degeneracy between the broadening parameter and the mixing ratio of sodium, because the strength of the line wings may be negated by assuming lower abundances of sodium (Figure 2).

From performing retrievals on low-resolution transmission spectra, we find that sub-Lorentzian wings are not necessary to explain the current state of the data (§3). This is a phenomenological, rather than theoretical, statement.

1.3 What are the limitations of interpreting low- versus high-resolution spectra?

There is a debate in the exoplanet literature concerning how to leverage the value of low- versus high-resolution data off each other to maximise the outcome of retrievals. For example, Brogi et al. (2017) studied a joint emission spectrum of HD 209458b with data from the HST Wide Field Camera 3 (WFC3), the Spitzer Space Telescope Infrared Array Camera (IRAC) and the VLT CRIRES spectrograph. These authors demonstrated that the posterior distributions of molecular abundances become narrower at the order-of-magnitude level when the high-resolution data are considered in the retrieval analysis of the low-resolution spectrum. However, as the computations were expensive they were unable to explore the entire parameter space and instead used low-resolution retrievals to guide the exploration. The conclusions of Brogi et al. (2017) on whether the dayside of HD 209458b is isothermal, as well as the retrieved range of metallicity values for the atmosphere, appear to be dependent on their choice of parametrization of the temperature-pressure profile.

As another example, Pino et al. (2018) combined low-resolution data from several HST instruments with high-resolution transmission spectra from HARPS to study the hot Jupiter HD 189733b, which includes the detection of the sodium doublet by Wyttenbach et al. (2015). Instead of performing a retrieval, they fixed the temperature-pressure profile to that reported by Wyttenbach et al. (2015), which was itself derived from an isothermal treatment. By paying attention to the line spread function of the different instruments, Pino et al. (2018) were able to fit the combined low- and high-resolution spectrum with a solar-composition model, i.e., the water, sodium and potassium abundances were held fixed. Pino et al. (2018) explicitly acknowledged that they did not explore the normalization degeneracy (Benneke & Seager, 2012; Griffith, 2014; Heng & Kitzmann, 2017).

In the current study, we ask a different question: what are the limitations encountered when separately interpreting low- versus high-resolution transmission spectra of the sodium doublet? This question is non-trivial, because high-resolution, ground-based spectroscopy is a “double differential” technique, where the absolute empirical normalization of the transmission spectrum is lost during data calibration and reduction (Wyttenbach et al., 2015, 2017; Pino et al., 2018). In other words, only the relative, and not the absolute, transit depths are measured. Quantitatively, it is unclear how the enhanced spectral resolution (relative to, e.g., HST transmission spectra), diminished wavelength range and lack of an empirical normalization play off or negate one another. Space-based HST transmission spectra cover a wider wavelength range and do possess an absolute empirical normalization, but may not encode enough information (due to the low spectral resolution) to completely break the normalization degeneracy (Fisher & Heng, 2018), which states that the theoretical normalization of the transmission spectrum is uncertain and leads to order-of-magnitude uncertainties in the chemical abundances (Benneke & Seager, 2012; Griffith, 2014; Heng & Kitzmann, 2017).

To this end, we apply our retrieval technique to the HST Space Telescope Imaging Spectrograph (STIS) transmission spectra curated by Sing et al. (2016), which includes HAT-P-1b, HAT-P-12b, HD 189733b, WASP-6b, WASP-17b and WASP-39b. The sodium doublet is not resolved in these STIS spectra, but is rather detected as a single, unresolved spike of two blended lines in absorption. As a single case study of high-resolution, ground-based spectra, we analyze the HARPS transmission spectrum of the hot Jupiter WASP-49b measured by Wyttenbach et al. (2017). The sodium doublet is resolved in this HARPS spectrum. Besides the deep transits associated with the sodium doublet (about 2% relative to the continuum), the spectral signature of the Rossiter-McLaughlin effect was not detected (Wyttenbach et al., 2017), which makes WASP-49b an ideal high-resolution case study.

Previously, Fisher & Heng (2018) demonstrated that HST WFC3 near-infrared transmission spectra encode enough information to partially break the normalization degeneracy. Specifically, the caveat is that the lower limits to the posterior distributions correspond to the maximum value assumed for the prior distribution of the reference pressure. The HST STIS transmission spectra of Sing et al. (2016) cover the optical/visible range of wavelengths, where it has previously been established (e.g., Lecavelier des Etangs et al. 2008; Fisher & Heng 2018) that if the optical/visible spectral slope may be uniquely attributed to molecular hydrogen alone then the normalization degeneracy is completely broken. However, the expectation is that if the optical/visible spectral slope is due to clouds or hazes, then the normalizaton degeneracy cannot be broken (Heng, 2016). In the current study, we interpret the STIS transmission spectra by assuming that absorption by sodium atoms, Rayleigh scattering associated with small cloud particles and gray absorption by large cloud particles are present. We find that the normalization degeneracy cannot be broken in the context of such a model interpretation using optical/visible data alone.

For the high-resolution spectra, we find that the theoretical normalization cannot be uniquely determined, but given the lack of an empirical normalization this has little to no effect on the posterior distributions of the retrieved atmospheric properties.

1.4 Is the degree of cloudiness of an exoplanetary atmosphere correlated with any exoplanet property?

Previously, Stevenson (2016) and Heng (2016) used WFC3 transmission spectra dominated by water and STIS observations of the sodium doublet, respectively, to study the degree of cloudiness of exoplanetary atmospheres and if it is correlated with any exoplanet property. Collectively, these studies find that hotter objects (Heng, 2016; Stevenson, 2016) with higher gravity (Stevenson, 2016) are more likely to have cloudfree atmospheres.

One of the objectives of the present study is to improve upon and revisit the study of Heng (2016), who designed a diagnostic to interpret the sodium and potassium lines observed in transmission spectra at low resolution, building on the theory developed in Heng et al. (2015). It is based on the reasoning that, since these are strong resonant lines, the line peaks are unaffected by the opacities of clouds or hazes, while the line wings are significantly affected, thus allowing for the peak-to-wing distance to be used as a diagnostic for the degree of cloudiness. The diagnostic was constructed to apply at two points: one at the line peak, one in the line wing. This two-point approach was previously applied to the high-resolution HARPS spectrum of WASP-49b (Wyttenbach et al., 2017).

A more accurate approach is to perform a retrieval that fits for the partially resolved shape of the sodium lines. From the fit, one can compute the distance between the line peak and wing, $\Delta R_{\rm fit}$ . One can then remove cloud opacity from the model and compute the corresponding distance in a cloudfree atmosphere, $\Delta R_{\rm fit,CF}$ . By taking the ratio of these quantities, we obtain the cloudiness index (Heng, 2016),

[TABLE]

such that cloudfree and cloudy atmospheres have $C=1$ and $C>1$ , respectively. The reason to use the best-fit model, rather than the actual data points, in the denominator of $C$ is to avoid situations where an imperfect fit to the measured transmission spectrum produces a best-fit model line peak that sits below the data point, which produces spurious values of $C<1$ . We further investigate if estimates of $C$ are affected by the assumption of LTE.

1.5 Layout of paper

In Section 2, we lay out our methodology, including the theory behind the two-level atom treatment, the cross section of the sodium lines, models of the transmission spectrum and treatments of the high-resolution data format. In Section 3, we study a suite of mock retrievals and perform retrieval analyses on the HST STIS and HARPS transmission spectra. In Section 4, we discuss the implications of our results, including a scrutiny of the information content of high-resolution spectra with respect to detecting NLTE effects. Table 1 contains the physical constants associated with the sodium doublet. Table 2 states the physical units and prior ranges of the retrieved parameters. Table 3 summarizes the values of the retrieved atmospheric properties for the sample of 7 exoplanets considered in the current study.

2 Theory & Methodology

2.1 Two-level sodium atom

The two-level atom is not a novel concept and has previously been described in monographs (e.g., Chapter 17.1 of Draine 2011). Here, we recast it in a form that is suitable for computation and application to sodium atoms.

The two-level atom has a ground state and an excited state. The number density of atoms in the excited state ( $n_{2}$ ) is described by the following evolution equation (e.g., equation [17.5] of Draine 2011),

[TABLE]

where $t$ denotes the time, $n_{1}$ is the number density of atoms in the ground state, $n_{e}$ is the number density of electrons, $n_{\rm total}$ is the total number density (which is the number density of hydrogen molecules in a hydrogen-dominated atmosphere), $A_{21}$ is the Einstein A-coefficient (with units of s*-1*) and $N_{\gamma}$ is the photon occupation number. The statistical weights are $g_{j}=2j^{2}$ , where $j=1,2$ . The rate coefficent for collisional de-excitation is given by equation (2.27) of Draine (2011),

[TABLE]

where $T_{4}\equiv T/10^{4}$ K and $T$ is the temperature. As explained by Draine (2011), the collision strength $\Omega_{21}$ is approximately independent of temperature for $T<10^{4}$ K and typically has values between 1 and 10. To maximise the NLTE effect, we approximate $\Omega_{21}\approx 1$ , which produces the highest possible value of the critical density.

Physically, besides the collisional excitation ( $n_{\rm total}n_{1}C_{12}$ ) and de-excitation ( $n_{\rm total}n_{2}C_{21}$ ) of the excited level as well as de-excitation by spontaneous emission ( $n_{2}A_{21}$ ), there is also photoabsorption ( $n_{1}N_{\gamma}g_{2}A_{21}/g_{1}$ ) and stimulated emission ( $n_{2}N_{\gamma}A_{21}$ ), as described in Chapter 17.1 of Draine (2011). However, we show in the bottom panel of Figure 1 that the photon occupation number,

[TABLE]

with $h$ being Planck’s constant, $k_{\rm B}$ being Boltzmann’s constant and $\nu$ being the frequency, is much less than unity, implying that the terms associated with $N_{\gamma}$ in equation (2) may be dropped. Note that the preceding expression for $N_{\gamma}$ assumes a Planckian radiation field.

Demanding a steady state ( $\partial n_{2}/\partial t=0$ ) yields

[TABLE]

where $E_{12}$ is the energy difference between the ground and excited levels and $f_{e}\equiv n_{e}/n_{\rm total}$ is the ionization fraction. For $f_{e}\ll 1$ , we obtain the standard NLTE expression for a two-level atom,

[TABLE]

To arrive at equation (6), we invoke the principle of detailed balance (e.g., Section 3.5 of Draine 2011),

[TABLE]

The critical density is defined as

[TABLE]

When $n_{\rm total}\gg n_{\rm crit}$ (i.e., collisions are important), we recover the Boltzmann distribution and the states obey LTE. At a pressure of 1 bar, $n_{2}/n_{1}$ is well described by its LTE limit for temperatures between $10^{3}$ – $10^{4}$ K (Figure 1). However, for lower pressures, the departures from LTE may become significant at higher temperatures.

For completeness, Appendix A describes an analytical expression for the electron density if all of the electrons are being sourced by collisional ionization of the sodium atom.

2.2 Cross section of sodium lines

The absorption cross section of the sodium atom is the product of the integrated line strength ( $S$ ) and the line shape ( $\Phi$ ),

[TABLE]

There are three different ways to write $S$ , depending on if the line shape is written in per wavenumber, per frequency or per wavelength units. We use per frequency units, such that $S$ has units of cm2 s*-1* and $\Phi$ has units of s.

2.2.1 Integrated line strength

The integrated line strength is (Penner, 1952)

[TABLE]

Following Appendix 2 of Goody & Yung (1989), a factor of $1/n_{\rm total}$ has been inserted to give $S$ the correct physical units, such that $\sigma_{\rm Na}$ has units of cm2. The energy density per unit frequency is given by

[TABLE]

where $c$ is the speed of light and

[TABLE]

is the Planck function in per frequency units. The Einstein A-coefficient is related to the oscillator strength, $f_{12}$ , by (e.g., equation [6.20] of Draine 2011)

[TABLE]

where $e$ is the elementary charge and $m_{e}$ is the mass of the electron.

It follows that the integrated line strength becomes

[TABLE]

There are several noteworthy aspects of the preceding expression. First, since $N_{\gamma}\ll 1$ and $n_{1}\approx n_{\rm total}$ , we obtain

[TABLE]

which is our final NLTE expression for the two-level atom. Second, in the LTE limit we obtain $S\approx\pi e^{2}f_{12}/m_{e}c$ , which is the same as equation (6.25) of Draine (2011). Third, since equation (10) is the same starting point for deriving the standard expression for the line strength in HITRAN / HITEMP, it implies that equation (15) is equivalent to equation (A5) of Rothman et al. (1998).

2.2.2 Line shape

The Voigt profile is standard knowledge in monographs (e.g., Chapter 6.5 of Draine 2011, Chapter 5.2 of Heng 2017). However, there is more than one way to define the Lorentz and Doppler widths. Here, we concisely restate the formalism to be explicit about the conventions adopted in the current study.

The Lorentz profile is

[TABLE]

where $\nu_{0}$ is the frequency at line center and $\Gamma_{\rm L}=f_{\rm broad}A_{21}/4\pi$ is the half-width at half-maximum (HWHM) of the Lorentz profile. This convention differs from, for example, Draine (2011) who defines it as the full-width at half-maximum. We use the dimensionless broadening parameter $f_{\rm broad}$ , which expresses the strength of the Lorentzian wings in terms of the natural width, as an agnostic way of including all sources of broadening. In particular, since the theory of pressure broadening is incomplete, the use of $f_{\rm broad}>1$ allows us to report retrieved broadening values in terms of the natural width, which will motivate future studies of the sodium line shape. It also allows us to isolate the pressure dependence of the cross section due to the NLTE term associated with the total density. The use of $f_{\rm broad}<1$ allows us to quantify the extent to which the line wings are sub-Lorentzian. We do not truncate the sodium line wings.

The Doppler profile is

[TABLE]

The thermal speed is $v_{\rm th}=\sqrt{2k_{\rm B}T/m}$ , while the Doppler shift associated with the thermal speed is $\sigma_{\rm th}=\nu_{0}v_{\rm th}/c$ . $\Gamma_{\rm D}=\sigma_{\rm th}\sqrt{\ln{2}}$ is the HWHM of the Doppler profile.

By convolving the Lorentz and Doppler profiles, the Voigt profile obtains,

[TABLE]

where the Voigt H-function is

[TABLE]

and $x\equiv(\nu-\nu_{0})/\sigma_{\rm th}$ . The damping parameter is defined as $a_{0}\equiv\Gamma_{\rm L}/\sigma_{\rm th}$ .

In theory, the Voigt H-function is evaluated using (Roston & Obaid, 2005; Zaghloul, 2007)

[TABLE]

where $\texttt{erfc}(a_{0})$ is the complementary error function (e.g., Chapter 10 of Arfken & Weber 1995) with $a_{0}$ as its argument.

In practice, the Voigt profile is evaluated using (Schreier, 1992)

[TABLE]

where $w(z)$ is the Faddeeva function (equation [7.1.3]111The Faddeeva function appears without being specifically attributed to Faddeeva in Abramowitz & Stegun (1970). of Abramowitz & Stegun 1970) evaluated at (Schreier, 1992)

[TABLE]

This representation of the Voigt function allows for faster and simpler computation of the line shape.

2.3 Transmission spectrum

2.3.1 LTE

At a given wavelength or frequency, the expression for the transit radius corresponding to an isothermal transit chord (Lecavelier des Etangs et al., 2008; de Wit & Seager, 2013; Bétrémieux & Swain, 2017; Heng & Kitzmann, 2017; Jordán & Espinoza, 2018) is

[TABLE]

where $H=k_{\rm B}T/mg$ is the isothermal pressure scale height, $m$ is the mean molecular mass, $g$ is the gravity, $\gamma=0.57721$ is the Euler-Mascheroni constant and $\kappa$ is the opacity (cross section per unit mass). The reference pressure, $P_{0}$ , corresponds to a reference transit radius, $R_{0}$ .

For our model atmospheres, the opacity is

[TABLE]

where $X_{\rm Na}$ is the volume mixing ratio of sodium and $\kappa_{\rm cloud}$ is a gray/constant cloud opacity. The assumption is the wavelength range probed is narrow enough that it is insensitive to changes in cloud opacity caused by small particles, or that the particles are large (compared to the wavelength).

Denoting the atomic mass unit by $m_{\rm amu}$ , the mean molecular mass is given by

[TABLE]

where $m_{\rm Na}=23m_{\rm amu}$ is the mass of the sodium atom. The mixing ratio of molecular hydrogen is determined by demanding that all mixing ratios sum to unity,

[TABLE]

where we have assumed that the helium mixing ratio follows cosmic abundance ( $X_{\rm He}=0.1X_{\rm H_{2}}$ ). When $X_{\rm Na}\ll 1$ , we have $X_{\rm H_{2}}\approx 0.91$ . A continuum of possibilities, from hydrogen-dominated to sodium-dominated atmospheres, is allowed.

2.3.2 NLTE

The NLTE expression describing the number densities of the excited and ground states of sodium explicitly contains the pressure $P$ , because $n_{\rm total}=P/k_{\rm B}T$ if we assume the ideal gas law. This means that we need to explicitly account for the variation of pressure across the sodium line profile. To keep the problem tractable, we assume hydrostatic equilibrium, $P=P_{0}\exp{[(R_{0}-R)/H]}$ . One of the assumptions involved in deriving equation (23) is worth emphasizing. Early in the derivation, one assumes that, in writing the optical depth as an integral involving the cross section and number density, one may bring the cross section out of the integral, i.e., assume that it is a constant with respect to the integration variable. Since the integral is carried out over the transit chord, this is equivalent to assuming that the cross section is constant across the chord at a given wavelength. It does not preclude the fact that the pressure probed by the transit chord may be rather different at different wavelengths. In that sense, equation (23) is only isobaric at a given wavelength and may be used to describe spectral lines where the pressure probed across the line (across different wavelengths) may vary markedly.

The opacity is given by

[TABLE]

where we have defined

[TABLE]

We assume that the particles associated with the gray cloud opacity track the behavior of the sodium atoms with regards to departures from LTE.

Substituting the expression for the opacity into equation (23) yields

[TABLE]

where we have defined

[TABLE]

This way of writing equation (29) means that the two effects of $P_{0}$ have been separated: its regular appearance in $R_{\rm LTE}$ versus its appearance in the NLTE term, which is related to the varying importance of collisions across the sodium line.

By defining

[TABLE]

equation (29) may be written as a quadratic equation,

[TABLE]

where $A\equiv n_{\rm crit}k_{\rm B}T/P_{0}$ and $x_{0}\equiv\exp{[(R_{\rm LTE}-R_{0})/H]}$ . By solving this quadratic equation, we obtain a closed-form solution for the transit radius,

[TABLE]

It is apparent that the dependence of the transit radius on $A$ , and hence the $P_{0}$ term responsible for the varying importance of collisions across the sodium line, is weak.

2.3.3 Values used for the reference transit radius and reference pressure

Following Fisher & Heng (2018), we estimate the reference transit radius using

[TABLE]

In Fisher & Heng (2018), $\bar{R}$ is the average transit radius in the WFC3 bandpass. Here, we take $\bar{R}$ to be the HST whitelight radius as stated in Sing et al. (2016).

Based on WASP-17b, Fisher & Heng (2018) estimated that the WFC3 and STIS bandpasses probe $\sim 10$ mbar. This implies that $P_{0}\sim 10$ bar. For the low-resolution HST STIS retrievals, we therefore set a tight, log-uniform prior of $1\leq P_{0}\leq 100$ bar.

For the high-resolution spectrum of WASP-49b, we match the spectral continuum to the HST whitelight radius of $1.198~{}R_{\rm J}$ (Wyttenbach et al., 2017), where $R_{\rm J}$ is the radius of Jupiter. This requires setting $R_{0}=0.752R_{\rm J}$ for $P_{0}=10$ bar.

2.4 Nested-sampling retrievals and Bayesian evidence

Our exploration of the parameter space of models is performed using a nested-sampling framework (Feroz & Hobson, 2008; Feroz et al., 2009, 2013). We use the open-source PyMultiNest package (Buchner et al., 2014), which was previously implemented in Fisher & Heng (2018). Nested sampling allows us to compute the Bayesian evidence and Bayes factor, which allows us to identify the best model given the quality of the data (Trotta, 2008). It also allows us to identify the family of models consistent with the data (e.g., Fisher & Heng 2018).

2.5 Equivalence with Heng (2016)

Let the transit radius at line center and wing be $R_{\rm c}$ and $R_{\rm w}$ , respectively. Let also $\Delta R\equiv R_{\rm c}-R_{\rm w}$ . It follows from equation (23) that

[TABLE]

where $\Phi_{\rm c}$ and $\Phi_{\rm w}$ are the Voigt profiles at line center and wing, respectively. If we approximate the Voigt profile as consisting of a Doppler core and Lorentzian wings, then we recover equation (10) of Heng (2016),

[TABLE]

Thus, the pointwise approach of Heng (2016), as designed for analyzing low-resolution spectra, is a limiting case of the continuous approach given by equation (23) that we are applying in the current study to high-resolution spectra.

2.6 Data format for un-normalized high-resolution spectra

Let the transit depth be $D=(R/R_{\star})^{2}$ , where $R_{\star}$ is the stellar radius. Ground-based, high-resolution data do not measure the transit depth, but rather the differential transit depth, which is usually denoted by ${\tilde{\Re}}$ and lacks an absolute empirical normalization. We wish to understand the relationship between $D$ and ${\tilde{\Re}}$ .

Generally, we have $D\ll 1$ . Thus, $1-D$ is a number that is almost unity across wavelength. All of the peaks of $D$ are converted into troughs for $1-D$ . Let the minimum value of $D$ , which corresponds to the continuum, be $D_{\rm min}$ . The quantity ${\tilde{\Re}}=(1-D)/(1-D_{\rm min})$ shifts the minimum value of the continuum to unity. The quantity that is typically used for measured spectra is ${\tilde{\Re}}-1$ (Wyttenbach et al., 2015, 2017), which shifts the continuum to zero.

3 Results

3.1 Does the normalization degeneracy affect our interpretation of the sodium doublet?

Before analyzing measured transmission spectra of the sodium doublet, we study a suite of mock retrievals to elucidate the degeneracies between the various model parameters. We also wish to understand the effect of the normalization degeneracy on spectra with a narrow wavelength range centered on the sodium doublet. To generate mock spectra, we assume 100 wavelength bins between 5880 and 5905 Å. Each data point is assumed to have an uncertainty of 0.005 Å on the transit radius. These numbers are guided by the measured HARPS transmission spectrum of WASP-49b (Wyttenbach et al., 2017). We use the measured gravity of WASP-49b of 689 cm2 s*-1* (Wyttenbach et al., 2017). For illustration, we assume $T=4500$ K, $X_{\rm Na}=10^{-5}$ , $f_{\rm broad}=100$ and $\kappa_{\rm cloud}=100$ cm2 g*-1* in our NLTE model of a mock transmission spectrum.

We study two types of retrievals of mock spectra. Un-normalized spectra are presented in the ${\tilde{\Re}}-1$ format previously described in §2.6. Normalized spectra are presented in transit radii ( $R$ ), where we match the continuum of the mock spectrum to the measured white-light radius of $1.198~{}R_{\rm J}$ for WASP-49b (Wyttenbach et al., 2017). For the choice of $P_{0}=10$ bar, this requires setting $R_{0}=0.752~{}R_{\rm J}$ .

The top panel of Figure 3 shows pairs of retrievals where we fix $P_{0}=10$ bar and $R_{0}=0.752~{}R_{\rm J}$ versus if we fix only $R_{0}=0.752~{}R_{\rm J}$ and allow $P_{0}$ to be part of the retrieval. Several striking trends are revealed even by examining the un-normalized spectrum (top left panel). First, the behavior of the temperature and broadening parameter are distinct enough from the other parameters that they may be accurately retrieved, regardless of whether $P_{0}$ is known. It demonstrates that our parametrisation of line broadening using $f_{\rm broad}$ effectively isolates it from other effects, including temperature variations.

Second, there is value added to the retrieval when the transmission spectrum is normalized (top right panel). The posterior distribution of sodium abundance spans about 5 orders of magnitude (for the full width at half-maximum or FWHM) for the un-normalized spectrum, but shrinks to about 2–3 orders of magnitude simply by matching the spectral continuum to the HST white-light radius. Our ignorance of $P_{0}$ for the normalized spectrum results in the cloud opacity being loosely unconstrained, but this does not affect our ability to accurately retrieve the temperature and broadening parameter.

Third, for the normalized spectrum (top right panel of Figure 3), both the degeneracy with cloudiness (i.e., any increase in the sodium abundance may be compensated by an increase in the gray cloud opacity) and the normalization degeneracy (i.e., any increase in the sodium abundance may be compensated by a decrease in $P_{0}$ ) are clearly present in the posterior distributions. These degeneracies clearly drive the somewhat large widths of the posterior distributions of both the gray cloud opacity and sodium abundance. When the spectrum is un-normalized (top left panel of Figure 3), the normalization degeneracy vanishes. Essentially, over the width of the posterior distributions of both $X_{\rm Na}$ and $\kappa_{\rm cloud}$ as determined by the cloudiness degeneracy, the reference pressure is unconstrained. Figure 4 shows the corresponding posterior distributions for the cloudiness index, which is robust to whether $P_{0}$ is treated as a fitting parameter and whether the retrieval is performed on un-normalized or normalized spectra.

In a real situation, the value of $R_{0}$ is unknown. As demonstrated by Fisher & Heng (2018), it suffices to make a reasonable guess for $R_{0}$ and then retrieve for $P_{0}$ . The bottom panels of Figure 3 elucidates the effect of good ( $R_{0}=0.85R_{\rm J}$ ), not-so-good ( $R_{0}=0.95R_{\rm J}$ ) and bad ( $R_{0}=1.05R_{\rm J}$ ) guesses. For the normalized spectrum, a bad guess results in a posterior for the cloud opacity that bumps up against its prior value. More insidiously, it also gives a posterior distribution of the sodium abundance that is narrower than the one for the good guess. A not-so-good guess also produces misleadingly narrow posterior distributions of the sodium abundance and cloud opacity.

The bottom left panel of Figure 3 shows that these concerns about guessing $R_{0}$ are mostly irrelevant for the un-normalized spectrum, because the posterior distributions of parameters are nearly identical for all quantities. The exception is the temperature, which is somewhat different for the bad guess. In other words, the retrieval of un-normalized spectra is unaffected by the normalization degeneracy, but at the expense of more uncertain posterior distributions compared to normalized spectra.

3.2 Can we distinguish between LTE and NLTE scenarios?

We now use mock retrievals to study a different question, which is whether a suite of nested-sampling retrievals are capable of distinguishing between LTE versus NLTE and cloudfree versus cloudy models? Using the same cloudy NLTE mock spectrum, we perform four retrievals using the cloudfree LTE model, cloudy LTE model, cloudfree NLTE model and cloudy NLTE model. Figures 5 and 6 show the outcomes of these retrievals. The first major outcome is that the shape of the sodium doublet, in un-normalized transmission spectra, does not encode enough information to distinguish between cloudfree and cloudy models based on the Bayesian evidence. It has the implication that the temperature is under-estimated if one uses a cloudfree model to interpret a cloudy transit chord.

A comparison of the Bayesian evidence also does not allow us to distinguish between NLTE and LTE models, which in our formulation have exactly the same number of parameters. A robust outcome of the retrievals is that** the value of the broadening parameter is accurately retrieved regardless of whether cloudy NLTE or LTE models are used**. However, when a LTE model is used to perform retrieval on the NLTE mock spectrum, the retrieved temperature is substantially under-estimated. In the example shown, we retrieve a temperature of about 2600 K compared to the true input value of 4500 K. This suite of mock retrieval suggest that retrievals performed on un-normalized HARPS transmission spectra of the sodium doublet are unable to discern between NLTE and LTE models.

In both the cloudy LTE and NLTE models, the width of the posterior distribution of the sodium abundance is determined by the degeneracy between the gray cloud opacity and the sodium abundance.

3.3 Retrieval analysis of HARPS data of WASP-49b

**The mock retrievals of Figure 5 assume idealized conditions, meaning that sources of contamination due to stellar activity, telluric lines, imperfect co-addition of multiple exposures, etc, are not considered. In analyzing the WASP-49b HARPS transmission spectrum of Wyttenbach et al. (2017), we assume that these issues have been addressed by the authors. **

Wyttenbach et al. (2017) recorded the peak of each sodium line within a narrow band with a width of 0.4 Å (see their Section 4.2). As discussed in Section 2.6, the reference transit depth, $D_{\rm min}$ , is theoretically taken as the minimum transit depth across the spectrum. However, due to the fluctuations in the continuum points in the data, Wyttenbach et al. (2017) set this reference value to the average measured continuum within two reference bands: blue (5874.94–5886.94 Å) and red (5898.94–5910.94 Å). We apply the same procedure to our models: within the same pair of reference bands, we compute the average value of the continuum. We then take $D_{\rm min}$ as this value when calculating ${\tilde{\Re}}$ , before shifting it to ${\tilde{\Re}}=1$ . Our fit to the data using this procedure is shown in Figure 7.

The top row of Figure 8 shows our retrieval analysis of the measured high-resolution, un-normalized transmission spectrum of WASP-49b (Wyttenbach et al., 2017). In our initial retrievals (not shown), we struggled to fit the deep line peaks with our model. This motivated us to extend our prior for $f_{\rm broad}$ to lower values, compared with those stated in Table 2. For this data, we use the range $[10^{-10},10^{4}]$ . As expected, we are unable to distinguish between the LTE versus NLTE interpretation from the computed Bayesian evidence (448.8 for LTE versus 448.4 for NLTE). On physical grounds alone, we favor the NLTE interpretation. The temperature is $7209^{+1763}_{-1892}$ K for the LTE interpretation versus $8415^{+1020}_{-1526}$ K for the NLTE interpretation, consistent with the lesson learned from our mock retrievals that the LTE model tends to predict a lower temperature. These retrieved temperatures are discrepant from the $2950^{+400}_{-500}$ K value reported by Wyttenbach et al. (2017), but we note that Wyttenbach et al. (2017) fitted Gaussians rather than Voigt profiles to the measured sodium doublets. The use of Gaussians is equivalent to fitting with a Doppler profile in the absence of Lorentzian wings. Our retrieved temperatures remain consistent with the expectation that radiative cooling by collisionally excited atomic hydrogen thermostats the temperature to $\sim 10^{4}$ K (Murray-Clay et al., 2009). The retrieved sodium volume mixing ratios are loosely constrained and somewhat insensitive to whether the LTE or NLTE interpretation is assumed ( $\log{X_{\rm Na}}=-5.64^{+2.55}_{-2.90}$ versus $-4.44^{+2.02}_{-2.51}$ ), as is the retrieved cloud opacity ( $\log{\kappa_{\rm cloud}}=-0.18^{+3.14}_{-2.77}$ versus $-0.90^{+2.82}_{-2.06}$ ). The broadening parameter is essentially unconstrained, but takes on $f_{\rm broad}<1$ values as the narrow wavelength range precludes full coverage of the sodium line wings.

The continuum of the measured transmission spectrum appears to possess a scatter with a half-width of about 0.005. To test the robustness of our retrieval results, we artificially increased the two line peaks by 0.005 each and re-ran our retrievals. The retrieval outcomes are shown in the bottom row of Figure 8. For the LTE retrieval, we obtained a reduced temperature of $T=6756^{+2048}_{-2030}$ K and a similar sodium abundance ( $\log{X_{\rm Na}}=-6.22^{+3.03}_{-2.86}$ ). For the NLTE retrieval, the temperature becomes $T=7925^{+1408}_{-1869}$ K and the sodium abundance becomes $\log{X_{\rm Na}}=-4.76^{+2.32}_{-3.01}$ . This pair of tests suggests that the retrieval outcomes are somewhat robust to the uncertainty in the continuum.

3.4 Retrieval analysis of low-resolution transmission spectra

We examine a sample of 6 low-resolution transmission spectra: HAT-P-1b, HAT-P-12b, HD 189733b, WASP-6b, WASP-17b and WASP-39b (Sing et al., 2016). Given the sparseness of the data, we assume a simple model that includes the sodium doublet and a spectral continuum sourced by clouds with both small and large particles of unspecified composition. Large cloud particles are represented by a constant opacity ( $\kappa_{\rm cloud}$ ). Following Sing et al. (2016), small cloud particles are assumed to contribute Rayleigh scattering that has a behavior, across wavelength, which is identical to that of molecular hydrogen, but with a magnitude that is offset by some constant, dimensionless factor $f_{\rm cloud}\geq 1$ . The cross section for Rayleigh scattering by molecular hydrogen is (Sneep & Ubachs, 2005),

[TABLE]

where $n_{\rm ref}=2.68678\times 10^{19}$ cm*-3* and the real part of the index of refraction is (Cox, 2000)

[TABLE]

Collectively, the opacity function is

[TABLE]

where $m_{\rm H_{2}}=2m_{\rm amu}$ is the mass of the hydrogen molecule.

The computed model spectrum is degraded to the measured one by taking the average value of the transit radius in each wavelength bin. Figures 9 and 10 show the outcomes of these 12 retrievals. The reference pressure is either unconstrained or only loosely constrained and its posterior distribution is largely prior-dominated. This is unsurprising, because the pressure probed by the transit chord is

[TABLE]

If the atmosphere is cloudfree ( $f_{\rm cloud}=1$ ) and hydrogen-dominated ( $X_{\rm H_{2}}\approx 0.91$ ), then the pressure can be uniquely determined. The transit-chord pressure and reference pressure are related by hydrostatic equilibrium. However, if the value of $f_{\rm cloud}$ is a priori unknown, then $P$ and hence $P_{0}$ cannot be uniquely determined.

At low spectral resolutions, the trend of decreasing $f_{\rm broad}$ producing strong lines (relative to the continuum) is negated by a higher abundance of sodium, and it is this degeneracy that dominates the width of the posterior distributions of $X_{\rm Na}$ . In all 12 retrievals, the temperature is tightly constrained as it functions like an independent “stretch mode” in the retrieval (Fisher & Heng 2018; see also Figure 2).

HD 189733b is singled out for discussion due to its historical importance in the literature (Pont et al., 2008, 2013; Sing et al., 2011). It is the only object among the sample of 6 objects that has a clearly cloudy STIS transit chord ( $f_{\rm cloud}\gg 1$ ), but it is important to note that the range of $f_{\rm cloud}$ values retrieved is prior-dominated by our choice of $1\leq P_{0}\leq 100$ bar. The temperature corresponding to the STIS transit chord is about 2800 K (assuming LTE) or about 4800 K (for NLTE), but these values are not inconsistent with those reported by Wyttenbach et al. (2015) and Heng et al. (2015). In all 12 retrievals, the values of $\kappa_{\rm cloud}$ retrieved are small, consistent with the continuum being non-gray.

Figure 11 shows the empirical trends (or lack thereof) between the various retrieved quantities. For the purpose of computing cloudiness index, the cloudfree scenario is obtained by setting $f_{\rm cloud}=1$ and $\kappa_{\rm cloud}=0$ . No clear trend between the cloudiness index ( $C$ ) and the equilibrium temperature, which is a proxy for the insolation or stellar irradiation, is seen, contrary to the tentative trend reported by Heng (2016). The retrieved $C$ values range from being consistent with unity to larger than unity, in broad agreement with the study of Barstow et al. (2017). NLTE interpretations yield a higher degree of cloudiness, because of the diminished strength of the sodium lines relative to LTE.

The uncertainties on the retrieved sodium abundances are large, but the general trend is that $X_{\rm Na}$ is consistent with being solar ( $\sim 10^{-6}$ ) or subsolar. The uncertainties on the retrieved broadening parameters are also large, but the retrieved $f_{\rm broad}$ values are broadly consistent with unity, suggesting that sub-Lorentzian wings are not needed to fit the data. Generally, $T/T_{\rm eq}\sim 1$ , where $T$ is the retrieved temperature of the transit chord, with the exception of HD 189733b where $T/T_{\rm eq}>2$ .

4 Discussion

Our retrieval outcomes in Figure 8 clearly demonstrate that ground-based, high-resolution spectra of the sodium doublet alone do not encode enough information to infer the pressure level being probed by the lines—not even at the order-of-magnitude level. This finding is consistent with the lessons learned from our mock retrievals in Figure 3, where we demonstrated that the normalization degeneracy and the lack of an empirical normalization prevents $P_{0}$ from being meaningfully constrained. The finding that ground-based high-resolution spectroscopy is incapable of accurately retrieving the pressures probed appears to contradict the work of Pino et al. (2018) for HD 189733b, but we note that this work did not address the normalization degeneracy and instead assumed a fixed value for $R_{0}$ that corresponds to $P_{0}=10$ bar. It is this assumption of fixing the reference radius that allows pressure levels to be inferred. The broader implication of this finding is that one cannot easily infer which part of the atmosphere one is probing (i.e., thermosphere, exosphere) if the sodium lines are analyzed in isolation. Future work should elucidate if ground-based, high-resolution spectra alone (with their lack of an empirical normalization), even across an extended wavelength range, encode enough information (e.g., in the sodium line wings) to provide precise information on the pressures probed by the transit chord.

In the current study, we have chosen to focus only on the optical part of the spectrum that contains the sodium lines in order to elucidate the limitations associated with such a restricted analysis. Future work should combine spectra from the optical and infrared, and elucidate if the near-infrared water lines probed by WFC3 also require a NLTE treatment as has been presented in the current study.

While the retrieved outcomes in the current study are consistent with Voigt profiles for the sodium line shapes, future work should elucidate a theory of $f_{\rm broad}=f_{\rm broad}(P)$ that reconciles and unifies the work of Burrows et al. (2000) and Allard et al. (2012).

We are grateful to Chris Hirata for guidance on the solution to the Saha equation for the sodium atom and to Jens Hoeijmakers and David Ehrenreich for useful discussions on the relative transit depth. We acknowledge financial support from the Swiss National Science Foundation, the European Research Council (via a Consolidator Grant to KH; grant number 771620), the PlanetS National Center of Competence in Research (NCCR), the Center for Space and Habitability (CSH) and the Swiss-based MERAC Foundation.

Appendix A Electron density via solution of Saha equation

For completeness, we provide the expression for $n_{e}$ if the electrons were entirely sourced by the collisional ionization of the sodium atom,

[TABLE]

To render the problem tractable, we retain the assumption that $n_{e}$ is a solution of the Saha equation, i.e., the ionization states are in LTE, which reads (e.g., Section 3.4 of Draine 2011)

[TABLE]

where $E_{\rm ion}$ is the ionization energy. $g_{{\rm Na}^{+}}$ and $g_{\rm Na}$ are the quantum degeneracies associated with Na+ and Na, respectively. If we assume that the electrons are solely provided by the ionization of the sodium atom, then we have $n_{e}=n_{{\rm Na}^{+}}$ and

[TABLE]

where $X_{\rm Na}\equiv n_{\rm Na}/n_{\rm total}$ is the volume mixing ratio of sodium. Previously, the preceding expression was used to describe the number density of electrons in protoplanetary disks for the purpose of studying the magneto-rotational instability, where it was assumed that all of the electrons were sourced by potassium, which has an ionization energy of $E_{\rm ion}=4.3407$ eV (Balbus & Hawley, 2000). It was also used to study Ohmic dissipation in hot Jupiters (Perna et al., 2010). In the current situation, we assume that the electrons are sourced only by sodium, which has $E_{\rm ion}=5.1391$ eV. We have $g_{{\rm Na}^{+}}/g_{\rm Na}=1/2$ because the singly-charged ion possesses a closed shell of electrons, but the neutral atom has a single valence electron in the $3s$ orbital.

Bibliography54

The reference list from the paper itself. Each links out to its DOI / PubMed record.

1Abramowitz & Stegun (1970) Abramowitz, M., & Stegun, I.A. 1970, Handbook of Mathematical Functions, 9th printing (New York: Dover)
2Allard et al. (2012) Allard, N.F., Kielkopf, J.F., Spiegelman, F., Tinetti, G., & Beaulieu, J.P. 2012, A&A, 543, A 159
3Arfken & Weber (1995) Arfken, G.B., & Weber, H.J. 1995, Mathematical Methods for Physicists, Fourth Edition (Academic Press)
4Balbus & Hawley (2000) Balbus, S.A., & Hawley, J.F. 2000, Space Science Reviews, 92, 39
5Barstow et al. (2017) Barstow, J.K., Aigrain, S., Irwin, P.G.J., & Sing, D.K. 2017, Ap J, 834, 50
6Benneke & Seager (2012) Benneke, B., & Seager, S. 2012, Ap J, 753, 100
7Bétrémieux & Swain (2017) Bétrémieux, Y., & Swain, M.R. 2017, MNRAS, 467, 2834
8Brogi et al. (2017) Brogi, M., Line, M., Bean, J., Désert, J.-M., & Schwarz, H. 2017, Ap JL, 839, L 2