Extensional tectonics are marked by shallow magma crystallization depths, whereas compressional tectonics are associated with deeper crystallization depths. This implies that variations in crystallization depths can be used to track changes in Earth’s dominant tectonic regimes through time. We therefore developed a new “pressure of crystallization” proxy based on the variation of the 176Lu/177Hf ratio in zircon. This ratio is controlled by zircon fractionation and residual garnet, and it can be used to monitor the evolution of a crystallizing magma ascending within the crust. The secular evolution of the 176Lu/177Hf ratio in zircon is characterized by cyclical oscillations that are broadly in tune with the δ18O record in zircon and with periods of continental collision and supercontinent amalgamation. The apparent mean depth of crystallization of zircon-bearing igneous rocks has decreased with time over the last ~3.0 b.y. This can be linked to shallowing of the primary crystallization depths and/or to the effect of time-integrated erosion in the geologic record. Prior to ca. 3.0 Ga, crystallization depth maxima and oscillations in apparent depth are less clear, perhaps suggesting that the nature of tectonic interactions was different in the Mesoarchean and earlier.

Magmas rise within the continental crust until they stall and solidify (Brown, 1994; Petford et al., 2000; Annen et al., 2006). Their crystallization tends to be shallower in extensional settings than in compressional settings due to contrasting normal stress gradients and conditions of magma propagation (Hogan and Gilbert, 1995; Vigneresse et al., 1999; Loucks, 2021). In geologically young terrains, the crystallization depth of magmas can be measured through geophysical approaches (e.g., Thybo and Nielsen, 2009). However, as a result of erosion and reworking of the crust through time, much of the ancient geologic record has been destroyed. As a consequence, there is little information with which to constrain the depths of crystallization and the periods of crustal extension and compression in the deep geologic past.

Zircon (ZrSiO4) crystallizes from a range of SiO2-rich magmas, and it is resistant to weathering and erosion processes, such that it is widely preserved in continental sediments across geologic time (Belousova et al., 2010; Voice et al., 2011). Many studies have explored the zircon archive to draw insights into Earth’s geodynamic regimes over time (e.g., Belousova et al., 2010; Voice et al., 2011; Dhuime et al., 2012; McKenzie et al., 2018; Balica et al., 2020). These include the use of U-Pb ages to infer the relative volume of magmatism (e.g., Schoene et al., 2012), Hf isotopes as a measure of mantle input and crustal reworking/recycling (e.g., Gardiner et al., 2016), trace-element ratios to infer mean crustal thickness (e.g., Tang et al., 2021), and oxygen isotopes as a proxy for crustal reworking (e.g., Dhuime et al., 2012; Spencer et al., 2014).

Variations in zircon 176Lu/177Hf have been documented to describe changes in magma conditions and metamorphic reactions at local scales (e.g., O’Brien and Miller, 2014; Rubatto, 2017), but global data sets of 176Lu/177Hf in zircon remain unexplored. Considering that variations in zircon composition may reflect the minerals crystallizing from the host magmas and residues during partial melting, we evaluated the relationship between zircon 176Lu/177Hf and pressure. We show that the 176Lu/177Hf ratio can be used to ascertain the depth of crystallization of zircon-bearing magmas, and we explore the implications for large-scale tectonic variations through time.

Figure 1 shows the variation of the geometric mean of 176Lu/177Hf in zircons from a worldwide data set of >120,000 analyses, subdivided into 0.1 b.y. time slices (Data Set S1 in the Supplemental Material1). The variation of the zircon 18O/16O ratios (expressed as δ18O relative to Vienna standard mean ocean water [VSMOW]) from an independent global database (Spencer et al., 2014) is plotted for comparison. It is striking that major troughs in 176Lu/177Hf correlate with peaks in δ18O, and these in turn broadly coincide with periods of supercontinent assembly. This observation was confirmed by the calculated anomaly in the running means of 0.1 b.y. intervals (inset Fig. 1), with a 0.5–0.6 b.y. cyclicity in δ18O peaks and 176Lu/177Hf troughs at 2.75–2.45 Ga, 2.1–1.9 Ga, 1.3–1.05 Ga, and 0.75–0.55 Ga (as indicated by coherent statistical variations between the two data sets; Tables S1–S5).

Higher values of δ18O in zircon can be linked to increased assimilation/reworking of supracrustal components and crustal thickening at times of major continental collisions (Spencer et al., 2014). The concomitant lower mean values of 176Lu/177Hf during periods of higher δ18O may similarly reflect crustal thickening and collisional tectonics.

Horizontal compression inhibits the ascent of magmas within the crust (Loucks, 2021), and crustal thickening expands the garnet stability field (e.g., Taylor et al., 2016; Rubatto, 2017). Given that Lu partitions preferentially into garnet, crustal thickening and an expanded garnet inventory in the crust may be influential factors. Zircon crystallization increases 176Lu/177Hf in the coexisting melts, and more fractionated magmas tend to have higher 176Lu/177Hf values. Thus, the 176Lu/177Hf of magmas and their zircons may be influenced by greater depths/pressures of crystallization in orogenic zones, and by the degree of fractionation as magmas ascend within the crust (Fig. 2A).

To assess the controls on 176Lu/177Hf in zircon from igneous magmas, we compiled studies with pressure estimates of zircon-bearing igneous rocks from which zircons were measured for 176Lu/177Hf isotopes. These studies generally relied on Al-in-hornblende to determine pressures of crystallization, while other petrologic/geochemical approaches were used to define pressures of magma generation (e.g., experimental petrology, thermodynamic modeling).

Zircon saturation-crystallization may occur from depths of magma generation to depths of final emplacement (Schoene et al., 2012), depending on the composition, temperature, and alkali content of the host crystal-rich magma reservoir (Boehnke et al., 2013; Laurent et al., 2020). At present, it is not possible to directly measure pressures of crystallization in zircon, but these can be estimated by analyzing mineral phases interpreted to be cogenetic with zircon crystallization (e.g., Barnes et al., 2019). Two or more minerals from a single paragenesis can show contrasting pressures of crystallization, and the reasons are often debated (e.g., Walker et al., 2013). To overcome this issue, our data set included a range of pressures from partial melting to emplacement depth estimates to better represent the complexity of transcrustal magmatic systems (Fig. 2A; Cashman et al., 2017).

The data spanned near-surface depths to ~2.0 GPa pressures and ranged in 176Lu/177Hf from 0.0005 (ultrahigh-pressure tonalites) to >0.002 (volcanic zircons). An exponential relationship of decreasing 176Lu/177Hf with increasing pressure is described in Equation 1:


where GPa(sample) is the pressure determined for crystallization, and 176Lu/177Hf(igneous zircon) is the geometric mean of many determinations for a sample or group of related samples from the same geologic unit. For an individual sample, the 176Lu/177Hf ratios can vary by ±25% in cogenetic zircons, whether volcanic or plutonic. As a result, ~10 or more analyses are required to estimate the mean 176Lu/177Hf and associated pressure for a given data point (Fig. 2B). One important caveat is therefore that Equation 1 is not appropriate for single-zircon pressure determinations.

Metamorphic zircons represent a modest component in the global geologic record (e.g., ~8% estimated by Balica et al., 2020). Similar to the approach for igneous zircons, we compiled data for metamorphic zircons using their age, 176Lu/177Hf, and phase equilibria–determined pressure-temperature conditions (Fig. 2B). Igneous zircons with 176Lu/177Hf < 0.0003 are scarce, and a large majority of metamorphic zircons have 176Lu/177Hf < 0.0005 (Fig. 2C). The overlap between igneous and metamorphic zircons is therefore small, and the 176Lu/177Hf ratio can be used for screening zircons of unknown origin when internal textures are ambiguous or unavailable (see the Supplemental Material).

The main sink of Hf in igneous rocks is zircon, whereas Lu is incorporated in a number of other minerals aside from garnet (e.g., amphibole, pyroxene, apatite, titanite). Lu is less compatible in zircon than Hf, and so the 176Lu/177Hf ratio is substantially lower in zircon than in the melt from which it crystallized, approximately by a factor of 5–30 (Data Set S2). Garnet has a much greater affinity for Lu, it can be abundant in the source of magmas in the middle and deep crust, and any melt generated in equilibrium with garnet is likely to have a lower 176Lu/177Hf ratio than the residue (Taylor et al., 2016; Rubatto, 2017; Gardiner et al., 2018). Other minerals such as amphibole, pyroxene, apatite, titanite, and ilmenite can influence 176Lu/177Hf in the melt, but in the presence of zircon and/or garnet, they have a minor effect on 176Lu/177Hf in the evolving magma or in zircon (see the Supplemental Material).

The 176Lu/177Hf ratio of a melt is sensitive to the melt fraction and to the proportion of garnet in the residue (Fig. 3A). Once this melt forms and starts to ascend within the crust, zircon saturation and fractionation become the dominant controls on 176Lu/177Hf in the evolving melt (assuming the lack of igneous garnet influence on the liquidus), by subtracting Hf from the system (Fig. 3B). The variation of 176Lu/177Hf with pressure is therefore thought to reflect higher amounts of garnet in the source at greater depths, combined with the increase in magma fractionation with decreasing depth (Dhuime et al., 2015). Other possible factors controlling Lu/Hf variation are temperature (Claiborne et al., 2017) and potential changes in magma types through time (Laurent et al., 2014). The impact of these two factors on the 176Lu/177Hf ratio in zircon is minor (see the Supplemental Material), and we conclude that garnet as a residual phase and zircon fractionation and increasing magma fractionation with decreasing depth are the major factors that control 176Lu/177Hf variations with depth (Fig. 2). Crystallization is the dominant control on 176Lu/177Hf through most of the depth range considered (Fig. 2B), and, as such, this is the focus of the following discussion.

We applied Equation 1 to geometric means of 176Lu/177Hf from the global zircon database to calculate the apparent mean depths of crystallization of zircon-bearing igneous rocks through time (Fig. 4, red curve). The deepest mean depths of crystallization (~30–35 km) occurred at ca. 2.6 and ca. 1.9 Ga. Shallower depths (~10–15 km) are observed for the Hadean and Paleoarchean/Mesoarchean, suggesting that changes in the depths of crystallization are not simply a function of time.

Strikingly, a marked cyclicity was observed over the last 3.0 b.y., but it is absent in the earlier record (Fig. 4). The broad coincidence of troughs in the crystallization depth curve (yellow stars, associated with the lowest 176Lu/177Hf ratios) with peaks in δ18O and periods of supercontinent assembly (Fig. 1) suggests a firstorder link with collisional processes, such as thickening and reworking of continental crust. This pattern implies either deeper exhumation through erosional and tectonic processes associated with continental collision and thickening (Blackburn et al., 2018), more efficient entrapment of magmas at deeper levels during global compressional episodes (e.g., Loucks, 2021), or a combination of both processes. Reciprocally, the peaks in the crystallization depth curve (green stars, associated with the highest 176Lu/177Hf ratios) reflect periods of maximum shallowing during tectonic extension, as in supercontinent breakup.

An unexpected feature of Figure 4 is that both the depths of deepest crystallization (yellow stars) and the peaks of shallowest crystallization (green stars) have progressively decreased since the Neoarchean, with slopes of ~0.015 and ~0.010 km/m.y., respectively. This may be interpreted as follows: (1) The red curve is a primary signal for mean crystallization depths, with deeper magmas being more abundant during older collisional events; or (2) the progressive apparent decrease of crystallization depths is due to preferential removal of uppercrustal zircons from the global zircon archive via erosion and physical destruction. The latter hypothesis does not necessarily involve changes in the primary mean depths of crystallization.

The primary mean depth of crystallization has decreased with time since 3 Ga (Fig. 4), and this may be linked to the secular cooling of the mantle (Herzberg et al., 2010) or to changes in average crustal thickness (Tang et al., 2021), although these are complex issues to address with little consensus (e.g., Balica et al., 2020; Tamblyn et al., 2022). In contrast, cumulative erosion is known to increase with time, and in young orogens in which the upper crust remains well preserved, deeper orogenic rocks are present in smaller volumes than in older orogens (e.g., Blackburn, et al., 2018). For example, the comparison between Tibet and older orogens with more deeply eroded roots (e.g., ca. 0.6 Ga Braziliano, ca. 1.1 Ga Grenvillian, and ca. 1.9 Ga Trans-Hudson collisional orogens) revealed a much greater exposure of middle to deep crust in the older orogens (e.g., Weller et al., 2021). In addition, continental collision involves significant reworking and exhumation of deeper crustal levels that may contain ancient reworked basement and continental margin sediments. This results in a progressive removal of upper-crust zircons and exposure of the grains that crystallized at greater depths.

It is likely that the 176Lu/177Hf global zircon record has been affected by both erosion and geothermal lithospheric conditions over billions of years, and a better understanding of the observed shallowing of maximum magma crystallization depths through Earth history may involve more detailed modeling of the interaction between erosion and thermal history through the cyclical generation and destruction of supercontinents. Nonetheless, it is striking that the cycles of magma crystallization depths started in response to global supercontinent activity between 3.0 and 2.5 Ga, which may be further evidence that plate tectonics became global at that time.

This project was funded by the European Research Council under the European Union’s Horizon 2020 research and innovation program (no. 817934). We thank the anonymous reviewers, Steven L. Goldstein, and Oscar Laurent for constructive and insightful comments that allowed us to clarify the manuscript. We also thank Urs Schaltegger for very helpful editorial guidance.

1Supplemental Material. Supplemental information, Figures S1–S5, Tables S1–S9, and Data Sets S1–S4. Please visit https://doi.org/10.1130/GEOL.S.21669482 to access the supplemental material and contact editing@geosociety.org with any questions.
Gold Open Access: This paper is published under the terms of the CC-BY license.