LETTER

doi:10.1038/nature12400

A rigid and weathered ice shell on Titan
D. Hemingway1, F. Nimmo1, H. Zebker2 & L. Iess3

a4

20Gravity (mGal)

–2–4

b80
60
40
200–20Elevation (m)Several lines of evidence suggest that Saturn’s largest moon, Titan,
has a global subsurface ocean beneath an outer ice shell 50 to200 kilometres thick1–4. If convection5,6 is occurring, the rigid portion of the shell is expected to be thin; similarly, a weak, isostaticallycompensated shell has been proposed7,8 to explain the observed
topography. Here we report a strong inverse correlation between
gravity3 and topography9 at long wavelengths that are not dominatedby tides and rotation. We argue that negative gravity anomalies(mass deficits) produced by crustal thickening at the base of theice shell overwhelm positive gravity anomalies (mass excesses) produced by the small surface topography, giving rise to this inverse
correlation. We show that this situation requires a substantiallyrigid ice shell with an elastic thickness exceeding 40 kilometres, andhundreds of metres of surface erosion and deposition, consistent
with recent estimates from local features10,11. Our results are therefore not compatible with a geologically active, low-rigidity ice shell.
After extrapolating to wavelengths that are controlled by tides androtation, we suggest that Titan’s moment of inertia may be evenhigher (that is, Titan may be even less centrally condensed) than iscurrently thought12.
Three solutions for Titan’s low-order gravity field3,12 have been obtained
via radio tracking of the Cassini spacecraft, with good constraints up tospherical harmonic degree 3. Topography data have been obtained
through a combination of radar altimetry and analysis of the overlapping portions of Cassini radar images13. Spherical harmonic models9,14
have been generated to represent the topography globally. Because ofdata gaps, the coefficients derived can depend on the maximum degree
allowed for the harmonic expansion, but are stable up to degree 6(Supplementary Information).
The degree-3 signal is useful because, although it is small comparedto those at degrees 2 and 4, it is not directly affected by tides, rotation or
tidal heating7. Figure 1 illustrates the negative correlation between thedegree-3 gravity (Fig. 1a) and topography (Fig. 1b) signals for one setof gravity and topography data. We employ an admittance analysis15,16,
which measures the gravity-to-topography ratio in a way that is not biased
by uncorrelated signals in the gravity data (Supplementary Information).
For example, there may be gravity anomalies sourced from the deeper
interior, but they should not show a strong correlation with surface
topography. The degree-3 admittance appears to be substantially negative, and a Monte Carlo analysis shows that this result is robust to theestimated uncertainties (Fig. 2).
Negative admittances are rare. They can occur in convecting systems
with strong viscosity stratification17, but it is unclear why Titan’s ice
shell would have such layering, and the shell is sufficiently thin that—asfor Enceladus18—convective features should have much shorter horizontal length scales than degree 3.
Titan’s long-wavelength surface topography probably arises fromsmall (,1%) variations in shell thickness7, which may be related totidal heating (Supplementary Information). Such variations in shell
thickness can give rise to negative admittances only if the shell isrigid and loading occurs primarily from below, leading to roots atthe base of the ice shell that are large compared with the isostatic case(Supplementary Information). Here we model this situation, taking–40–60–80

Figure 1 | Titan’s degree-3 gravity and topography. a, Gravity field derived
from the potential coefficients of the SOL1a gravity field representation of ref. 3(multi-arc analysis, 3 3 3 gravity field). b, Spherically referenced topography
based on degree-6 harmonic expansion9 (Supplementary Information). The
two signals display a strong negative correlation (with correlation coefficientc 5 –0.61) and give rise to an admittance of –32 mGal km21, based on a Monte
Carlo analysis accounting for the uncertainty in the two signals (SupplementaryInformation). Maps are shown in Mollweide projection, centred on the antiSaturnian point (180u longitude).

into account loading from both the top (for example, sedimentation orerosion at the surface) and bottom (for example, freezing or melting at
the base of the ice shell), and calculate the resulting deformation using
the equations of flexure on a thin elastic shell19–21.
If the observed surface topography h and the amount of top loading
dt are specified, and if top and bottom loads are in phase, it can beshown (Supplementary Information) that the degree-l admittance isgiven by1
32 0dtl1{
lz1dd76 B
h z tCð1Þ
Z(l)~4pGrcrust 41{@A 1{5Ch2lz1Rwhere d is the mean ice-shell thickness, R is Titan’s mean radius andrcrust is the crustal (that is, ice shell) density. The quantity dt is theamount of material added or subtracted at the top of the shell: negativevalues of dt indicate removal of material (erosion) from pre-existinghighs and deposition at lows. Here C is the compensation factor atdegree l; when C 5 1, the topography is fully compensated (isostatic)
and admittance is necessarily positive. When C , 1, the load is partly

1Department of Earth and Planetary Sciences, University of California Santa Cruz, 1156 High Street, Santa Cruz, California 95064, USA. 2Departments of Geophysics and Electrical Engineering, Stanford
University, Stanford, California 94305, USA. 3Dipartimento di Ingegneria Meccanica e Aerospaziale, Università La Sapienza, 00184 Rome, Italy.

5 5 0 | N AT U R E | VO L 5 0 0 | 2 9 AU G U S T 2 0 1 3©2013 Macmillan Publishers Limited. All rights reserved

LETTER RESEARCH
Degree-3 admittance (mGal km–1)

Topography–40–200Gravity

a 400

100+2080–27 ± 12 (–0.49)

SOL1bSOL2Figure 2 | Admittance estimates for nine sets of gravity and topography
data. Three distinct gravity models3 and three distinct topography models9(truncated at harmonic degrees 4, 5 and 6) were used to produce a total of nine
separate admittance estimates. Each estimate is based on a Monte Carloanalysis in which the admittance was computed for each of 100,000 distinct setsof gravity and topography coefficients, distributed according to the 1s
uncertainties in the model coefficients. For each Monte Carlo analysis,
diamonds are plotted at the mean admittance estimate, and 1s error bars areshown to represent the distribution. To the right of the error bars are themean 6 standard deviation of the admittance estimates in milligals perkilometre, followed by the correlation coefficients c in parentheses.

–20

150–40100–6050–80

0–1,500–1,000–500
Top load, dt (m)

0500

b 600

–100

0500

–200
00–2–400

400

–600

300

–800
–1,000200

Top load, dt (m)supported by elastic flexure of the ice shell and admittance may be
negative for sufficiently small values of C or sufficiently negative valuesof dt/h. In the case of degree 3, and assuming reasonable ice-shellproperties (Supplementary Information), C is so large that admittancecannot be negative unless dt/h is negative (that is, when erosion has
occurred in areas of positive topography).
The black line plotted in Fig. 3a indicates the admittance corresponding to Fig. 1 and shows that more than 100 m of surface erosion isrequired even for very large (,400 km) elastic thicknesses. If the elastic
thickness is less than 40 km, more than ,1 km of erosion is required.
Figure 3b plots combinations of d, T and dt that satisfy the observedadmittance. For an ice-shell thickness of 200 km, for example, morethan 200 m of erosion are required. As the elastic layer decreases inthickness, there is an increase in both the implied amount of erosion,and the absolute uncertainty in that estimate (Supplementary Information). Using Z(3) 5 –32 6 11 mGal km21, the uncertainty in theadmittance estimate produces approximately 630% uncertainty inthe erosion estimate (where one galileo is defined as one centimetre
per second squared).
The survival of large-impact basin rims and other topographic features10,11
rules out erosion amplitudes greater than ,1 km, implying a shell elastic
thickness greater than 40 km (Fig. 3b). This high rigidity could be theresult of a cold (ammonia-rich) subsurface ocean22,23, a clathrate-richshell4 or a low heat flux from the interior (Supplementary Information).
The addition of a rigid shell would slightly reduce the k2 tidal Love
number, but not enough to conflict with the measured value3 (Supplementary Information). The implied elastic thickness rules out a vigorously
convecting shell, limits the potential for widespread cryovolcanism22,24,
and permits the survival of lateral shell thickness variations (Supplementary Information). To generate the observed topography, ourmodel requires shell thickness variations at degree 3; possible sources
include tidal heating in a laterally heterogeneous shell25, and redistribution of material via non-Newtonian lateral flow26 (Supplementary Information). A rigid conductive shell is also only weakly dissipative, potentially
helping to explain Titan’s high present-day orbital eccentricity4,7.
Our lower bound on the extent of erosion and deposition (200 m)
over global length scales is compatible with estimated local amounts
based on impact crater degradation10,11, but somewhat larger than estimates from fluvial incision27. The implied vertical erosion/deposition
rate is of the order of a metre per million years, comparable to the lower
end of aeolian deposition rates measured on Earth28; erosion may be

000

SOL1a–16 ± 12 (–0.34)

200–4–32 ± 11 (–0.61)
Degree ≤ 6

SOL2

2000

–15 ± 8 (–0.33)

SOL1aSOL1b–6–8 ± 9 (–0.21)

40250

–800

–19 ± 8 (–0.45)
Degree ≤ 5SOL260300
Admittance (mGal km–1)–22 ± 10 (–0.43)

SOL1b–1,6–1,4 00–1,20000–1,000–13 ± 10 (–0.30)

Shell thickness, d (km)Degree ≤ 4

SOL1a

Elastic shell thickness, T = d (km)350

–24 ± 9 (–0.51)–1,200100–1,400
0050

100

150

200

250

300

350

400

–1,600

Elastic thickness, T (km)Figure 3 | Model predictions of admittance and erosion. a, Degree-3
admittance predicted by our model for various combinations of elastic
thickness T and top load dt, here assuming that T equals the shell thickness d.
Larger values of T correspond to smaller values of C (that is, less compensation).
Negative top load indicates erosion at topographic highs and deposition at
topographic lows. The solid black line indicates the admittance corresponding
to Fig. 1 (–32 mGal km21); the dashed lines indicate the 1s uncertainty on thatestimate. Whereas the rheology of ice implies that T < 0.5d (SupplementaryInformation), assuming that T 5 d leads to conservative estimates of erosion.
b, Top load amplitude, dt, required to produce the observed degree-3
admittance for various combinations of total shell thickness d and elastic
thickness T. Dashed white lines highlight the likely range of Titan’s ice shell
thickness, 50–200 km. Both a and b were generated assuming a degree-3
topographic amplitude h of 66 m.

occurring by physical comminution, dissolution or sublimation. Because
Titan’s topography is high at the equator, we predict maximum erosion
occurring at the equator and sediment transport, via fluvial or aeolianprocesses, predominantly towards the poles. An alternative possibility
is transport in the vapour phase, if the mobile material has a sublimation temperature close to that of the surface temperature.
The degree-2 gravity signal is dominated by tides and rotation, andhas been used to determine Titan’s moment of inertia12. If the processes identified here are also operating at degree 2, then the measured
degree-2 gravity signal is an underestimate of that caused by tides androtation (Supplementary Information). As a result, Titan’s fluid Love
number may be larger than the value estimated directly from theobserved gravity field, indicating that Titan may be even less centrally
condensed than previously thought (Supplementary Information).
This result reinforces the need to understand how so large a body could
have accreted without undergoing more complete differentiation12,29.
2 9 AU G U S T 2 0 1 3 | VO L 5 0 0 | N AT U R E | 5 5 1©2013 Macmillan Publishers Limited. All rights reserved

RESEARCH LETTER
If the ice shell is sufficiently rigid (that is, the elastic layer exceeds
half of the total shell thickness), then the admittance at degree 4 should
also be negative. Assuming a fully elastic shell (that is, T 5 d) of thickness 200 km, we predict the degree-4 gravity field amplitude and thedegree-4 admittance to be 1.4 mGal and –5 mGal km21, respectively(Supplementary Information). Future Cassini gravity fly-bys will improvethe determination of the degree-4 gravity field by a factor of two andwill therefore provide a partial test of this prediction. If a non-negative
degree-4 admittance is observed, that may be an indication that asmaller portion (less than half) of the ice shell is elastic.
Received 8 May; accepted 20 June 2013.
1.

Béghin, C., Sotin, C. & Hamelin, M. Titan’s native ocean revealed beneath some
45 km of ice by a Schumann-like resonance. C. R. Geosci. 342, 425–433 (2010).
2. Bills, B. G. & Nimmo, F. Rotational dynamics and internal structure of Titan. Icarus
214, 351–355 (2011).
3. Iess, L. et al. The tides of Titan. Science 337, 457–459 (2012).
4. Tobie, G., Lunine, J. I. & Sotin, C. Episodic outgassing as the origin of atmosphericmethane on Titan. Nature 440, 61–64 (2006).
5. Mitri, G. & Showman, A. P. Thermal convection in ice-I shells of Titan andEnceladus. Icarus 193, 387–396 (2008).
6. Tobie, G., Grasset, O., Lunine, J. I., Mocquet, A. & Sotin, C. Titan’s internal structure
inferred from a coupled thermal-orbital model. Icarus 175, 496–502 (2005).
7. Nimmo, F. & Bills, B. G. Shell thickness variations and the long-wavelength
topography of Titan. Icarus 208, 896–904 (2010).
8. Choukroun, M. & Sotin, C. Is Titan’s shape caused by its meteorology and carbon
cycle? Geophys. Res. Lett. 39, 1–5 (2012).
9. Zebker, H. A. et al. Titan’s figure fatter, flatter than its gravity field. AGU Fall Meet.
abstr. P23F–01 (2012).
10. Neish, C. D. et al. Crater topography on Titan: implications for landscape evolution.
Icarus 223, 82–90 (2013).
11. Moore, J. M., Howard, A. D. & Schenk, P. M. Bedrock denudation on Titan: estimatesof vertical extent and lateral debris dispersion. Lunar Planet. Sci. Conf. XXXXIIII,abstr. 1763 (2013).
12. Iess, L. et al. Gravity field, shape, and moment of inertia of Titan. Science 327,1367–1369 (2010).
13. Stiles, B. W. et al. Determining Titan surface topography from Cassini SAR data.
Icarus 202, 584–598 (2009).
14. Zebker, H. et al. Size and shape of Saturn’s moon Titan. Science 324, 921–923(2009).
15. McKenzie, D. The relationship between topography and gravity on Earth andVenus. Icarus 112, 55–88 (1994).

16. Wieczorek, M. A. Gravity and topography of the terrestrial planets. Treat. Geophys.
10, 165–206 (2007).
17. Richards, M. A. & Hager, B. H. Geoid anomalies in a dynamic Earth. J. Geophys. Res.
89, 5987–6002 (1984).
18. Roberts, J. H. & Nimmo, F. Tidal heating and the long-term stability of a subsurface
ocean on Enceladus. Icarus 194, 675–689 (2008).
19. Kraus, H. Thin Elastic Shells (Wiley, 1967).
20. Turcotte, D. L., Willemann, R. J., Haxby, W. F. & Norberry, J. Role of membranestresses in the support of planetary topography. J. Geophys. Res. 86, 3951–3959(1981).
21. McGovern, P. J. et al. Localized gravity/topography admittance and correlation
spectra on Mars: implications for regional and global evolution. J. Geophys. Res.
107, 5136 (2002).
22. Moore, J. M. & Pappalardo, R. T. Titan: an exogenic world? Icarus 212, 790–806(2011).
23. Grasset, O., Sotin, C. & Deschamps, F. On the internal structure and dynamics ofTitan. Planet. Space Sci. 48, 617–636 (2000).
24. Lopes, R. M. C. et al. Cryovolcanic features on Titan’s surface as revealed by theCassini Titan Radar Mapper. Icarus 186, 395–412 (2007).
25. Běhounková, M., Tobie, G., Choblet, G. & Čadek, O. Tidally-induced melting eventsas the origin of south-pole activity on Enceladus. Icarus 219, 655–664 (2012).
26. Nimmo, F. Non-Newtonian topographic relaxation on Europa. Icarus 168,205–208 (2004).
27. Black, B. A., Perron, J. T., Burr, D. M. & Drummond, S. A. Estimating erosional
exhumation on Titan from drainage network morphology. J. Geophys. Res. 117,E08006 (2012).
28. Patterson, D. B., Farley, K. A. & Norman, M. D. He-4 as a tracer of continental dust: a
1.9 million year record of aeolian flux to the west equatorial Pacific Ocean.
Geochim. Cosmochim. Acta 63, 615–625 (1999).
29. O’Rourke, J. G. & Stevenson, D. J. Stability of ice/rock mixtures with applications toTitan. Lunar Planet. Sci. Conf. XXXXII, abstr. 1629 (2011).
Supplementary Information is available in the online version of the paper.
Acknowledgements We thank the Cassini radar science team, M. Manga, D. Stevenson,
R. Pappalardo and W. McKinnon for their suggestions. Portions of this work weresupported by NASA grants NNX13AG02G and NNX11AK44G.
Author Contributions F.N. initiated the effort. D.H. and F.N. developed the loading
models and analysed the results. L.I. led the development of the gravity field models.
H.Z. synthesized the topography models. All authors discussed the results andimplications and commented on the manuscript.
Author Information Reprints and permissions information is available atwww.nature.com/reprints. The authors declare no competing financial interests.
Readers are welcome to comment on the online version of the paper. Correspondenceand requests for materials should be addressed to D.H. (djheming@ucsc.edu).

5 5 2 | N AT U R E | VO L 5 0 0 | 2 9 AU G U S T 2 0 1 3©2013 Macmillan Publishers Limited. All rights reserved