# Methane seasonal cycle at Gale Crater on Mars consistent with regolith adsorption and diffusion

John Moores, Raina Gough, German Martinez, Pierre-Yves Meslin, Christina

Smith, Sushil Atreya, Paul Mahaffy, Claire Newman, Christopher Webster

► To cite this version:

John Moores, Raina Gough, German Martinez, Pierre-Yves Meslin, Christina Smith, et al.. Methane seasonal cycle at Gale Crater on Mars consistent with regolith adsorption and diffusion. Nature Geoscience, 2019, 12 (5), pp.321-325. 10.1038/s41561-019-0313-y. hal-04767526

HAL Id: hal-04767526

https://hal.science/hal-04767526v1

Submitted on 2 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
# Methane seasonal cycle at Gale Crater on Mars consistent with regolith adsorption and diffusion

John E. Moores<sup>1*</sup>, Raina V. Gough<sup>2</sup>, German M. Martinez<sup>3</sup>, Pierre-Yves Meslin<sup>4</sup>, Christina L. Smith<sup>1</sup>, Sushil K. Atreya<sup>3</sup>, Paul R. Mahaffy<sup>5</sup>, Claire E. Newman<sup>6</sup> and Christopher R. Webster<sup>7</sup>

A strong, repeatable seasonal cycle in the background methane mixing ratio has been observed at the Gale Crater landing site of the Mars Science Laboratory rover with the Tunable Laser Spectrometer of the Sample Analysis at Mars instrument. However, as of yet, no physical process has been proposed that can explain both the timing and amplitude of the observations. Here we show that a one-dimensional numerical model considering adsorption onto and diffusion through the regolith can reproduce the variation, including a phase lag, if the regolith is impregnated with methane from a prior plume or supplied from below by micro-seepage. Combining the model results with geological constraints, we estimate that the amount of microseepage at Gale is at most $3 \times 10^{-5}$ tonnes km$^{-2}$yr$^{-1}$. Gale's unique dynamical environment makes such seeps easier to detect in surface sampling measurements. Over most of the Martian surface, atmospheric mixing is stronger or atmospheric transport more effective, and we expect the amplitude of the seasonal cycle to be smaller for the same strength of seep.

Studying the seasonal cycles of volatile materials can be especially helpful to deduce the processes operating on a planetary body in the present day. These observations do not need to be global in nature to be useful. For instance, on Earth, the seasonal cycle of carbon uptake by forests from the atmosphere can be observed as a variation in observed CO$_2$ from atmospheric sampling conducted in Hawaii<fnref n="1" />. On Mars, the variation in total pressure observed at Gale Crater describes the annual cycle of CO$_2$ ice deposition at the poles and can also elucidate atmospheric transport² and secular changes in the Martian climate³. As such, the seasonal cycle in methane detected in the Tunable Laser Spectrometer of the Sample Analysis at Mars (SAM-TLS) enrichment runs<fnref n="1" />,⁵ is of particular interest, as methane has a relatively short atmospheric lifetime of several hundred years on Mars⁶ and candidate processes for its production have geological and astrobiological significance<fnref n="7" />.

Both the mean mixing ratio and amplitude of the methane seasonal cycle are small compared with past measurements⁸⁻¹² and to the 7 ppbv methane spikes previously reported at Gale Crater¹³,¹⁴. However, in a relative sense, the seasonal cycle is large with a variation of nearly a factor of three between the lowest and highest values observed over the Martian year. This large amplitude excludes many environmental processes that might be responsible for the variation⁴. For instance, while ultraviolet radiation can be roughly correlated with methane concentration, the seasonal variation in ultraviolet flux at the surface of Gale⁴,¹⁵ is only ~20% and no physical model of increasing methane release from the surface as a direct result¹⁶ of increased ultraviolet explains how this small change in the flux could bring about such a large relative change in methane concentration. Several other mechanisms for explaining the variation were proposed in previous work⁴ and subsequently ruled out, leaving no convincing mechanism for producing the seasonal cycle observed.

A potential solution is offered by processes that can amplify small seasonal variations through exponential dependencies, and that correlate with insolation. An example of such a process is adsorption

onto the regolith when combined with diffusion into and out of the regolith. The temperature of the regolith correlates well with the incoming solar flux¹⁷ and the rates of adsorption and desorption of methane onto mineral surfaces have been previously observed to have an Arrhenius dependency¹⁸. This is not the first time that methane adsorption on the Martian regolith has been examined to infer a seasonal cycle. Another study¹⁹ has used a global circulation model and laboratory-derived constants¹⁸ to determine the seasonal variation of methane across the planet due to adsorptive transfer to and from the regolith. At Gale's latitude, this seasonal variation was found to be less than a few percent¹⁹ consistent with the low adsorption coefficient measured for this gas, extrapolated from laboratory experiments¹⁸,¹⁹ to Mars ground temperatures.

However, this study¹⁹ did not explore a source of methane within the subsurface in their simulations of the seasonal cycle, nor thermophysical parameters different from laboratory experiments. In particular, microseepage is not unreasonable, nor unexpected<fnref n="7" />, especially in sedimentary environments like that of Gale²⁰. At these locations, it is hypothesized that organic carbon from early in Mars' history may be sequestered below the surface where it may then be thermally converted to methane and slowly released<fnref n="7" />. It is also possible that there are ancient methane clathrates below the annual temperature wave²¹,²² that are slowly degrading. Furthermore, serpentinization reactions may produce methane at depth in the current era<fnref n="7" />. Regardless of the origin, as this ancient methane rises towards the surface, it will need to diffuse through the regolith, adsorbing and desorbing onto mineral grains as it travels.

## Model results compared with SAM-TLS data

A one-dimensional numerical model described in the Methods was run with parameters selected to provide the best fits to the SAM-TLS data. The largest effect on the quality of the fit is exerted by three thermophysical regolith parameters: (1) the adsorption enthalpy, $ΔH$, (2) the uptake to evaporation coefficient ratio, $γ/η$, which appears as a pre-factor in the adsorption equilibrium coefficient,

<footnote n="1">¹Centre for Research in Earth and Space Science, York University, Toronto, Ontario, Canada. ²University of Colorado, Boulder, CA, USA. ³University of Michigan, Ann Arbor, MI, USA. ⁴Université Paul Sabatier, Toulouse, France. ⁵NASA-Goddard, Greenbelt, MD, USA. ⁶Aeolis Research, Pasadena, CA, USA.</footnote>

<footnote n="7">⁷Caltech-JPL, Pasadena, CA, USA. *e-mail: jmoores@yorku.ca</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1** | Full diffusive-adsorptive model result for methane in the near-surface atmosphere³³ at Gale. These results (black line with one s.e.m shown in cyan) correspond to a subsurface seep of 2.8 × 10⁻¹⁶ kg m⁻² s⁻¹, an EADT of 30 sols, a methane enthalpy of 31.5 kJ mol⁻¹ and a γ/η ratio of 10⁰. SAM-TLS measurements of methane with one s.e.m error bars are shown while the rover travelled along the floor of Gale Crater (purple circles) and after it began its ascent of Aeolus Mons (yellow circles). The quality of the fit of the model to the purple circles is described by the inset using the same error bars as in the main figure. The dashed blue line is a surface-only Arrhenius relationship. SEM, standard error of mean.</figcaption></figure>

and (3) the microseep mass flux. An example of a well-fit solution is shown in Fig. 1 with a comparable Arrhenius curve with three free parameters shown for reference. For this fit, the predicted value of the equilibrium absorption constant (equation (7)) was calculated for equal uptake and evaporation coefficients. This, in turn, required modelling the adsorption enthalpy as 31.5 kJ mol⁻¹, higher than terrestrial laboratory measurements run on Mars-analogue soils¹⁸ of between 16 and 20 kJ mol⁻¹ but consistent with estimates²³ of the methane enthalpy required to produce the 7 ppbv plumes¹³ through deliquescence-induced regolith desorption. The higher values needed in the model could suggest that methane on the Martian surface is bound more tightly than observed in the laboratory. For this particular case, the strength of the seep was modelled to be 2.8 × 10⁻¹⁶ kg m⁻² s⁻¹.

The strength of the seep is not an intuitive parameter. However, each value of the seep strength can be related through conservation of mass (equation (9)) to a single value of the near-surface effective atmospheric dissipation timescale (EADT) of methane. We define the EADT as the exponential timescale over which methane is removed from the portion of the atmosphere in continuous direct contact with the surface within the model, as given by the nighttime planetary boundary layer (PBL) thickness²⁴. The EADT is therefore different from the photochemical lifetime and the mixing timescale for a point source emission. When considering the EADT, the reader should note that the EADT is sensitive to the unknown concentration of methane in the air that mixes with the PBL above Gale. If exterior air mixing with Gale air is rich in methane, the EADT will be enhanced in the model whereas exterior air that is poor in methane will result in EADT values that are suppressed. For the seep strength modelled in Fig. 1, the EADT is 30 sols.

The fits obtained with the model compare favourably to the SAM-TLS observations. Computing the reduced chi squared ($\chi_v^2$) statistic (inset of Fig. 1, see ref. ²⁵ for a full definition of $\chi_v^2$) gives a value of 0.91, much closer to the $\chi_v^2=1$ value corresponding to an ideal fit, compared with previous work⁴ in which the best value was $\chi_v^2=2.07$. Based on the $\chi_v^2$ statistic, the model has a goodness of fit

probability of 0.65. Note that values of the probability greater than 0.1 are considered excellent fits²⁵. This statistical methodology²⁵ was chosen to allow direct comparison with the fits found in previous work¹. The greatest divergence between the model and observations is seen in the second half of the year when dynamical modelling²⁶,²⁷ and observation of clouds and dust within the crater²⁸,²⁹ suggest that regional mixing becomes more important.

The points marked in yellow are excluded, with the exception of the point at solar longitude $L_s=265.3°$, which is during the transition between roving on the floor of Gale Crater and the beginning of the ascent of Aeolus Mons. These points were excluded based on observations from other instruments that this change in terrain corresponds to a change in the dynamical state of the overlying atmosphere²⁶. A potential cause of this change is an increase in the daily minimum PBL depth, consistent with models that predict enhanced atmospheric mixing on the slopes of Aeolus Mons²⁶ and as demonstrated by an observed increase in dust devil activity³⁰. While there are, as yet, too few observations to be well fit by a modified model, if the PBL depth is increased by 70%, the new observations can be reconciled with the diffusive-adsorptive model. Note that a more extensive sensitivity analysis (see Methods) indicates that no matter which set of points are selected, an adsorption enthalpy of close to 31.5 kJ mol⁻¹ is preferred.

Finally, while the surface-based Arrhenius curve may appear to capture a great deal of the variation (blue dashed line in Fig. 1), the difference between it and the subsurface diffusive-adsorptive model is significant. Computing the goodness of fit for a surface-only model¹ given by $(3 \times 10^4)^T \exp(-32 \text{ kJ mol}^{-1}/RT_{\text{surf}}) + 0.1 \text{ ppb v}$, where $R$ is the universal gas constant and $T$ is the temperature in K (see equation (3))⁴ produces a fit of $\chi_v^2 = 2.50$ and therefore the probability is 0.0328, which is lower than 0.1, suggesting the fit should be rejected²⁵. The major reason for this difference is the phase lag seen between $L_s=100°$ and 200°. In a purely surface-based Arrhenius model, methane emission from the surface follows the temperature closely whereas for the diffusive-adsorptive model, time is required for this methane to travel through the subsurface, giving rise to
<figure><img src="image_3.png" /><figcaption>**Fig. 2** | An example of the sensitivity analysis showing the variation in the quality of fit with changing values of $\Delta H$ and $\gamma/\eta$ ratio. The best fits are those with the lowest values of the $\chi^2$ parameter, here coloured in yellow. These define a line of best fit in this EADT slice, or a plane in a three-dimensional representation that includes seep mass flux/EADT (see Methods). Previous terrestrial laboratory work¹⁸ is shown as a white 'X' with one s.e.m.</figcaption></figure>

the phase lag. Furthermore, the surface-based Arrhenius model requires both an unreasonably large value of the pre-exponential factor, indicating an unphysical amount of regolith in direct contact with the atmosphere (corresponding to tens of kilometres depth) as well as a constant background of 0.1 ppbv. Neither of these adjustments are required in the diffusive-adsorptive model.

## Uniqueness of the solution

A full sensitivity analysis of the result can be found in the Methods. However, the defining feature of this analysis is that the best fits to the SAM-TLS data lie along a plane in a three-dimensional space defined by (1) adsorption enthalpy, $\Delta H$, (2) uptake to evaporation coefficient ratio, $\gamma/\eta$, and (3) seep mass flux or by extension EADT. A representative example of this behaviour can be seen in Fig. 2, which shows a slice along a plane of constant seep mass flux/EADT. In this figure, the line of best-fitting values for the model lies along a line in the $\Delta H-\gamma/\eta$ plane given by:

$$ \Delta H = -4,300 \log \left( \frac{\gamma}{\eta} \right) + 31,000 \text{ J mol}^{-1} $$

This line is valid for lifetimes above 17 sols, the adsorption timescale of atmospheric methane along the line of best fit and consistent with the timescale for absorption of a plume into an initially methane-free surface layer. For values of the EADT lower than the mean absorption timescale of approximately 17 sols, the model is challenged to reproduce the observed variation except at significantly higher values of $\Delta H$ or $\gamma/\eta$ that imply faster kinetics of adsorption and desorption. Without these faster kinetics, the methane in the atmosphere disappears faster than it can interact with the surface.

## Setting limits on methane seepage at Gale

If microseepage occurs constantly and evenly everywhere across the surface of Mars, then methane can dissipate only by oxidation in the lower atmosphere or photolysis above as air transported in from outside the crater would contain an equal methane concentration

compared to air within the crater. This would produce an EADT value equal to the photochemical lifetime of methane at Mars, estimated at 329 yr (ref. ⁶), and the corresponding seepage strength would be $7.2 \times 10^{-20} \text{ kg m}^{-2} \text{ s}^{-1}$. At the opposite end of possibility, if microseepage of methane is restricted to Gale Crater only, with methane-free air outside of the crater, then the atmospheric lifetime would be given by the dissipation time for a point source, of order 1 sol (refs. ²⁶,³¹) and the seepage strength would be $8.4 \times 10^{-15} \text{ kg m}^{-2} \text{ s}^{-1}$. These limits may be narrowed by the application of the numerical model from the Methods and comparing the seepage rates with geological constraints on the emitting area.

### Upper limit for methane microseepage.
The changes in measured methane in different parts of Gale, for instance, on the floor of the crater as opposed to the slopes of Aeolis Mons, are suggestive of incomplete horizontal mixing. As such, a minimum EADT on the order of tens of sols is likely and the model gives relatively better fits for such values. This effect is observed in our sensitivity analysis where for atmospheric lifetimes below 17 sols, the shape of the annual methane pressure curve begins to deteriorate and progressively greater values of the adsorption enthalpy are required to fit the SAM-TLS data. This change forces the exchange timescale to shorten in lock step with the EADT. Therefore, the seep strength is likely no greater than $10^{-15} \text{ kg m}^{-2} \text{ s}^{-1}$, or approximately $3 \times 10^{-5} \text{ tonnes km}^{-2} \text{yr}^{-1}$. This is significantly lower than the strength of terrestrial microseepage⁴ of $1-1,000 \text{ tonnes km}^{-2} \text{yr}^{-1}$.

Based on the SAM-TLS results⁴, the mean atmospheric concentration can be no more than ~0.41 ppbv and when this is combined with the photochemical lifetime of 329 yr (ref. ⁶), only ~11 tonnes of methane would need to be supplied by microseepage annually for the atmosphere to remain in equilibrium. For the upper limit of the seepage rate at Gale found above, a total of $3.8 \times 10^5 \text{ km}^2$ or 0.26% of the surface of Mars would be needed to supply the needed methane. Thus, if Gale is a typical seep, then only a very small fraction of the Martian surface need be active to supply the planetary methane budget.

### Lower limit for methane microseepage and the potential for fossil plumes.
Above an EADT of 17 sols, the behaviour shown in the sensitivity analysis (see Methods) stabilizes up to the photochemical lifetime of methane on Mars⁶. Furthermore, models run explicitly with a microseepage rate of zero also can yield matches to SAM-TLS results⁴ up to a $\chi^2$ value of 1.58 (see Methods). However, such zero microseepage values require substantial subsurface adsorbed methane and would degrade over time as the methane is released back into the atmosphere and dispersed. If not provided from below by microseepage, it is possible that such adsorbed methane could be the result of previous plumes of methane emitted at Gale or nearby that would have adsorbed into the regolith down to a substantial depth and whose slow re-release into the atmosphere is being controlled by adsorption and diffusion. This possibility may partially explain the rapid disappearance³² of some of the previously observed methane plumes⁸,¹³ on Mars and would be evinced by a variation in the amplitude of the seasonal cycle at different distances from large plume sources. Such behaviour could be visible to column methane instruments, such as the Trace Gas Orbiter's Nadir and Occultation for MArS Discovery instrument (TGO-NOMAD)³³.

### Geological constraints.
In the sensitivity analysis, the best fits are obtained for EADT on the order of a few tens of sols, near the upper range of microseepage. However, the rate at which the fits degrade at higher values of the EADT is gradual and it is therefore difficult to exclude these alternative fits, given the uncertainties in the modelled values and in the SAM-TLS measurements (for example, Fig. 1). Indeed, a zero-seepage model cannot be explicitly ruled out
based on quality of fit alone. A stronger constraint can be levied by considering the geological environments that may be producing microseeps of methane. A recent review⁷ of many potential production environments and processes that could produce methane observed that none of these processes are active over the entire surface of the planet. For instance, less than 50% of the planet has the correct mineral types to favour the Fisher-Tropsch-type (Sabatier) reaction. Furthermore, even if source materials are universally available at depth and the reactions to convert those source materials into methane are equally effective over the surface of the planet, the pathways required to allow this subsurface methane to escape are not ubiquitous. While there are major fault systems that could release subsurface methane⁷, these appear to cover a relatively small part of the surface⁷. Thus, it would seem reasonable to assume that somewhere between a few percent and a few tens of percent of the Martian surface is capable of outgassing methane through microseepage. This would suggest that if Gale is typical, methane outgassing rates should be of the order of ~$10^{-6}$ tonnes km$^{-2}$ yr$^{-1}$, consistent with an EADT of order ~100 sols.

## Implications for future measurements

Gale is not the only repository of sedimentary materials on Mars and there may be a wide variety of locations emitting subsurface methane at a variety of rates. How easily this seepage may be detected from the ground depends not only on the seepage rate, but also on the local thickness of the PBL. This arises because the observed concentration of methane in the modelled atmosphere above Gale Crater depends on the thickness of the atmosphere that is in constant contact with the surface and can exchange material on a short timescale. The minimum nighttime depth of the PBL, which has previously been modelled in global circulation models and mesoscale models²⁴,²⁶, gives the thickness of this continuously exchangeable layer. The mass balance described in equation (9) has a dependency on the PBL height in which the same rate of seepage will produce a larger concentration of methane in a thinner atmospheric layer. It has previously been noted that the PBL is especially suppressed at Gale Crater compared with other locations on Mars where PBL height can be up to five times greater²⁴,²⁶. This behaviour at Gale has since been confirmed through observation of the line-of-sight optical depth using the navigation cameras of the Curiosity rover²⁸,²⁹. As such, for the same seepage rate of methane, the methane concentrations observed elsewhere should be substantially lower.

## Online content

Any methods, additional references, Nature Research reporting summaries, source data, statements of data availability and associated accession codes are available at https://doi.org/10.1038/s41561-019-0313-y.

Received: 20 May 2018; Accepted: 24 January 2019;
Published online: 4 March 2019

## References

1. IPCC Climate Change 2014: Synthesis Report (eds Core Writing Team, Pachauri, R. K. & Meyer, L. A.) (IPCC, 2014)

2. Haberle, R. M. K., Melinda, A., De la Torre, M., Kass, D. M. & Mars Science Laboratory Science Team. Detection of Northern Hemisphere transient baroclinic eddies in REMS pressure data at Gale Crater Mars. In *American Astronomical Society, DPS Meeting 48* id.220.13 (2016).

3. Haberle, R. M. et al. Secular climate change on Mars: an update using one Mars year of MSL pressure data. In *American Geophysical Union Fall Meeting 2014* abstr. P51B-3947 (2014).

4. Webster, C. R. et al. Background levels of methane in Mars' atmosphere show strong seasonal variations. *Science* **360**, 1093–1096 (2018).

5. Webster, C. R. & Mahaffy, P. R. Determining the local abundance of Martian methane and its $^{13}$C/$^{12}$C and D/H isotopic ratios for comparison with related gas and soil analysis on the 2011 Mars science laboratory (MSL) mission. *Planet. Space Sci.* **59**, 271–283 (2011).

6. Atreya, S. K., Mahaffy, P. R. & Wong, A.-S. Methane and related species on Mars: origin, loss, implications for life and habitability. *Planet. Space Sci.* **55**, 358–369 (2007).

7. Oehler, D. Z. & Etiope, G. Methane seepage on Mars: where to look and why. *Astrobiology* **17**, 1233–1264 (2017).

8. Mumma, M. J. et al. Strong release of methane on Mars in northern summer 2003. *Science* **323**, 1041–1045 (2009).

9. Krasnopolsky, V. A., Maillard, J. P. & Owen, T. C. Detection of methane in the Martian atmosphere: evidence for life? *Icarus* **172**, 537–547 (2004).

10. Formisano, V., Atreya, S., Encrenaz, T., Ignatiev, N. & Giuranna, M. Detection of methane in the atmosphere of Mars. *Science* **306**, 1758–1761 (2004).

11. Geminale, A., Formisano, V. & Giuranna, M. Methane in Martian atmosphere: average spatial, diurnal, and seasonal behavior. *Planet. Space Sci.* **56**, 1194–1203 (2008).

12. Geminale, A., Formisano, V. & Sindoni, G. Mapping methane in Martian atmosphere with PFS-MEX data. *Planet. Space Sci.* **59**, 137–148 (2011).

13. Webster, C. R. et al. Mars methane detection and variability at Gale Crater. *Science* **347**, 415–417 (2015).

14. Webster, C. R. et al. Low upper limit to methane abundance on Mars. *Science* **342**, 355–357 (2013).

15. Vicente-Retortillo, A. et al. Determination of dust aerosol particle size at gale crater using REMS UVS and mastcam measurements. *Geophys. Res. Lett.* **44**, 3502–3508 (2017).

16. Moores, J. E., Smith, C. L. & Schuerger, A. C. UV production of methane from surface and sedimenting IDPs on Mars in light of REMS data and with insights for TGO. *Planet. Space Sci.* **147**, 48–60 (2017).

17. Martinez, G. M. et al. The modern near-surface Martian climate: a review of in-situ meteorological data from Viking to Curiosity. *Space Sci. Rev.* **212.1-2**, 295–338 (2017).

18. Gough, R. V., Tolbert, M. A., McKay, C. P. & Toon, O. B. Methane adsorption on a Martian soil analog: an abiogenic explanation for methane variability in the Martian atmosphere. *Icarus* **207**, 165–174 (2010).

19. Meslin, P.-Y., Gough, R., Lefèvre, F. & Forget, F. Little variability of methane on Mars induced by adsorption in the regolith. *Planet. Space Sci.* **59**, 247–258 (2011).

20. Grotzinger, J. P. et al. A habitable fluvio-lacustrine environment at Yellowknife Bay, Gale Crater, Mars. *Science* **343**, 1242777 (2014).

21. Stevens, A. H., Patel, M. R. & Lewis, S. R. Modelled isotopic fractionation and transient diffusive release of methane from potential subsurface sources on Mars. *Icarus* **281**, 240–247 (2017).

22. Lasue, J., Quesnel, Y., Langlais, B. & Chassefière, E. Methane storage capacity of the early Martian cryosphere. *Icarus* **260**, 205–214 (2015).

23. Hu, R., Bloom, A., Gao, P., Miller, C. E. & Yung, Y. L. Hypotheses for near-surface exchange of methane on Mars. *Astrobiology* **16**, 539–550 (2016).

24. Newman, C. E. et al. Winds measured by the rover environmental monitoring station (REMS) during the Mars Science Laboratory (MSL) rover's Bagnold Dunes Campaign and comparison with numerical modeling using MarsWRF. *Icarus* **291**, 203–231 (2017).

25. Press, W. H., Teukolsky, S. A., Vetterling, W. T. & Flannery, B. P. *Numerical Recipes in C: The Art of Scientific Computing* 2nd edn (Cambridge Univ. Press, New York, 1992).

26. Rafkin, S. C. R. et al. The meteorology of Gale Crater as determined from rover environmental monitoring station observations and numerical modeling. Part II: interpretation. *Icarus* **280**, 114–138 (2016).

27. Tyler, D. Jr. & Barnes, J. R. Mesoscale modeling of the circulation in the Gale Crater region: an investigation into the complex forcing of convective boundary layer depths. *Mars* **8**, 58–77 (2013).

28. Moores, J. E. et al. Observational evidence of a suppressed planetary boundary layer in northern Gale Crater, Mars as seen by the Navcam instrument onboard the Mars Science Laboratory rover. *Icarus* **249**, 129–142 (2015).

29. Moore, C. A. et al. A full Martian year of line-of-sight extinction within Gale Crater, Mars as acquired by the MSL Navcam through sol 900. *Icarus* **264**, 102–108 (2016).

30. Lemmon, M. T. et. al. Dust devil activity at the Curiosity Mars rover field site. In *48th Lunar and Planetary Science Conference* Contribution No. 1964, id.2952 (2017).

31. Moores, J. E. et al. Transient atmospheric effects of the landing of the Mars Science Laboratory rover: the emission and dissipation of dust and carbazic acid. *Adv. Space Res.* **58**, 1066–1092 (2016).

32. Lefèvre, F. & Forget, F. Observed variations of methane on Mars unexplained by known atmospheric chemistry and physics. *Nature* **460**, 720–723 (2009).

33. Robert, S. et al. Expected performances of the NOMAD/ExoMars instrument. *Planet. Space Sci.* **124**, 94–104 (2016).
**Acknowledgements**

J.E.M. and C.L.S. acknowledge funding for this work provided by the Canadian Space Agency's Mars Science Laboratory Participating Scientist Program.

**Author contributions**

J.E.M. led all aspects of the manuscript including writing, code development and simulation. All members of the team participated in the writing and revision of the manuscript and in discussion about the scientific interpretation of the results. C.L.S. performed the sensitivity analysis. R.V.G. and P.-Y.M. provided valuable insights into the absorption of methane on mineral grains. G.M.M. provided valuable insights on the REMS dataset. C.E.N. provided valuable insights on numerical modelling of the martian atmosphere. S.K.A., P.R.M. and C.R.W. provided valuable insights on the SAM-TLS methane measurements.

**Competing interests**

The authors declare no competing interests.

**Additional information**

**Supplementary information** is available for this paper at https://doi.org/10.1038/s41561-019-0313-y.

**Reprints and permissions information** is available at www.nature.com/reprints.

**Correspondence and requests for materials** should be addressed to J.E.M.

**Publisher's note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© The Author(s), under exclusive licence to Springer Nature Limited 2019
Methods

**Model description.** The diffusive–adsorptive model. The model employed is based on a previous diffusive–adsorptive code$^{14}$, which itself uses the mathematical convention$^{15}$ whereby the Knudsen diffusion of a gaseous species, in this case methane, through small passageways and also undergoing adsorption and desorption may be expressed as:

$$D_{\text{eff}} = \frac{D_{\text{Kn}}}{1 + \psi \frac{\partial \theta}{\partial p}} \qquad (1)$$

where

$$D_{\text{Kn}} = \frac{d}{3} \sqrt{\frac{8R_{\text{U}}T}{\pi M_{\text{CH4}}}} \qquad (2)$$

$$\psi = R_{\text{CH4}} T \rho_{\text{reg}} A_s m_{\text{ML}} \qquad (3)$$

and

$$\frac{\partial \theta}{\partial p} = \frac{k_{\text{eq}}}{(1 + k_{\text{eq}} p)^2} \qquad (4)$$ $$[1 + (k_{\text{eq}} p k_d t (1 + k_{\text{eq}} p) - 1) \exp(-k_d (1 + k_{\text{eq}} p)t)]$$

Here the effective diffusivity of the gas, $D_{\text{eff}}$ is calculated as the Knudsen diffusivity, $D_{\text{Kn}}$ modified by the capacity of the regolith to take up adsorbate, $\psi$, and the slope of the adsorption isotherm that relates the fractional coverage of the surface, $\theta$, and the partial pressure of the gas, $p$. The expression for Knudsen diffusivity, equation (2), is taken from the kinetic theory of gases and depends on the temperature, $T$, the diameter of the pores in subsurface regolith, $d$, taken to be $1 \times 10^{-5}$ m for consistency with previous models$^{19}$, the universal gas constant, $R_U$, is 8.314 $\text{J mol}^{-1} \text{K}^{-1}$ and the molecular weight of methane, $M_{\text{CH4}}$, is 0.016 kg mol$^{-1}$. While this formulation is substantially less complex than previous work$^{19}$ it offers the advantage of depending on fewer unknown factors, especially the tortuosity. $\psi$ is itself a function of the specific gas constant of methane, $R_{\text{CH4}}$, the temperature of the regolith at a particular depth, $T$, the bulk density of the regolith, $\rho_{\text{reg}}$, which is selected as $1,300 \text{ kg m}^{-3}$ for consistency with previous models$^{19}$ and the mass per unit surface area of a single adsorbed monolayer of methane$^{19}$, $m_{\text{ML}}$, $1.38 \times 10^{-7} \text{ kg m}^{-2}$. The specific surface area of the regolith, $A_s$, was allowed to vary between the commonly used laboratory value for JSC-1 Martian simulant of$^{36}$ $1.06 \times 10^5 \text{ m}^2 \text{ kg}^{-1}$ and the low $A_s$ Viking soils$^{37}$ with $A_s$ values near $1.7 \times 10^4 \text{ m}^2 \text{ kg}^{-1}$. However, for most runs, the $A_s$ was selected to be $3.0 \times 10^4 \text{ m}^2 \text{ kg}^{-1}$, an upper limit as inferred$^{38}$ from examining variation in soil hydrogen measurements by the ChemCam instrument. Finally, equation (4), which describes the slope of the adsorption isotherm, takes into account both the equilibrium change in surface coverage with changes in pressure, through $k_{\text{eq}}$, the equilibrium rate constant, which we define as:

$$k_{\text{eq}} = \frac{k_a}{p k_d} = \frac{R_a / (1 - \theta) p}{R_d / \theta} \qquad (5)$$

as well as the kinetics of how quickly that change is achieved through $k_d$, the adsorption rate constant and $k_a$, the adsorption rate constant (see equation (8)). Note that in equation (5), $R_a$ and $R_d$ are the absolute rates of adsorption and desorption. For methane at these temperatures, surface coverage is estimated to be of order $10^{-10}$ and therefore the simpler Langmuir model shown in equation (5) can be used instead of the Brunauer–Emmett–Teller (BET) formulation. Rearranging equation (5) to solve for the surface coverage when $R_d$ and $R_a$ are equal yields:

$$\theta_{\text{eq}} = \frac{k_{\text{eq}} p_{\text{CH4}}}{1 + k_{\text{eq}} p_{\text{CH4}}} \qquad (6)$$

For the pressures encountered for the seasonal cycle at Gale, $p_{\text{CH4}} \le 10^{-6}$ Pa and $k_{\text{eq}}$, which describes the ratio between the adsorption rate constant and the desorption rate constant, can be expressed by adapting two equations, 12 and 13, from previous work$^{18}$ to a pressure basis:

$$k_{\text{eq}} = \frac{\gamma}{\eta} \frac{\nu h}{4 M L_{\text{CH4}}} \left( \frac{1}{k_B T} \right)^2 \exp(\Delta H / RT) \qquad (7)$$ $$= 6.073 \times 10^{-6} T^{-3/2} \frac{\gamma}{\eta} \exp(\Delta H / RT)$$

Note that the enthalpy, $\Delta H$, in equation (7) may vary for different materials and the value of the leading constant presented assumes that all adsorption sites are equally accessible to methane molecules. Here, $\nu$ is the mean molecular speed, $M_{\text{CH4}}$ is the number of methane molecules per $\text{m}^2$ of adsorptive surface, $h$ is

Planck's constant and $k_B$ is Boltzmann's constant. Measurements made$^{18}$ near 150 K indicate a value for $\Delta H$ of approximately $18 \text{ J mol}^{-1}$, comparable to previous work with methane on a variety comparable mineral surfaces$^{39–41}$ that showed a range of $\Delta H$ from 16 to $25 \text{ J mol}^{-1}$ depending on the material on which methane was adsorbed. As will be seen later, the regolith surface on Mars may differ from laboratory analogue materials and so the numerical implementation of equation (7) allows the value of $\Delta H$ and the ratio between the uptake coefficient, $\gamma$, and the evaporation coefficient, $\eta$, to be varied as fitting parameters both at the surface and in the subsurface. $M_{\text{CH4}}$, although fixed here at a given value, is not as important in the control of atmospheric variations.

Finally, to allow for computation of equation (4), the desorption rate constant, $k_d$, was derived from equation 12 of previous work$^{18}$:

$$k_d = \gamma \frac{k_B T}{h} \exp(-\Delta H / RT) \qquad (8)$$

where $\gamma$ is the uptake coefficient as determined experimentally$^{18}$.

**Modelling temperature in the subsurface.** Temperatures as a function of time and depth are calculated in the model using the heat diffusion equation. The model temperature at the surface, used as the upper boundary condition for the model, is described by the measured rover environmental monitoring station (REMS) ground temperatures acquired over the first two full Martian years following the landing of Curiosity$^{12}$ and begin and end at $L_s=158^\circ$. These two Martian years were seen to be very comparable in terms of the seasonal cycle in ground temperature and appear to be characteristic of conditions on the relatively flat floor of Gale Crater surrounding Aeolus Mons. In the numerical implementation of the model, the bottom of the regolith is taken to be 30 m depth, which is well below the expected annual temperature wave and is assumed to remain constant at the average annual temperature$^{17}$ of 233 K. We select this depth purely as a numerical convenience as it is below the annual thermal skin depth. The bulk density of the regolith is taken to be$^{19}$ 1,300 kg m$^{-3}$, the thermal conductivity of the regolith is taken to be $0.5 \text{ W m}^{-1} \text{K}^{-1}$ and the constant pressure heat capacity, $C_p$, is taken to be $246 \text{ J m}^{-2} \text{K}^{-1}$. These values are chosen to be appropriate for a thermal inertia of $\sim 400 \text{ J m}^{-2} \text{K}^{-1} \text{s}^{-1/2}$ consistent with measurements made at Gale Crater$^{12}$.

The diffusive–adsorptive model is very sensitive to temperature, as such, errors for partial pressure and methane surface coverage outputs were calculated based on the uncertainty in temperature. Typical $1\sigma$ errors in REMS ground temperature sensor measurements are of the order $\pm 3$ K while the day-to-day $1\sigma$ variability in the recorded REMS temperatures are $\pm 2.16$ K for the maximum temperature, $\pm 5.42$ K for average temperatures and $\pm 10.42$ K for minimum temperatures$^{17}$. This trend is explained by the decreasing signal-to-noise ratio of the sensor in colder conditions. The model outputs are influenced most strongly by the methane surface coverage/pressure equilibrium that, in turn, depends on the calculation of $k_{\text{eq}}$ in equation (7). If we take the day-to-day variability as more representative of the physical variability of temperatures at Gale and average the maximum and average temperature errors, this gives a typical temperature error of $\pm 3.79$ K. Propagating this error through equation (7) suggests a $1\sigma$ error in the output pressure of $\pm 25\%$. This value is indicated in Fig. 1 as a cyan band surrounding the model solution.

**Modelling the seep and atmosphere boundary conditions.** The diffusive–adsorptive code requires boundary conditions. At the base of the regolith stack, the seep may be represented either as (1) a constant pressure of methane or (2) as a specific mass flux into the regolith stack. Given that there are no temperature variations at this level in the regolith, both are equivalent. However, since the rate of escape of methane from a seep is the more fundamental parameter and a more useful boundary condition to set directly, we use the mass flux rate of methane into the subsurface, or 'seep strength', at 30 m depth to describe the lower boundary condition.

At the surface, we approximate the atmosphere as a single reservoir that receives methane from the subsurface. At each time step, gaseous methane diffuses into the atmosphere from the regolith stack. Additionally, the equilibrium surface coverage of methane will change with temperature. In the model, only the top regolith grid point is permitted to directly communicate with the atmosphere. As such, the adsorption state of this top point changes in response to the atmospheric pressure of methane with methane movement in the subsurface described by the diffusive–adsorptive model. In the case where the atmospheric pressure of methane is smaller than the equilibrium pressure at the top of the regolith, methane is permitted to desorb following the kinetic rate (see equation (7)). In the case where the atmospheric pressure is larger than the equilibrium pressure at the top of the regolith, methane is permitted to adsorb onto the regolith as defined by the relative kinetics (see equation (7)).

Under these conditions, the model was seen to accumulate methane over time as mass seeped into the system with concentrations in the regolith and atmosphere increasing at an equal rate to the seep strength, when integrated. To model a repeatable seasonal cycle year to year, it was therefore necessary to remove methane from the system through the atmosphere. This was accomplished by prescribing a lifetime to near-surface atmospheric methane, termed the EADT, which resulted in
mass loss proportional to the amount of methane in the near-surface atmosphere at each time step. The EADT implicitly combines several effects, including mixing, advection and eventual photolysis by Lyman-α or other destructive processes. As such, it is not anticipated that the EADT value in the model will compare to overall methane lifetimes previously calculated⁶. Instead, the EADT is kept as an independent parameter that can be varied and each value of which is matched to a particular seep strength through a mass balance:

$$ \left. \frac{dm}{dt} \right|_{\text{seep}} = \left. \frac{dm}{dt} \right|_{\text{PBL}} = \frac{-\rho_{\text{CH4}} h_{\text{PBL}}}{L_{\text{CH4}}} \quad (9) $$

Here $ρ_{CH4}$ is the density of methane within the diurnal minimum PBL, $h_{PBL}$, $L_{CH4}$ is the EADT of that methane within the diurnal minimum PBL, the thickness of which was allowed to vary seasonally²¹.

**Competitive adsorption with water vapour.** The model also takes into account competitive adsorption with water vapour in the topmost layer as part of the calculation of the upper boundary condition on methane from the previous section. Based on previous work³,³⁶,⁴³, water is assumed able to displace methane molecules in the upper metre of regolith as a fraction of total sites occupied by water. Furthermore, diffusion of water into the subsurface over short timescales occurs relatively slowly³⁴ and it will be assumed that there is no net change in the adsorbed water content of the regolith below 1 m. The model assumes that water also achieves equilibrium with the atmospheric mixing ratios of water vapour derived from the REMS relative humidity sensor¹⁷, and that this equilibrium is instead representative of the coldest time of the day when adsorption is most favoured. Furthermore, it should be noted that as a result of the difference in adsorption enthalpies of methane and CO₂, which lead to a three orders of magnitude difference in their adsorption timescales¹³, co-adsorption of methane and CO₂ was not considered in the calculations described here. Experimental data on the co-adsorption of CH₄, CO₂ and H₂O at relevant pressures and temperatures for Mars would help to better understand this system.

**Initializing and running the model.** A schematic of the model is shown in Supplementary Fig. 1. To initialize the model, the seep rate was propagated vertically throughout the regolith and the atmospheric methane content was set to 0.45 ppbv. The model was then allowed to run until a repeatable seasonal cycle was obtained. The amount of time required for this steady state varied considerably depending on the choice of enthalpy and the uptake to evaporation coefficient ratio under investigation. Each time step of the model proceeds by first solving the heat diffusion equation to determine the temperature at each model grid point in the subsurface. These temperatures are calculated once per sol using the average surface temperature for that sol of the Martian year, averaged over the first two Martian years following landing. Next, the diffusive-adsorptive model is run to propagate the methane pressure, surface coverage and mass flux in time at each model grid point. The seep rate specifies the boundary condition of mass flux at depth while the surface boundary condition is established though atmospheric exchange of methane. The next section (Model sensitivity analysis) describes the range of runs performed and the sensitivity of the results obtained with the model to the parameters selected.

**Model sensitivity analysis.** Sensitivity parameters in the diffusive-adsorptive model. To understand the sensitivity of the model to changes in various parameters and to examine the quality of fit at any one EADT slice, the two variables with the largest effect on the quality of fit, namely ΔH and the leading coefficient of equation (7) through the uptake to evaporation coefficient ratio (γ/η) were varied. Changes in EADT/seep rate were also considered. These plots show the quality of fit as a function of the leading coefficient and ΔH for different values of the EADT. Example EADT slices are shown in Supplementary Fig. 2, which shows that each EADT slice has a line of best fit (shown in yellow) trending from the upper left to lower right.

There is no strong preference for any particular EADT with the exception of excluding the 1 and 10 sol cases where the line of best fit lies above where it is found in the other EADT slices. In these cases, the EADT is lower than the anticipated adsorption timescale. In fact, the range of good fits is seen to move to higher values of the enthalpy for which the adsorption timescale is smaller as the system attempts to equilibrate over shorter timescales.

Changes in water vapour enthalpy were found to have a relatively minor effect on the result. Similarly, changes in the specific surface area were also found to be a minor factor, which could be incorporated with the variation in the leading coefficient, as the equilibrium adsorption constant is often multiplied by the specific surface area when calculating the diffusivity, for instance, in equation (1).

**Testing alternative hypotheses for changes with rover movement.** The base model (termed 'hypothesis 2' below) assumes that all points where the rover was travelling over the floor of Gale Crater are equally valid points of comparison. These are the purple points in Fig. 1. There are, however, alternative groupings of the included and excluded points that are reasonable to consider, with reference to the colours of Fig. 1:

(1) All points are equally valid measurements of the methane cycle at Gale. This requires that the atmosphere homogenizes over spatial scales comparable to the entire rover traverse and that changes in concentration with altitude are small. All points are purple points.

(2) The three points acquired after climbing began have a different PBL height and should be excluded.

(3) The point near $L_s=330°$ is part of a plume event and is not a valid measurement, which requires exclusion.

(4) The point near $L_s=300°$ is artificially low due to the greater inter-crater mixing during this time of the year. Other datasets also record enhanced cloud activity and dust at this time of year.

(5) The first point near $L_s=220°$ is included, as the rover only just began climbing, while the points near $L_s=300°, 10°$ (during climb) and $L_s=330°$ (part of plume) are excluded.

The results of this analysis for all five hypotheses are given as Supplementary Fig. 3. No matter which hypothesis we selected, the line of the best fit remained a feature of all figures produced and remained in the same location.

**Zero-flux models.** At longer values of the EADT, the sensitivity plots shown in Supplementary Fig. 2 indicate that good fits stabilize and good fits may be obtained for models that explicitly set the seepage flux to zero and the EADT to the maximum reasonable timescale of 329 yr, consistent with the photochemical lifetime of methane⁶. Supplementary Fig. 4 shows an example of a good fit at zero flux and Supplementary Fig. 5 shows a corresponding sensitivity analysis. Good fits are not obtainable in the zero-flux case except with extremely long values of EADT that are greatly in excess of a Mars year.

A critical parameter in the zero-flux case is the initial methane loading of the regolith. Sufficient methane must be present in the subsurface whether emplaced from below by microseepage or from above by adsorption from high concentration plume events. The regolith can potentially store considerably more methane than does the atmosphere; however, the mean value of the seasonal cycle will decrease over time as methane escapes the system. Therefore, to have a repeating cycle, a stabilizing source of methane is needed. However, the methane cycle at Gale appears to be relatively repeatable year-to-year, changing only when the rover began to climb Aeolis Mons. Such repeatability points towards a more constant source and removal mechanism.

Furthermore, it is unclear whether a plume would be able to push sufficient methane sufficiently deep into the subsurface to produce the same surface cycle that is observed¹. A key feature of the pattern observed in this work is a phase lag due to the interaction between diffusing methane and the thermal wave. This phase lag in emitted methane is a key difference between the cycle seen here and simple surface models, such as the surface models discussed in Fig. S39 of past work¹, and explains in part the better fits that we observe. If a plume is only adsorbed on a thickness of regolith thinner than the annual temperature skin depth, the regolith will not be sufficiently charged. As such, while zero-flux models cannot be ruled out, the microseepage model for regolith supply from below seems more likely.

## Code availability

The code used to generate the figures shown in this paper and in the Supplementary Information are available at www.yorku.ca/jmoores/MarsMethane1D.tar.gz and there are no restrictions on availability.

## Data availability

All data used in this study, REMS and SAM-TLS values, are available on the Planetary Data System (https://pds.nasa.gov). For methane data, the reader is further directed to ref. ⁴ where reduced tabulated values are available. All other data used in the production of this paper are listed in the text, and the figures were generated from this data using the code described under Code availability.

## References

34. Moores, J. E., Smith, P. H. & Boynton, W. V. Adsorptive fractionation of HDO on JSC MARS-1 during sublimation with implications for the regolith of Mars. *Icarus* **211**, 1129–1149 (2011).

35. Chevrier, V., Ostrowski, D. R. & Sears, D. W. G. Experimental study of the sublimation of ice through an unconsolidated clay layer: Implications for the stability of ice on Mars and the possible diurnal variations in atmospheric water. *Icarus* **196**, 459–476 (2008).

36. Pommerol, A., Schmitt, B., Beck, P. & Brissaud, O. Water sorption on martian regolith analogs: thermodynamics and near infrared reflectance spectroscopy. *Icarus* **204**, 114–136 (2009).

37. Ballou, E. V., Wood, P. C., Wydeven, T., Lehwalt, M. E. & Mack, R. E. Chemical interpretation of Viking Lander 1 life detection. *Nature* **271**, 644–645 (1978).

38. Meslin, P.-Y. et al. Soil diversity and hydration as observed by ChemCam at Gale Crater, Mars. *Science* **341**, 1238670 (2013).
39. Golden, T. C. & Sircar, S. Gas adsorption on silicalite. *J. Colloid. Interface Sci.* **162**, 182–188 (1994).

40. Pereira, P. R., Pires, J. & de Carvalho, M. B. Adsorption of methane and ethane in zirconium oxide pillared clays. *Sep. Purif. Technol.* **21**, 237–246 (2001).

41. Zhang, S. Y., Talu, O. & Hayhurst, D. T. High pressure adsorption of methane in NaX, MgX, CaX, SrX and BaX. *J. Phys. Chem.* **95**, 1722–1726 (1991).

42. Vasavada, A. R. et al. Thermophysical properties along Curiosity's traverse in Gale Crater, Mars, derived from the REMS ground temperature sensor. *Icarus* **284**, 372–386 (2017).

43. Zent, A. P., Howard, D. J. & Quinn, R. C. H₂O adsorption on smectites: application to the diurnal variation of H₂O in the Martian atmosphere. *J. Geophys. Res.* **106**, 14667–14674 (2001).