The Structure of Titan’s N2 and CH4 Coronae

In this study, we analyze the structures of Titan’s N2 and CH4 coronae using a large data set acquired by the Ion Neutral Mass Spectrometer (INMS) instrument on board Cassini. The N2 and CH4 densities measured from the exobase up to 2000 km imply a mean exobase temperature of 146 K and 143 K, respectively, which is lower than the mean upper atmospheric temperature by 4 and 7 K. This indicates that on average, Titan possesses a subthermal rather than suprathermal corona. A careful examination reveals that the variability in corona structure is not very likely to be solar driven. Within the framework of the collisionless kinetic model, we investigate how the CH4 energy distribution near the exobase could be constrained if strong CH4 escape occurs on Titan. Several functional forms for the CH4 energy distribution are attempted, assuming two representative CH4 escape rates of s−1 and s−1. We find that the double Maxwellian and power-law distributions can reproduce the shape of the CH4 corona structure as well as the imposed CH4 escape rate. In both cases, the escape rate is contributed by a suprathermal CH4 population on the high-energy tail, with a number fraction below 5% and a characteristic energy of 0.1–0.6 eV per suprathermal CH4 molecule. The coexistence of the subthermal CH4 corona revealed by the INMS data and substantial CH4 escape suggested by some previous works could be reconciled by a significant departure in the exobase CH4 energy distribution from ideal Maxwellian that enhances escape and causes a noticeable redistribution of the corona structure.


Introduction
A corona is the uppermost part of an atmosphere that extends outward from the exobase and gradually fades into interplanetary space. The corona structure on terrestrial planets, including Venus, Mars, and Titan, has been extensively studied over the past few decades. On Venus, the atomic H corona was first discovered by the Mariner 5 ultraviolet (UV) photometer (Barth et al. 1967). Using these data, Anderson (1976) was able to fit the atomic H corona with a thermal population at 275K and a suprathermal population at 1020K. The presence of the suprathermal population was later confirmed by Chaufray et al. (2012) based on Lyα observations made with the Spectroscopy for Investigation of Characteristics of the Atmosphere of Venus (SPICAV) instrument on board Venus Express (VEx). Such a suprathermal population was thought to be contributed by momentum transfer from solar wind (SW) energetic protons to atomic H in the Venusian upper atmosphere (Hodges 1999). A hot O corona on Venus was also revealed by 130.4nm observations made with the UV spectrometer on board Pioneer Venus Orbiter (PVO) (Nagy et al. 1981) and Venera 11 and 12 (Bertaux et al. 1981), and was thought to be produced by dissociative recombination of + O 2 , the dominant constituent of the Venusian ionosphere (Johnson et al. 2008). However, the presence of this hot O corona was not confirmed by more recent VEx SPICAV observations (Lichtenegger et al. 2009).
On Mars, the atomic H and O coronae were observed by the UV spectrometer on board Mariner 6, 7, and 9 (Anderson 1974), the Spectroscopy for Investigation of Characteristics of the Atmosphere of Mars (SPICAM) on board Mars Express (MEx) (Chaufray et al. 2008), the Alice far-UV imaging spectrograph on board Rosetta (Feldman et al. 2011), and the Imaging Ultraviolet Spectrograph (IUVS) on board Mars Atmosphere and Volatile Evolution (MAVEN) Chaufray et al. 2015;Deighan et al. 2015). It has been suggested that both coronae are composed of a thermal population and a suprathermal population, although the existence of a hot atomic H population was questioned by Rosetta Alice observations (Feldman et al. 2011). The most recent measurements made by MAVEN IUVS revealed a transition from thermal O corona  to suprathermal O corona at around an altitude of 600 km . The mechanisms responsible for the production of hot atoms on Mars are analogous to those occurring on Venus (e.g., Lee et al. 2015). For instance, the MAVEN IUVS observations of the O corona are best reproduced by including a suprathermal population with a characteristic energy of 1.1 eV, consistent with a scenario driven by dissociative recombination of ionospheric + O 2 .
On Titan, reliable measurements of corona N 2 and CH 4 made by the Ion Neutral Mass Spectrometer (INMS) on board Cassiniextend up to an altitude of 2000km (De La Haye et al. 2007a), whereas measurements of corona H 2 extend up to at least 7000km (Cui et al. , 2011. Existing studies revealed no remarkable deviation in exobase temperature of H 2 from the mean temperature in the upper atmosphere (Cui et al. , 2011, implying that a hot H 2 population does not exist. For N 2 and CH 4 , the analysis of the TA, TB, and T5 INMS data by De La Haye et al. (2007a) revealed the presence of suprathermal coronae for both constituents with a typical temperature enhancement of ∼20-70K over the mean state below the exobase. The above distinction calls for certain heating mechanisms to be operative for N 2 and CH 4 , but not for H 2 near Titan's exobase. De La Haye et al. (2007b) later evaluated the role of exothermic ion-neutral chemistry, but found it inadequate to produce the suprathermal N 2 and CH 4 populations. In addition, the atomic H corona on Titan was sampled by the Cassini Hydrogen Deuterium Absorption Cell (HDAC) instrument, and existing studies using the HDAC data did not support the presence of a hot atomic H population (Hedelt et al. 2010).
The present work is aimed at a thorough investigation of Titan's N 2 and CH 4 coronae, which becomes timely with the accumulation of significantly more INMS data since the publication of the early work by De La Haye et al. (2007a) and the continuous improvement in INMS data analysis over the past 10 years (Magee et al. 2009;Cui et al. 2009Cui et al. , 2012Teolis et al. 2015). The paper is organized as follows. We start with a brief description of INMS data reduction in Section 2. This is followed by an analysis of neutral temperature in Titan's upper atmosphere in Section 3, as an extension of numerous existing studies (e.g., Müller-Wodarg et al. 2008;Westlake et al. 2011;Snowden et al. 2013). We then derive in Section 4 the characteristics of Titan's corona using the traditional kinetic model of Chamberlain (1963) to evaluate the possible existence of suprathermal N 2 and CH 4 populations. The structure of the CH 4 corona is of special importance because it is highly relevant to CH 4 escape from Titan, an issue that has caught extensive research interest over the past ten years and is still under debate Strobel 2009;Tucker & Johnson 2009;Bell et al. 2010Bell et al. , 2011Cui et al. 2012;Schaufelberger et al. 2012). The link between these two issues is the topic of Section 5. Finally, we discuss and provide concluding remarks in Section 6.

Retrieval of the N 2 and CH 4 Abundances from the INMS Data
The analysis presented in this study relies exclusively on the Cassini INMS data in the closed source neutral (CSN) mode available at the Planetary Plasma Interactions (PPI) node of the NASA Planetary Data System (PDS) public archives (http://ppi.pds.nasa.gov). Thirty Cassini flybys with Titan, from TA on 26 October 2004 to T107 on 10 December 2014, are included, with continuous INMS measurements of Titan's upper atmosphere and corona from the closest approach (CA) up to about 2000 km above the surface. Only the inbound data of these flybys are used because surface chemistry occurring on the walls of the INMS antechamber may seriously contaminate the outbound data of various constituents (Vuitton et al. 2008). Meanwhile, all measurements made with ram angle, i.e., the angle between the INMS aperture normal direction and the spacecraft velocity direction, greater than 75°are excluded since the configuration of the instrument in these cases was not favored for collecting ambient neutral particles. Note that for a few recent flybys (T86, T91, T92, T95, T98, and T100), only the region of Titan's atmosphere below the exobase was sampled by the INMS. Accordingly, these flybys are not considered here even though they are useful for characterizing the neutral upper atmosphere.
The approach used to retrieve the N 2 and CH 4 abundances from the raw INMS counts is based on the well-established algorithm described in Cui et al. (2012). The calculation of the ram enhancement factor is improved by taking into account the realistic neutral gas transmission through the instrument, including flow pathways that have not been recognized in previous works (Teolis et al. 2015). The N 2 and CH 4 sensitivity values used to convert count rates into number densities are lower than the values in Cui et al. (2009) andCui et al. (2012) by a factor of 1.55 to account for INMS detector fatigue (Teolis et al. 2015). With the above updates, the empirical enhancement factor of around 3 that has previously been imposed to force equality between the INMS densities and the Cassini Ultraviolet Imaging Spectrograph (UVIS) densities is no longer required (Koskinen et al. 2011).
The choice of the maximum altitude of 2000km quoted above is a crucial constraint for the purpose of this study, ensuring that the residual N 2 and CH 4 gases within the INMS antechamber do not substantially affect the density determination in Titan's corona (Cui et al. 2009). For instance, the mean CH 4 density at 2000km above the surface is1 10 5 cm −3 (see Section 4). Using a representative ram enhancement factor of 27 and an appropriate sensitivity value of´-4.4 10 4 cm 3 s −1 , we estimate the count rate contributed by ambient CH 4 to be 1200s −1 at the same altitude. The typical count rate over the inbound portion of a Titan flyby is 10s −1 for the residual CH 4 gas within the INMS antechamber (see Figure36 of Cui et al. 2009), ensuring that the contribution from this residual gas to the observed count rate is no more than 1%.
An example of INMS data analysis is provided in Figure 1 for the T83 flyby on 22 May 2012, with CA at a latitude (LAT) of 73°N, a longitude (LON) of 128°E, a solar zenith angle (SZA) of 71°, and a Saturn local time (SLT) of 13.7hr. Panel (a) shows the N 2 densities obtained from the primary counter (C1) of mass channel 28, (  . The actual N 2 densities are therefore determined by combining the count rates in all the three cases according to the weighting functions shown in panel (b). The choice of these weighting functions smoothly connects different regimes of counter saturation (Cui et al. 2012) and avoids artificial discontinuity in the N 2 density data that has been easily produced in previous data analysis works (e.g., Müller-Wodarg et al. 2008;Cui et al. 2009). The N 2 density profile shown in panel (a) appears to be asymmetric about CA, with a clear enhancement over the outbound portion of the flyby when > t 400 s CA . This asymmetry is indicative of wall surface chemistry (Vuitton et al. 2008) and not an actual horizontal variation in the ambient atmosphere . The CH 4 densities (not shown in the figure) are derived exclusively from measurements made with the C1 counter of mass channel 16, with slight saturation near CA corrected following the nonlinear algorithm outlined in Cui et al. (2012) and also with a clear signature of wall surface chemistry over the outbound portion of the flyby.

Characterization of the Mean Thermal State in Titan's Upper Atmosphere
To characterize Titan's corona structure, it is essential to derive the mean thermal state below the exobase a priori. Here for each flyby, the N 2 density data are used to calculate the isothermal neutral temperature in Titan's upper atmosphere, denoted as T n , assuming hydrostatic equilibrium. An example is presented in Figure 1(c) for inbound T83, with the best-fit hydrostatic equilibrium model indicated by the black solid line with » T 148 n K. The best-fit values of T n for all flybys are listed in Table 1 for reference.
The analysis of the N 2 density profiles reveals a highly variable neutral temperature ranging from 113 to 197K, with a globally averaged temperature of 150K. This is in good agreement with previous results also obtained from the INMS data, but including a smaller set of Titan flybys (e.g., Cui et al. 2009;Westlake et al. 2011;Snowden et al. 2013). The lowest temperature is encountered in T43 and the highest temperature in T104, the latter being a dayside Titan flyby that has not been analyzed in any existing works. For some of the common flybys analyzed in this and previous studies, we note that the neutral temperature differs by a non-negligible amount. This is mainly due to the update in INMS data reduction, especially in terms of the ram enhancement calculation, and also due to the different altitude coverage used for model fitting. However, the crucial results in terms of both the globally averaged value and the flyby-to-flyby variability remain the same.
Despite the large variability in neutral temperature, no reliable horizontal and diurnal trends could be identified, as implied by Figure 2, where we show the variations in neutral temperature with LAT, LON, SLT, and SZA, all referred to CA, and with the 10.7cm solar radio flux. The color style of the circles drawn in the figure indicates the absolute value, whereas the size of the circles is chosen to be proportional to the 10.7cm solar radio index at 1au in solar flux unit (SFU, i.e., 10 −19 erg cm −2 s −1 Hz −1 , see Table 1). The diurnal variation is our particular interest. Taking into account the extended nature of Titan's atmosphere (e.g., Müller-Wodarg et al. 2000), we define the dayside subsample as including all flybys with SZA at CA below 110°, and the nightside subsample as SZA at CA beyond 10°. Analysis of the INMS data then reveals an average dayside temperature that is slightly lower than the average nightside temperature by about 4K, indicating that the variability in the thermal structure of Titan's upper atmosphere is unlikely to be driven by solar energy deposition (Snowden et al. 2013). A multivariate regression analysis is also performed to parameterize the correlation between the derived neutral temperature and various geophysical parameters. The correlation is weak in most cases, and the strongest correlations occurs between T n and F10.7, with a correlation coefficient of 0.35.

Characterization of Titan's N 2 and CH 4 Coronae
In this section, the density profiles of N 2 and CH 4 above Titan's exobase are modeled following the simple kinetic approach developed by Chamberlain (1963), with the corona density, N, of either constituent given by where z is the altitude, N e is the exobase density, ζ is the partition function that depends on the Jeans parameter, λ, defined as l = , with G being the gravitational constant, k b the Boltzmann constant, M T Titan's mass, R T Titan's radius, m the molecular mass, and T e the weighting functions used to derive the actual N 2 densities free of counter saturation, for inbound T83. Panel (c): INMS N 2 and CH 4 density profiles in Titan's upper atmosphere (with altitude indicated on the left) and corona (with altitude indicated on the right), also for inbound T83, superimposed by the best-fit hydrostatic equilibrium model for N 2 below the exobase with a neutral temperature of = T 148 n K (solid black) and the best-fit Chamberlain models in the corona with an exobase N 2 temperature of 177K (solid red) and an exobase CH 4 temperature of 159K (dashed red). exobase temperature. l e indicates the Jeans parameter at Titan's exobase, placed at 1500km throughout this study (Westlake et al. 2011). The functional form of z l ( ) can be found in Chamberlain (1963) and is not repeated here. In practice, N e is obtained directly from the INMS measurements, whereas T e is treated as the only free parameter to be constrained by data-model comparison. For illustrative purpose, the best-fit corona models for T83 are given by the red lines in Figure 1(c), with an inferred exobase temperature of 177K for N 2 and 159K for CH 4 , respectively. The same approach has also been used in previous analyses of Titan's H 2 corona (Cui et al. , 2011, atomic H corona (Hedelt et al. 2010), as well as N 2 and CH 4 coronae (De La Haye et al. 2007a).
Applying the Chamberlain model to the data acquired over all available flybys, we find a best-fit exobase N 2 temperature from 105 to 177K with a globally averaged value of 146K, and a best-fit exobase CH 4 temperature from 100 to 173K with a globally averaged value of 143K. These mean temperatures are slightly lower than the mean temperature of Titan's upper atmosphere by 4K and 7K, respectively (see Table 1). Interestingly, these observations imply that on average, Titan possesses a subthermal rather than suprathermal corona, in contrast to the early finding of De La Haye et al. (2007a). In that work, both the inbound and outbound portions of TA and T5 were analyzed, as well as the outbound portion of TB, and suprathermal N 2 and CH 4 coronae were identified for all cases, except for inbound T5. Since the publication of this early result, it has been recognized that the outbound densities of many constituents within Titan's upper atmosphere and corona are problematic because of contamination by wall surface chemistry (Vuitton et al. 2008). Therefore we speculate that to some extent, the early conclusion drawn by De La Haye et al. (2007a) reflects instrumental bias rather than physical reality. For the two cases analyzed in both works, i.e., inbound TA and T5, we confirm the early result of suprathermal N 2 and CH 4 coronae for the former and subthermal N 2 and CH 4 coronae for the latter, although the exact values of temperature change across the exobase do not perfectly agree between the two works. For inbound TA, we find a temperature enhancement of 5K for N 2 and 24K for CH 4 , to be compared with the respective values of 24K and 21K from De La Haye et al. (2007a). For inbound T5, these authors reported a common temperature decrement of 13K for both constituents, to be compared with our values of 7K for N 2 and 21K for CH 4 . The discrepancies quoted above are primarily a result of improved INMS data reduction over the years (see Section 2 for details). , obtained from Chamberlain model fits to the N 2 and CH 4 data above the exobase, and the corresponding 10.7cm solar radio index at 1au. The globally averaged values are provided in the last row.
Our analysis reveals a highly variable corona structure on Titan for both N 2 and CH 4 , showing evidence for suprathermal corona with temperature enhancement across the exobase for some flybys, but subthermal corona with temperature decrement for the others. The Chamberlain fits to two representative CH 4 density profiles above Titan's exobase are shown in Figure 3, one for inbound T28 in panel (a), a typical example with temperature enhancement over the exobase, and the other for inbound T29 in panel (b), a typical example with temperature decrement. In each panel, an alternative Chamberlain model using the neutral temperature below the exobase but the same exobase CH 4 density is superimposed with the blue solid line for comparison. A similar change in temperature across the exobase is also seen in the N 2 density data for each flyby (see Table 1). For N 2 , the largest temperature enhancement is 31K encountered in T107 and the largest temperature decrement is 42K in T104, both dayside flybys. For CH 4 , the largest temperature enhancement of 29K and decrement of 44K are encountered in the same flybys. A careful examination reveals that in all the 30 flybys, only 5 (T39, T57, T59, T83, and T107) show clear evidence of a suprathermal N 2 corona at more than 3σ significance level, and 5 more (TA, T28, T40, T58, and T87) at 2-3σ. In contrast, 12 out of the remaining 20 flybys (T5, T23, T25, T27, T29, T36, T42, T48,  T56, T61, T65, and T104) show clear evidence of a subthermal N 2 corona at more than 3σ, and 3 more (T32, T50, and T71) at 2-3σ. Similar statistics applies to CH 4 as well. Only 6 (TA, T28, T39, T57, T59, and T107) show clear evidence of a suprathermal CH 4 corona at more than 3σ, and 2 more (T83 and T87) at 2-3σ. Twelve of the remaining 22 flybys (T5, T23,  T25, T27, T29, T36, T42, T48, T61, T65, T71, and T104) show clear evidence of a subthermal CH 4 corona at more than 3σ, and two more (T40, T50) at 2-3σ. In general, the above observations suggest that the presence of a suprathermal corona is sporadic rather than stable on Titan, whereas the presence of a subthermal corona is a more common feature.
We further show in Figure 3(c) the relation between the neutral temperature derived in Section 3 and the exobase N 2 (red) and CH 4 (blue) temperatures derived in this section. The figure reveals that the corona structures for both constituents are highly correlated, with the temperature change for N 2 across the exobase coincident with the temperature change for CH 4 in most cases. The only exceptions are T18 and T40. The INMS T18 data imply a temperature decrement of 1K for N 2 but a temperature enhancement of 4K for CH 4 , both insignificant by no more than 2σ. The situation is reversed for T40, giving a temperature enhancement of 7K for N 2 but a temperature decrement of 6K for CH 4 , with a significance level of 3.5σ and 2σ, respectively. The general coincidence in the N 2 and CH 4 corona structures calls for a common heating  Chamberlain fits to two representative CH 4 density profiles above Titan's exobase, one for inbound T28, a typical example with temperature enhancement over the exobase, and the other for inbound T29, a typical example with temperature decrement. In each panel, the best-fit model is given by the solid black line, while the model using the neutral temperature below the exobase is given by the solid blue line. Panel (c): correlation between the neutral temperature derived in Section 3 and the exobase N 2 (red) and CH 4 (blue) temperatures derived in Section 4, with the solid line in the middle representing no temperature change across the exobase. mechanism that is operative for both N 2 and CH 4 near Titan's exobase.
In Figure 4 we show the variations of the best-fit exobase N 2 and CH 4 temperatures with LAT, LON, SLT, and SZA, all referred to CA, and with the 10.7cm solar radio flux. The drawing convention is the same as in Figure    various geophysical parameters, as shown schematically in Figure 5. The best correlation occurs between exobase N 2 temperature and solar radio flux with a correlation coefficient of only 0.57. The above results suggest that the structural variations of Titan's N 2 and CH 4 coronae are unlikely to be solar driven. Further support could be gained by comparing the corona structures in Figures 3(a) and (b). The two flybys shown in the figure sample roughly the same nightside region of Titan's atmosphere over the northern hemisphere and separated by only one Titan day. This implies that the structure of Titan's corona may vary substantially over relatively short timescales, which is more likely caused by the highly variable plasma environment in the vicinity of Titan (e.g., Arridge et al. 2011a).
In reaching the above conclusions, we have implicitly assumed that the N 2 and CH 4 density profiles measured by the INMS, especially above Titan's exobase, properly reflect the vertical trends and that any horizontal variations could be reasonably ignored. Indeed, Cui et al. (2011) have shown that the observed variation in Titan's H 2 corona is more likely indicative of variation in the temporal domain rather than in the spatial domain. We expect a similar situation for Titan's N 2 and CH 4 coronae as well based on the physical argument that any horizontally spatial variation seen below the exobase tends to be smoothed out when propagating upward into the corona, but the temporal variation could be reserved as long as the timescale for the variation is longer than the typical dynamical timescale of 10 3 s. Meanwhile, our key conclusion that Titan possesses on average a subthermal rather than suprathermal corona is free of systematic uncertainties caused by horizontal variations.

Implications for CH 4 Escape from Titan
The structure of Titan's CH 4 corona is highly relevant to the issue of CH 4 escape from the satellite, since both depend on the CH 4 energy distribution function (EDF) near the exobase (De La Haye et al. 2007a). The former is obtained by integrating the exobase CH 4 EDF over the full momentum space occupied by realistic particle trajectories that cross the exobase with upward velocities, and the latter obtained by integrating the same EDF over a restricted momentum space with kinetic energy greater than the local escape energy. The Chamberlain model used in Section 4 implicitly assumes an ideal Maxwellian for the exobase CH 4 EDF, given by where E is the particle energy, and other parameters have been defined above. Such a functional form gives a Jeans escape rate of no more than 10 20 s −1 according to the range of exobase CH 4 temperatures quoted in Section 4. Adopting a non-Maxwellian EDF, the escape rate could be much higher.  Li et al. (2014) shows that their eddy diffusion (red line) crosses the CH 4 binary diffusion (black line) at around 950 km, giving a high homopause for this study as well.
Placing the homopause at 1000km, Bell et al. (2014) inferred a much lower CH 4 escape rate of1.2 10 25 s −1 . Second, Bell et al. (2010) proposed an alternative mechanism for CH 4 loss via aerosol trapping in Titan's upper atmosphere, from which they obtained an even lower CH 4 escape rate at the level of 10 21 s −1 . Third, the direct simulation Monte Carlo (DSMC) calculations made by Johnson (2009) andSchaufelberger et al. (2012) suggested that the INMS CH 4 measurements could be well reproduced with a negligible CH 4 escape rate that is compatible with the Jeans escape rate. Meanwhile, it has been argued that a torus composed of carbon-group ions would be produced if the CH 4 escape rate is as high as 10 27 s −1 , but such a torus has not been observed by the Cassini Plasma Spectrometer (CAPS) and Magnetosphere Imaging Instrument (MIMI) (Johnson 2009;Arridge et al. 2011b). Clearly, existing studies predict a wide range of CH 4 escape rates from negligible to several 10 27 s −1 , with the lower limit being compatible with the Chamberlain model presented in Section 4. We allow here other likely choices for the CH 4 escape rate, and examine their implications on the exobase CH 4 EDF based on the observed CH 4 corona structure in a globally averaged sense. Two representative cases are considered: one with a CH 4 escape rate of1.2 10 25 s −1 from Bell et al. (2014), and the other with an escape rate of2.2 10 27 s −1 that reflects the maximum possible value. The latter is determined from the INMS CH 4 data below the exobase following the approach of Cui et al. (2012).
In both cases, the imposed CH 4 escape rates are much higher than the Jeans value. This clearly calls for a departure in CH 4 EDF at the exobase from ideal Maxwellian, in that some of the CH 4 molecules at relatively low energy are accelerated and populate the high-energy tail. Such a characteristics is compatible and not intuitively in conflict with the existence of subthermal CH 4 corona on Titan. Specifically, if the departure in EDF is large enough, it would lead to a noticeable redistribution of the corona structure manifested as a density decrement in the lower part of the corona and an accompanying density enhancement in the upper part. This could be understood with the fact that more energetic particles ejected from the exobase tend to populate more distant regions of a planetary corona. Accordingly, the departure in EDF has the apparent effect of reducing the exobase temperature if the region of interest is constrained to the lower part of the corona, which is likely the case studied here since this is the same region where the INMS is able to gather sufficient count rates from ambient CH 4 and allows reliable determination of its density. An analogous modeling of the upper part of the corona alone well above 2000km likely predicts a much higher exobase CH 4 temperature instead, reflecting a transition from subthermal to suprathermal somewhere beyond the altitude range considered here. Unfortunately, validating such a speculation is not possible with the INMS data since the CH 4 abundances cannot be retrieved reliably at such high altitudes (see Section 2). For comparison, we note that such a transition in atomic O corona does occur on Mars at around 400km above the exobase (e.g., Chaufray et al. 2015;Deighan et al. 2015).
There are several ways that the departure in CH 4 EDF at the exobase could be invoked, which we discuss in turn below.

k-distribution
The κ-distribution is frequently used as a straightforward replacement for ideal Maxwellian for systems in stationary states out of thermal equilibrium (Livadiotis & McComas 2013). It was also used in the early INMS study of De La Haye et al. (2007a) to parameterize the CH 4 EDF at Titan's exobase with enhanced high-energy tail. Following Tucker et al. (2016), the κ-distribution is written as where Γ is the gamma function, κ is a nondimensional index reflecting the contribution from nonthermal particles, and k ( )k T 3 2 b represents the mean particle energy. k ( ) ( ) f E e reduces to ideal Maxwellian in the limit of k  +¥. With N e retrieved directly from the INMS data, κ and k T are treated as two free parameters to be constrained by comparing the corona density structure measured by the INMS and that calculated based on the formula outlined in De La Haye et al. (2007a) using Equation (3) as the exobase CH 4 EDF.
The best-fit model for CH 4 corona structure is shown by the solid black line in Figure 6(a), superimposed on the globally averaged INMS data. The uncertainties shown in the figure are evaluated following the Monte Carlo procedure described in Cui et al. (2012). The solution is obtained by searching for the minimum normalized c 2 within the two-dimensional parameter space, as demonstrated in Figure 6(b). The best-fit parameters are = k T 152 K and k = 30, corresponding to a globally averaged CH 4 escape rate of2.4 10 21 s −1 , which is higher than the Jeans value, but still falls short by several orders of magnitude to match the range of CH 4 escape rate predicted by existing fluid model calculations (e.g., Yelle et al. 2008;Strobel 2009;Cui et al. 2012;Bell et al. 2014). We further apply an additional constraint that the imposed CH 4 escape rates of1.2 10 25 s −1 and2.2 10 27 s −1 be reproduced. For the former, the allowed range of models falls on the solid red line in the figure inset, b1, indicating that κ should not exceed 8 over a wide range of k T up to 190K. The pair of model parameters that best matches the INMS CH 4 densities above the exobase is found to be = k T 153 K and k = 6.4. The corresponding model corona structure is indicated by the solid red line in Figure 6(a), providing an unsatisfactory fit to the INMS data. To reproduce the alternative CH 4 escape rate of 2.2 10 27 s −1 , κ should not exceed 3 over a range of k T up to 300K, as indicated in the figure inset, b2. The pair of model parameters that best matches the observed INMS CH 4 densities above the exobase is = k T 210 K and k = 2.2, with the model corona structure (given by the solid blue line in Figure 6(a)) providing an even poorer fit to the INMS data. A careful examination of a sequence of model runs clearly demonstrates that a high κ value is required to reproduce the shape of the corona structure, whereas a low κ value is required to reproduce a CH 4 escape rate much higher than the Jeans value. In practice, it is impossible to reconcile the two requirements with the κ-distribution, and a more flexible functional form of EDF that incorporates more free parameters has to be invoked instead.

Double Maxwellian Distribution
One of the simplest functional forms of EDF with such a flexibility is the double Maxwellian distribution that treats the CH 4 gas near Titan's exobase artificially as the combination of two populations, both under thermal equilibrium, but physically separated. The tenuous nature of the region near Titan's exboase renders the collisional coupling between CH 4 molecules weak, and the coexistence of two separate best-fit model assuming a κ-distribution for the exobase CH 4 EDF given by the solid black line. This model predicts a low CH 4 escape rate of2.4 10 21 s −1 . The best-fit models with CH 4 escape rates of1.2 10 25 s −1 and 2.2 10 27 s −1 are shown by the solid red and blue lines, respectively, both providing unsatisfactory fits to the data. Panel (b): distribution of normalized c 2 with model parameters, k T and κ (see text), for the case when no constraint on the CH 4 escape rate is imposed. The best-fit solution is indicated as a red star for clarity. In the figure insets, the solutions for two models with imposed CH 4 escape rates of1.2 10 25 s −1 (b1) and2.2 10 27 s −1 (b2) are indicated with similar drawing convention.
populations is therefore not entirely unrealistic provided that the timescale for extra heating is shorter than the mean collision timescale. The double Maxwellian distribution is formulated as where subscripts 1 and 2 denote the cold and hot populations, respectively, with + N N e1 e2 fixed by the INMS measurements made near the exobase. The inclusion of more free parameters (three as compared to two in Section 5.1) allows a rigorous determination of a best-fit solution that is able to reproduce both the measured shape of the CH 4 corona structure and the imposed CH 4 escape rates reasonably well.
The best-fit models are given by the solid lines in Figures 7(a) and ( K for the hot population. It is clear that in both cases, the CH 4 abundance in the corona is dominated by the cold population, at least over the altitude range probed by the INMS, whereas the CH 4 escape rate is almost exclusively contributed by the hot population. These facts reveal a possible scenario of (2-5)% of the CH 4 molecules near Titan's exobase being accelerated to a characteristic energy of eV by some unknown mechanism, presumably nonthermal. Although this energy is lower than the CH 4 escape energy of 0.37eV on Titan, the spread of CH 4 molecules around the characteristic energy allows a considerable portion of these hot molecules to gain sufficient energy and overcome Titan's gravitational potential.

Power-law Distribution
Finally, we consider a model for exobase CH 4 EDF with a power-law tail superimposed on an ideal Maxwellian above some threshold energy, E th . This EDF, also used by De La Haye et al. (2007a), is formulated as  where T p is a model parameter that characterizes a cold CH 4 population (to be distinguished from T e in the Chamberlain model, see Section 4), and A and E th are two additional parameters that characterize the fraction and mean energy of a hot CH 4 population on the tail. It is easily obtained from Equation (5)b that the fraction of hot CH 4 molecules is a -( ) A 2 3 3 2 , where a = ( ) E k T th b p , and the mean energy is E 3 th . The power index is fixed as 2.5 in this study for illustrative purposes, as done in De La Haye et al. (2007a) and Tucker et al. (2016). This choice ensures that the CH 4 escape flux and energy density be convergent. Strictly speaking, other values of the power index are physically allowed, but we do not take them into account partly because introducing the power index as an additional parameter does not provide a unique solution that reproduces both the shape of the corona and the imposed escape rates.
Again, N e is directly obtained from the INMS CH 4 data, and the other three parameters (T p , A, α or E th ) are constrained by data-model comparison, with the additional requirement that the globally averaged CH 4 escape rates of1.2 10 25 s −1 oŕ 2.2 10 27 s −1 be reproduced. The best-fit model parameters for an imposed CH 4 escape rate of1.2 10 25 s −1 are found to be = T 152 p K, =´-A 7.0 10 3 , and a = 18, as compared to = T 131 p K, A=1.4, and a = 19 for an imposed CH 4 escape rate of2.2 10 27 s −1 . A common threshold energy of = E 0.2 eV th is found for both escape rates. Similar to the double Maxwellian case, CH 4 escape, if present on Titan at a rate above 10 25 s −1 , should be exclusively contributed by suprathermal CH 4 molecules on the power-law tail, whereas less energetic CH 4 molecules are responsible for the corona structure sampled by the INMS. We note that the best-fit model parameters for the power-law distribution lead to a hot CH 4 fraction of no more than 2% at the exobase and a respective mean energy of 0.6eV, in broad agreement with the values predicted by the double Maxwellian model (see Section 5.2).

Discussions and Concluding Remarks
In this study, we analyze the N 2 and CH 4 distributions in Titan's upper atmosphere and corona using a large data set acquired by the INMS instrument on board Cassini. The density measurements made above Titan's exobase are described by the collisionless Chamberlain model, predicting a mean exobase temperature of 146K for N 2 and 143K for CH 4 . These values, when compared with the mean N 2 temperature of 150K in the upper atmosphere, reveal that on average, Titan possesses subthermal rather than suprathermal N 2 and CH 4 coronae, in contrast to the conclusion made by De La Haye et al. (2007a) using a small INMS data set and subject to instrumental uncertainties not recognized in early INMS analysis works (cf. Teolis et al. 2015, and references therein). A careful examination of the flyby-to-flyby variability reveals that the occurrence rate of suprathermal N 2 and CH 4 coronae is no higher than 20% at more than 3σ significance level, with the maximum temperature enhancement in N 2 and CH 4 across Titan's exobase found to be 31K and 29K, respectively, both occurring during the dayside flyby of T107. In contrast, the occurrence rate of subthermal N 2 and CH 4 coronae is more than twice as high, with the maximum temperature decrement found to be 42 and 44K for another dayside flyby of T104. The formation of suprathermal N 2 and CH 4 coronae on Titan reported by De La Haye et al. (2007a) appears to be sporadic, in analogy to the sporadic CH 4 escape from the same satellite claimed by Cui et al. (2012).
The structure of Titan's CH 4 corona is highly relevant to the issue of CH 4 escape from the satellite, since both depend on the CH 4 EDF near the exobase. In this study, we investigate such a distribution with two imposed scenarios of CH 4 escape, one with a rate of1.2 10 25 s −1 obtained by Bell et al. (2014), and the other with a rate of2.2 10 27 s −1 based on the updated INMS sample and the approach of Cui et al. (2012). The latter value characterizes the maximum level of CH 4 escape reported in previous works (e.g., Yelle et al. 2008;Strobel 2009). Several functional forms for the CH 4 EDF at Titan's exobase are attempted in this study, all characterized by an enhanced high-energy tail, including the κ-distribution, double Maxwellian distribution, and power-law distribution. These are presumably induced by nonthermal processes such as sputtering and exothermic chemistry in Titan's upper atmosphere (e.g., Cravens et al. 1997;Shematovich et al. 2003;Michael et al. 2005). The collisionless kinetic approach of De La Haye et al. (2007a) is then used to model the observed CH 4 corona structure and infer the optimal values of the free parameters for each case. We find that in a globally averaged sense, both the double Maxwellian and power-law distributions are able to reproduce the CH 4 corona structure and meanwhile match the imposed CH 4 escape rates. In both cases, the CH 4 escape rate is exclusively contributed by a suprathermal population on the high-energy tail, with a number fraction of no more than 5% and a characteristic energy of 0.1-0.6eV per particle.
The above conclusion is compatible and not intuitively in conflict with the existence of a subthermal CH 4 corona revealed by the INMS data, since an enhancement in high-energy tail of the EDF has to be accompanied by a depletion of less energetic particles. This likely leads to a noticeable redistribution of the corona structure in terms of a density decrement in the lower part of the corona and an accompanying density enhancement in the upper part. In consequence, Chamberlain modeling of Titan's lower corona, where reliable determination of the ambient CH 4 density is possible, naturally predicts a temperature decrement in CH 4 across the exobase. Further support of such a scenario is gained from the fact that of all the 15 flybys showing clear evidence of strong CH 4 escape from Titan (Cui et al. 2012), the presence of subthermal CH 4 corona is found three times more often than that of suprathermal CH 4 corona.
In most cases, the density at a given location within a planetary corona is primarily contributed by ballistic particles with relatively low energies ejected from the exobase, whereas the escape rate is entirely contributed by particles with sufficiently high energies to overcome the gravitational potential of the parent body. This implies that the density distribution within the corona is not a fair diagnostic of the exobase EDF, and some additional constraint has to be included in kinetic modeling. One of the simplest constraints could be imposed by requiring that the escape rate from diffusion model fitting be reproduced from a kinetic point of view by integrating the EDF over all energies exceeding the escape energy at the exobase. This makes the crucial difference in method between the present work and that of De La Haye et al. (2007a). We have shown in Section 5 that the above constraint, though simple, is efficient for a crude parameterization of the high-energy portion of the CH 4 EDF at Titan's exobase. In both De La Haye et al. (2007a) and this work, an abrupt transition from strongly collisional to collisionless is assumed at the exobase. In the more recent work of Tucker et al. (2016), the effect of finite collisions near the exobase is evaluated via the implementation of a DSMC model. The difference in CH 4 density from the exobase up to at least 2000km is insignificant between the DSMC model and the collisionless model (see Figure3 of Tucker et al. (2016)), and accordingly, we expect the results presented here to be valid.
Several representative functional forms of the CH 4 EDF at Titan's exobase, obtained with the imposed escape rate of 2.2 10 27 s −1 , are compared in Figure 8, from which we could evaluate their respective impacts on CH 4 escape as well as CH 4 redistribution in the corona. On the one hand, for both the ideal Maxwellian distribution, given by the solid black line, and the κ-distribution with k = 30, given by the dashed red line, the sharp decline in EDF above the escape energy implies a small fraction of energetic CH 4 molecules capable of leaving Titan's gravitational potential. On the other hand, for the κ-distribution with k = 2.2, given by the solid red line, and the double Maxwellian distribution, given by the solid blue line, the enhancement in EDF at the high-energy tail is substantial such that a relatively high CH 4 escape rate could be reproduced. In contrast, whether the observed corona density structure could be correctly modeled relies on the shape and magnitude of the EDF at energies that are well below the escape energy. Figure 8 indicates that for three of the models shown in the figure (the ideal Maxwellian distribution, the double Maxwellian distribution, and the κ-distribution with k = 30), the low-energy portions of the corresponding EDFs agree accurately. All these models reproduce the observed CH 4 corona structure below 2000km well, despite the substantial difference in EDF above around 0.1eV. This clearly demonstrates that the CH 4 corona structure sampled by the INMS is not appropriate for constraining the CH 4 escape rate. For comparison, the remaining model (κ-distribution with k = 2.2) is characterized by an EDF significantly depleted below 0.05eV, leading to an underestimate of CH 4 corona density below 1750km and an overestimate above, as revealed by Figure 6(a).
We ought to emphasize that the collisionless kinetic model used in this study is an empirical model and not able to provide a unique solution to the exobase CH 4 EDF. Two models presented in Section 5 with entirely different functional forms of exobase EDF (double Maxwellian and power law) appear to be acceptable, and we cannot exclude the possibility that other models work equally well. Despite this, the two acceptable models predict very similar values for the fraction and characteristic energy of the suprathermal CH 4 population. These values should have important implications on the nature Figure 8. Comparison of several representative functional forms of the CH 4 EDF at Titan's exobase, obtained with an imposed CH 4 escape rate of2.2 10 27 s −1 . These functional forms include the Maxwellian distribution in solid black and the κ-distribution in dashed red that reproduce the corona structure but not the escape rate, the κ-distribution in solid red that reproduces the escape rate but not the corona structure, as well as the double Maxwellian distribution in solid blue that reproduces both the corona structure and the escape rate. The vertical dashed line on the bottom axis marks the CH 4 escape energy at the exobase. of possible nonthermal mechanisms if strong CH 4 escape does occur on Titan, which we discuss below.
A number of nonthermal mechanisms were investigated in previous works (Johnson et al. 2008, and references therein). The DSMC modeling of Shematovich et al. (2003) predicted a total nitrogen escape rate of3.6 10 25 s −1 in the form of N 2 and N, due to sputtering by incident magnetospheric ions and pick-up ions. A similar result was obtained by Michael et al. (2005), also based on DSMC modeling. Multiplying this value by the CH 4 to N 2 cross section and density ratios near Titan's exobase gives a crude estimate of the CH 4 escape rate aś 2.6 10 25 s −1 due to sputtering, comparable to one of the two imposed values of the CH 4 escape rate in Section 5. We caution that the early DSMC results for Titan's atmospheric escape due to sputtering were based on limited pre-Cassini knowledge of the varying plasma environment of this satellite, and an updated analysis of sputter-induced escape is clearly deserved, taking advantage of the extensive measurements made over the past 10 years by several in situ plasma instruments on board Cassini(e.g., Arridge et al. 2011a). However, we do not expect a substantial difference in the mean energy of escaping N 2 and CH 4 molecules since the physical processes leading to atmospheric sputtering are relatively well known (Johnson 1994). For N 2 , Shematovich et al. (2003) reported a mean energy of 2eV for sputtering by incident magnetospheric ions and 7eV for sputtering by pick-up ions. The mean energy of escaping CH 4 molecules is likely higher (since CH 4 is lighter than N 2 ) and is therefore an order of magnitude too high compared to the characteristic energy of suprathermal CH 4 molecules obtained here.
Exothermic chemistry in Titan's upper atmosphere, initiated by either incident solar EUV/X-ray photons or precipitating charged particles from Saturn's magnetosphere, may also play an important role in driving atmospheric escape from Titan. This has been evaluated using either a simplified fluid model (Cravens et al. 1997) or a more sophisticated kinetic model (De La Haye et al. 2007b). The former predicted a CH 4 escape rate of2 10 25 s −1 , and the latter predicted an escape rate about a factor of 2 lower. Cravens et al. (1997)  4 . The exothermic energies released from these reactions are 7.9eV and 0.9eV, respectively, of which 0.5eV and 0.6eV are imparted to CH 4 . This is in fair agreement with the characteristic energy found in Sections 5.2 and 5.3. Detailed models of Titan's upper atmosphere are clearly required to ensure that the fraction of suprathermal CH 4 molecules inferred from realistic ion-neutral chemistry is also in agreement with the range we find.
Finally, the observed variations of Titan's corona do not reveal any systematic trend with SZA or solar radio index, essentially ruling out the possibility that these variations are driven by solar energy deposition. The solar-driven scenario has also been excluded for interpreting the thermal structure of Titan's upper atmosphere by various authors (Westlake et al. 2011;Snowden et al. 2013;Snowden & Yelle 2014). The drastic change in corona structure from T28 to T29, two consecutive flybys sampling roughly the same regions of Titan's atmosphere but separated by one Titan day (see Table 1), reveals a variability on relatively short timescales that could only be driven by Titan's highly variable plasma environment (Arridge et al. 2011a). Therefore we expect that charged particle precipitation from the ambient plasma likely plays an important role in determining the EDF near Titan's exobase and the corona structure as well, although the nature of the controlling mechanisms still remains to be determined by future studies.