On Synchronous Supereruptions

The Youngest Toba Tuff (YTT) supereruption from Toba Caldera in Sumatra at ca. 74,000 years BP is the largest volcanic event recorded in the Pleistocene. Intriguingly, recent radioisotopic dating of the near antipodal Los Chocoyos (LCY) supereruption from the Atitlán caldera in Guatemala finds an identical age within uncertainties to that of YTT. This opens the question of whether these synchronous supereruptions may be a coincidence or could be a consequence of each other? Using the known eruptive record from the past 2 Myr, we find that the likelihood of having two near antipodal supereruptions (>1,000 km3 tephra volume) within centuries (<400 years), as suggested by volcanic proxies and annual counting layer chronology in the ice core records, is very small (0.086%), requiring a non-random cause and effect. Considering this analysis, we speculate that one potential physical mechanism that could explain the temporal relationship between these supereruptions is that seismic energy released during YTT eruption focused on the antipodal region, where concentrated stresses ultimately promoted the eruption of the perched LCY magma system (or vice versa). This supereruption “double-whammy” may thus be the more compelling source of the significant environmental impacts often attributed individually to the YTT supereruption. Improving the existing age information of YTT and LCY, and a better understanding of caldera collapse events will enable further testing of the hypothesis that synchronous supereruptions do not result by pure chance.


INTRODUCTION
Supereruptions are low-frequency but large-magnitude and unpredictable events capable of explosively ejecting ≥1 × 10 15 kg (≥450 km 3 dense rock equivalent, DRE) of silica-rich magma at geologically instantaneous timescales (magnitude (M) ≥ 8; Pyle, 2015). The recorded and expected disruptive impacts of such supereruptions range from local to global in scale (Miller and Wark, 2008;Self, 2015;Brenna et al., 2021) and may have devastating consequences to human civilization (Sparks et al., 2005). Moreover, their rarity hinders meaningful eruptions forecasting. Adding to this challenge is that whereas the trigger mechanisms of small-scale eruptions are relatively well understood (e.g., recharge, overpressure from magmatic volatile exsolution), what drives magmatic systems to generate supereruptions is still much debated (Wilson et al., 2021). It has been suggested that supereruptions and smaller-scale eruptions may be triggered by fundamentally different processes that influence their frequency and erupted magma volume (Mason et al., 2004;Deligne et al., 2010;Caricchi et al., 2014;de Silva and Gregg, 2014;Gregg et al., 2015). As such, supereruptions are extreme, enigmatic, and unpredictable natural "black swan" events that society needs to understand and anticipate (Aven, 2013).
In the last 2 Myr, 13 supereruptions are known globally (Crosweller et al., 2012) with an estimated recurrence interval of ca. 150 kyr, a timescale shorter than the frequency of large meteorite impacts (>0.6 Myr, ≥20 km crater diameter; Bland, 2005) that could have potentially similar environmental consequences (Rampino, 2002). If only the well-preserved eruption record of the last ca. 100 kyr is considered, a much shorter recurrence interval of ca. 17 kyr is obtained (Rougier et al., 2018). This is because established eruption databases are incomplete due to preservation biases (Crosweller et al., 2012), especially going back in geologic times. Thus, these recurrence FIGURE 1 | Spatial and geochronological information for the YTT and LCY projected over climate and volcanic proxy signals from the northern and southern hemisphere ice core records. (A) Location map of the Toba (Sumatra) and Atitlán calderas (Guatemala), their respective antipodes (hexagons) and approximate tephra distribution. (B) Synchronization of the YTT and LCY radioisotopic ages (±1σ) with the NGRIP oxygen isotope and sulfate concentration records around the Greenland Interstadial 20 (GI-20) and the Greenland Stadial (GS-20) as well as the sulfur isotopic compositions from the EPICA Dronning Maud Land (EDML, Antarctica), and the EPICA dome C (EDC, Antarctica) ice core records (Crick et al., 2021). Dashed gray lines indicate the candidate sulfate anomalies for the YTT supereruption in the NGRIP and Antarctic ice cores (Supplementary Figure S1). We only consider that bipolar S spikes are candidates for an equatorial supereruption. intervals could be considered maxima, and a temporal coincidence of supereruptions is not a priori unlikely. Synchronous, paired or clustered, large (M7 to M8) eruptions have been proposed within various volcanic regions (e.g., de Silva et al., 2006;Gravley et al., 2007), but synchroneity of eruptions ≥ M8 on a global scale is hitherto unknown. The discovery of two apparently synchronous recent supereruptions, the ca. 74 ka Youngest Toba Tuff (YTT;~5,300 km 3 DRE; Costa et al., 2014), Sumatra, and Los Chocoyos (LCY), Guatemala (~730 km 3 DRE; Kutterolf et al., 2016;Cisneros de León et al., 2021), has implications for the global record of supereruptions and warrants an evaluation of the randomness of paired eruptions at the colossal scale. This study contributes to this discussion and lays out future research directions.

BACKGROUND
The Toba caldera (100 × 30 km) in Sumatra and Atitlán caldera (18 × 12 km) in Guatemala are both the outcome of continental subduction-related magmatism (Rose et al., 1987;Chesner et al., 1991). These multi-cyclic calderas have produced large silicic eruptions from which the YTT and LCY stand out. The YTT supereruption represents the largest event in the Quaternary period, discharging more than 8,600 km 3 tephra (M9.1; Costa et al., 2014; Figure 1A), whereas the LCY represents the largest eruption in the American continent in the past~100 kyr (~1,220 km 3 ; Cisneros de León et al., 2021). According to deep-sea tephra settling model estimations, the total duration of the YTT was~9-14 days, whereas the LCY lasted for 20-27 days (Ledbetter and Sparks, 1979). The tephra record for the YTT lacks evidence for an initial Plinian-style eruption but distal deposits were attributed to co-ignimbrite deposits generated in ring-caldera fractures (Chesner, 1998). In contrast, the LCY presents proximal and distal tephra fallout deposits at the base of its sequence with characteristics suggesting an ultra-Plinian initial eruption (>40 km high eruptive column) followed by thick non-welded ignimbrite deposits Rose et al., 1987). The YTT supereruption is considered to be the consequence of progressive thermal maturation of its magmatic system and surrounding crust by prolonged magma influx over timescales of hundreds of thousands to millions of years (Liu et al., 2021). Zircon geochronology from the LCY also points to the protracted accumulation of silicic magmas prior to eruption (>80 kyr; Cisneros de León et al., 2021). The generation of silicic magmas associated with explosive eruptions is usually attributed to extensive mineral fractionation and melt extraction/ segregation from crystal-rich (mush) domains (Bachmann and Bergantz, 2008). The evidence of protracted zircon crystallization with diverse physicochemical signatures as well as variable glass compositions in the YTT magmas suggests mush-derived melt heterogeneity in an extensive upper crustal magma reservoir (Reid and Vazquez, 2017;Tierney et al., 2019).
The potential release of significant amounts of sulfur gases during the YTT supereruption has been linked to a major global climatic downturn that may have challenged the survival of modern humans (Ambrose, 1998;Figure 1B). Although this hypothesis is debated (Oppenheimer, 2002;see Supplementary Material), new evidence is tilting the opinion back towards a significant environmental impact around the time of the YTT supereruption (Black et al., 2021;Osipov et al., 2021). The LCY supereruption has not been studied for direct climate proxies but modelling connotes a strong potential in forcing climate change, especially by previously unconsidered sulfur coupled to halogen loading to the atmosphere Kutterolf et al., 2015;Brenna et al., 2020;Brenna et al., 2021).

Constraints on the Timing of YTT and LCY Supereruptions
The exact location of the SO 4 2− spike related to the YTT in ice cores remains highly ambiguous due to the lack of volcanic glass shards in both northern and southern hemisphere archives (Oppenheimer, 2002;Robock et al., 2009;Williams, 2012). Nevertheless, a prominent sulfate anomaly occurring in bipolar ice-core records has been tentatively correlated to the YTT (e.g., T2 sulfate spike, Figure 1B and Supplementary Figure S1; Svensson et al., 2013). However, eight other significant volcanic-derived sulfate anomalies from unknown sources (T1-T9; Figure 1B and Supplementary Figure S1; Svensson et al., 2013) also occur within uncertainty of the currently accepted radioisotopically determined eruption ages for the YTT between 73.9 ± 0.3 ka BP (1σ; 40 Ar/ 39 Ar in sanidine; Storey et al., 2012) and 75.0 ± 0.9 ka BP (1σ; 40 Ar/ 39 Ar in sanidine and biotite; Mark et al., 2014). An updated 40 Ar/ 39 Ar age for the YTT at 73.7 ± 0.3 ka (Mark et al., 2017) has been recently calculated by integrating the age data from Storey et al. (2012) and Mark et al. (2014) (two different laboratories) based on a new optimization model from Niespolo et al. (2017). Moreover, recent (U-Th)/He zircon dating of YTT has shown indistinguishable age results (75.7 ± 4.1 and 75.5 ± 4.7 ka, 1σ; Mucek et al., 2021) from those of 40 Ar/ 39 Ar ages, demonstrating the consistency of both dating methods, although (U-Th)/He zircon ages are less precise than 40 Ar/ 39 Ar ages.
A recent study indicates that the LCY supereruption (Cisneros de León et al., 2021) is a potential source for one of these significant sulfate spikes. The age of the LCY was initially estimated from δ 18 O stratigraphy at 84 ± 5 ka BP (Drexler et al., 1980) and remained radioisotopically untested for several decades due to the lack of sanidine and the ubiquitous presence of altered biotite that precluded reliable 40 Ar/ 39 Ar dating (Rose et al., 1999). However, recent radioisotopic dating applying (U-Th)/He in zircon has produced the first radiometric age for LCY of 74.8 ± 1.8 ka BP (1σ; Cisneros de León et al., 2021). This age is strikingly close to that of YTT (overlapping within 1σ error). Adding another dimension to these potential connections is the fact that Atitlán caldera is nearly antipodal (opposite side of Earth's sphere) to Toba caldera ( Figure 1A).

Timespan Between YTT and LCY
Based on the size and magnitude of the YTT and LCY supereruptions, significant deposition of sulfate on bipolar ice sheets could be expected, especially given their tropical vent locations (YTT = 35-3,500 Mt; Chesner and Luhr, 2010;Costa et al., 2014;LCY = 523 ± 94 Mt;Brenna et al., 2020). However, these sulfur masses are only minima due to the unknown sulfur mass that was partitioned into a fluid phase prior to the eruption and potentially vented during the onset of eruption. Thus, we propose that two of the nine sulfate signals are likely genetically related to the YTT and LCY ( Figure 1B and Supplementary Figure S1). A potential relative time difference between the two supereruptions can be estimated by counting the icedeposition annual layers (Svensson et al., 2013; Figure 1B). This yields a maximum of ca. 2,000 years (T1 to T9 spikes) and a minimum of ca. 25 years (T5 to T6 spikes; Svensson et al., 2013). We note that sulfate spikes T1 to T4 show relatively large-magnitude sulfur massindependent fractionation (S-MIF) isotopic signatures ( Figure 1B; Crick et al., 2021), which are indicative of large eruptions from tropical locations, the plumes of which reached altitudes at or above the ozone layer in the stratosphere. If only such spikes are considered, the potential timespan between the YTT and LCY can be further constrained to ca. 400 and 87 years, respectively; orders of magnitude shorter than the estimated recurrence interval of supereruptions (~15-100 kyr). If we consider the two peaks with the largest sulfate and S-MIF signals in EPICA Dronning Maud Land (EDML; peaks T1 and T2) as derived from the LCY and YTT supereruptions, the relative time between the two eruptions is 87 ± 6 years. Because the amplitude of the sulfate peaks is variable across the different ice cores due to spatial variation in the deposition, the relative ordering differs slightly between different datasets (see Supplementary Figure S1 for other records). The temporal proximity of the YTT and LCY raises the question of whether there is a causal relationship between these two geologically concurrent events.

METHODS
We evaluate whether the temporal clustering of large eruptions is purely random using the timing of >400 km 3 bulk volume (>M7) Quaternary eruptions (LaMEVE; Crosweller et al., 2012) that produced well-preserved deposits in the geological record. We chose the >400 km 3 threshold to increase the sample size to n = 28 for statistical analysis to avoid bias from having two coeval supereruptions (LCY and YTT) out of 13 in the past ca. 2 Myr (~10%). To assess any temporal eruption clustering in the geological record spanning the last ca. 2 Myr we calculated the coefficient of variation value (CV: the ratio of the standard deviation and the mean value for the time between two successive volcanic eruptions) for the reported eruption record (n = 28). Furthermore, we used a Monte Carlo simulation to generate 50,000 different possible synthetic eruption histories after the reported eruption record and their 1σ uncertainties (for details see Supplementary Material).

Supereruption Clustering
We find that the median value of the CV is~1.035 indicating an approximately random distribution (Supplementary Figure S2) while the median value of the mean time between eruptions is 76.3 kyr (28 eruptions in 2.054 Myr). Using the CV values obtained from the synthetic sequences of n = 28, we find that our >400 km 3 bulk volume LaMEVE distribution lies within the 5-95th percentile for a random distribution (inset Figure 2). Thus, the LaMEVE dataset does not display any significant nonrandomness/clustering at the 95% confidence limit. This conclusion is enhanced by the clear difference in the CV value between the LaMEVE dataset and synthetic eruption histories with either periodically spaced eruptions or close eruption pairs (~5% of the average time between eruption groups, Supplementary Figure S2). We describe additional results of our statistical analysis in the Supplementary Material.
Although our statistical analysis shows the randomness of large eruptions in the LaMEVE dataset as a whole, the probability of two supereruptions occurring within only 80-400 years may not be random by itself. Among 50,000 synthetic histories with random spacing between eruptions and volumes sampled from our LaMEVE dataset, we find that only 1.73% of the synthetic histories have an eruption pair that matches the YTT-LCY characteristics (Inset Figure 2). Moreover, even if we assume that the LaMEVE database is only complete for the last 100 kyr as suggested by Rougier et al. (2018) (6 eruptions with >400 km 3 in last 100 kyr, recurrence time of ca. 17 kyr), there is still only a 4.2% probability of a YTT-LCY type eruption pair (Supplementary Figure S6). Thus, the statistical likelihood for two closely spaced supereruptions is small. In order to incorporate the near-antipodal nature of YTT-LCY in our statistical analysis, we assign, for each synthetic history, an eruption location randomly from amongst the events in the LaMEVE catalog (see Figure 3). We find a very low probability of 0.086% (for the full LaMEVE dataset time period, Supplementary Figure S7) when taking into account the synthetic eruption pairs with a time separation of 80-400 years and a <3,000 km distance between the antipodal location of the first eruption in the eruption pair and the second eruption's location (comparable to Toba and Atitlán source calderas as shown in Figure 3).

Speculation About Physical Processes for Synchronous Supereruption Initiation
Our statistical model demonstrates that the probability of quasisynchronous supereruptions within millennial timescales at near antipodal locations is extremely low (0.086%), indicating that a "pure-chance" scenario as highly unlikely. A short recurrence interval between two supereruptions is unique, even considering present dating uncertainties (mainly due to the ±1.8 ka uncertainty for LCY). At present, it can only be speculated that a physical connection is strictly limited to dynamic farfield stresses since any static stress change associated with caldera collapse and eruptive venting will be restricted to the near-source region (Manga and Brodsky, 2006). In the absence of better constraints, we present a working hypothesis that accounts for the near antipodal location of Toba and Atitlán calderas and the known effects of antipodal seismic focusing, as explained in detail below. We present this as a challenge to the community hoping that future studies might shed light on this hypothesis. Geological effects including extensive crustal fracturing and surface disruption have been reported in antipodal locations after major meteorite impacts on Mercury, Mars, Pluto, and the Moon resulting from spherical focusing of impact-generated seismic energy (Schultz and Gault, 1975;Watts et al., 1991;Williams and Greeley, 1994;Denton et al., 2021). On Mars, volcanism has been linked to the antipodal focusing of seismic energy after meteorite impacts by creating a deeply antipodal fractured region that acted as magma pathways into the surface (Williams and Greeley, 1994). The idea behind antipodal deformation following shock events is that spherical seismic waves propagate through the planet and focus on the antipode, producing enhanced near-surface stress gradients leading to deformation features (Hood and Artemieva, 2008;Bowling et al., 2013). On Earth, antipodal effects from meteorite impacts have been potentially associated with triggering or enhancing volcanic activity (Hagstrum, 2005;Meschede et al., 2011;Richards et al., 2015).
Over the past few decades, there have been several observations illustrating interactions (eruption or unrest) between nearby magmatic systems (<200 km; Linde and Sacks, 1998;Biggs et al., 2016). The primary proposed physical mechanism is the effect of evacuation of the magma reservoir on other magmatic systems through either direct elastic stress changes or poro-visco-elastic response of crustal and asthenospheric mush zones related to these other reservoirs (Gonnermann et al., 2012;Tarasewicz et al., 2012;Bégué et al., 2014;Walter et al., 2014;Sulpizio and Massaro, 2017). Large-magnitude tectonically-generated earthquakes have also been associated with antipodal seismic focusing (O'Malley et al., 2018). There is clear observational evidence of focusing of seismic surface waves near Atitlán caldera from large earthquakes in Sumatra (e.g., Dec 2004, M9.1 and April 2012, M8.7; Supplementary Figure  S8) although this focusing is not exactly antipodal due to Earth's internal heterogeneity (e.g, subduction zones, mid-ocean ridges). An estimate for the total elastic energy released during the YTT supereruption is on the order of 10 19 [J] (Gudmundsson, 2014), which is the same order of magnitude as the largest instrumentally recorded earthquake, the M9.5 Valdivia (Chile) earthquake (Gudmundsson, 2016). The YTT seismic estimate should be considered an integrated estimate since the total energy released from the caldera collapse could have occurred over hours to days rather than an instantaneous rupture event (Stix and Kobayashi, 2008). However, it is noteworthy that there is potential for constructive wave interferences from a long-lived event (or multiple events) especially close to an antipodal location that may enhance the ground motion and local deformation (Supplementary Figure S8). However, the potential causal relationship between seismic energy released and initiation of volcanic eruptions hundreds of kilometers away remains poorly constrained (Hill et al., 2002;Seropian et al., 2021). It has been documented for only 0.4% of historical eruptions, though this probability may increase to 10% when considering a 2-year window between a leading large earthquake and a subsequent explosive eruption (Sawi and Manga, 2018). Causal effects are further supported by a temporal link between large magnitude earthquakes and volcanic activity at a global scale as proposed for the M9.1 Sumatra earthquake FIGURE 2 | Cumulative number of eruptions (>400 km 3 bulk tephra volume) through the last ca. 2.2 Myr (LaMEVE). The color and size of the symbols represent the magnitude of the eruptions. Black bars through symbols are 1σ age uncertainties. The upper inset plot shows the probability density curves for the coefficient of variation (CV). Values of CV > 1 indicate clustering of eruptions and CV < 1 periodic eruptions. The lower inset histogram shows the number of eruption histories (among 50,000 synthetic eruptions histories assuming that eruptions are randomly distributed) that contain paired eruptions within 80-400 years and at least one supereruption (>1,000 km 3 ) divided by the number of eruption pairs. A paired supereruption with YTT-LCY characteristics would be represented by "1.5" or "2.0" (either 1 or 2 eruption pairs) on the x-axis. On the other hand, if only one of the two closely spaced eruptions is a supereruption, it would be represented by the "0.5" or "1" (either 1 or 2 eruption pairs) bin in the x-axis. The numbers on each histogram show the percentage probability of being in that bin based on the synthetic eruptive histories.
Frontiers in Earth Science | www.frontiersin.org April 2022 | Volume 10 | Article 827252 5 (Hill-Butler et al., 2020). The peak ground motion and dynamic stresses induced by passing seismic waves (particularly long-period surface waves) have been linked to the onset of multiple magmatic processes, affecting the host-rock, magma chamber, or associated hydrothermal system that could ultimately lead to an eruption (Davis et al., 2007;Seropian et al., 2021).
Large supereruption-feeding magma systems are thought to be buffered and "perched" at a critical threshold for extended periods and then externally triggered (de Silva and Gregg, 2014;Wilson et al., 2021). Long residence in a melt-present buffered state for both the YTT and LCY supervolcanic magmatic systems is suggested by their protracted zircon crystallization records. Notably, both the LCY and YTT exhibit strikingly similar thermochemical histories based on the crystallization ages of zircon ( Figure 4A) which is sensitive to changes in magma chemistry and temperature. Magma accumulation timescales inferred from zircon ages are on the order of tens of thousands of years prior to both supereruptions, with a coincident maximum at ca. 96 ka. This suggests that the main phase of silicic magma differentiation and segregation of a melt-dominated magma body for YTT and LCY likely occurred within a similar time window of ca. 20 kyr before eruption. The maximum apparent rates of crystallization in magma bodies as indicated by zircon age maxima in the probability-density curve ( Figure 4A) have been shown to occur over a wide range of timescales prior to eruption in different volcanic systems and are still not well understood (4-247 kyr; Simon et al., 2008;Wilson and Charlier, 2009). A rough correlation has been observed between older pre-eruption zircon age maxima (tens of kyr) and magmatic systems showing lower temperatures and higher water contents (Tierney et al., 2019). The exact significance of coincident zircon age spectra and crystallization maxima in YTT and LCY magmatic systems requires further investigation, yet enhanced zircon crystallization is usually preceded by significant recharge events (Klemetti and Clynne, 2014). However, if a coincident recharge event in both systems occurred shortly before the ca. 96 ka crystallization maxima, excursions in zircon trace elements (e.g., elevated Ti-in-zircon) are expected, which is presently unsupported. Consequently, whereas a similar temporal evolution based on the zircon ages is suggestive of synchronous magmatic evolution towards a critical state, the available data do not permit a non-unique interpretation.
At present, the eruption order between YTT and LCY is unresolved given existing dating uncertainties. Nevertheless, we propose that the YTT preceding LCY is the more likely scenario since the larger YTT supereruption would have produced stronger seismic perturbations (Zürn and Widmer, 1996). The passage of long-period (tens of seconds) Rayleigh seismic waves generated during the YTT catastrophic eruption and caldera collapse through a crystal-mush-dominated LCY reservoir may have affected the system's stability ultimately culminating in a supereruption on a decadal-century timescale ( Figure 4B). At these seismic wave frequencies (10-20 s periods), the surface Rayleigh waves have peak depth sensitivities of 5-30 km (upper crust; O'Donnell et al., 2019) and thus are ideally suited to influence crustal magma reservoirs. We posit that the enhanced attenuation of the seismic surface waves in a FIGURE 3 | Large volcanic eruptions (>400 km 3 bulk tephra volume) from the LaMEVE database over the last 2.1 Myr showing vanishingly small chance of supereruptions occurring at antipodal locations. Eruption magnitudes (bulk tephra volume) are represented by the size and color of the symbol. The antipodal location for each eruption is shown as blue circles. It is noteworthy that the main potentially antipodal relationship between two supervolcano eruptions also close in time is the YTT and LCY eruption pair (~2,200 km between the Atitlán caldera and the antipodal location of the Toba caldera). Note that Toba has sourced two Quaternary supereruptions, with YTT represented by the larger circle (id. Taupo).
Frontiers in Earth Science | www.frontiersin.org April 2022 | Volume 10 | Article 827252 6 granular mushy LCY magma reservoir could have increased the local melt permeability due to mush dilation, enabling rapid assembly of an eruptive magma body on decadal timescales (Winkler and Nur, 1982;O'Donovan et al., 2016).
Alternatively, the surface wave attenuation could have enhanced local pore pressure leading to small-scale fracturing/ micro-seismicity, enhanced nucleation and growth of bubbles (Linde and Sacks, 1998;Hill et al., 2002;Nishimura, 2017), and overall promoting liquefaction (transformation of a solid where it behaves fluid dynamically as a liquid; Weinberg et al., 2021) in the crystalline mush (Müller et al., 2010;Sharma, 2012;Waymel et al., 2018;Zaccherini et al., 2020;Chapman et al., 2021;Wang et al., 2021;Yuan et al., 2021). Since mush domains behave as permeable materials, these spatio-temporal changes in melt and stress changes as well as crystal reorganization through decompaction, rotation/local shearing, or dissolutionreprecipitation of crystals lead to liquefaction and accumulation of eruptible magma (Léopoldès et al., 2013;Giacco et al., 2018;Brum et al., 2019;Dörfler and Kenkmann, 2020;Weinberg et al., 2021;Li et al., 2022). For the LCY system (or conversely YTT), the timescales necessary for the magmatic system to reach a critically perched state leading to a supereruption will depend on several parameters including crustal thermal state, magma reservoir(s) size and spatial distribution, volatile content, crustal stress state, and material properties in addition to magma composition. Nevertheless, since permeable melt flow is a slow process (compared to flow in open fractures), a lag in triggering a supereruption from another is to be expected. Thus, a decadal to a centennial time delay between the YTT and LCY eruption pair is physically plausible ( Figure 4B). The reverse triggering of YTT by LCY would also be possible through similar physical effects of LCY associated long-period seismic waves.
Assuming that the seismic energy localization effect for the antipodal magma reservoir is valid, we posit that specific requirements still need to be fulfilled in order to trigger a (super)eruption doublet, including a sufficient number of magma bodies at a critical state and a sufficiently large magnitude of the sequence of volcanic-seismic events itself. This represents the "advance-clock" hypothesis wherein the seismic energy helped synchronize two eruptions that may still have occurred independently but at different times (Gomberg et al., 1997;Sawi and Manga, 2018). Failing to fulfill any of the above-mentioned requirements will preclude a correlated eruption. This scenario can also potentially explain why known supereruptions lack corresponding antipodal-triggered eruptions in the LaMEVE catalog (except YTT-LCY). Since there is significant uncertainty about the physical processes associated with the seismic source during a large caldera collapse eruption (e.g., instantaneous versus prolonged duration; efficiency of body versus surface wave generation), a detailed wave propagation model is beyond the scope of the present work. We can, however, outline one consequence of our proposal of supereruptions triggering other eruptions, namely that the YTT or LCY event may have also primed smaller volcanic systems, especially close to their respective locations. However, these less-significant eruptions are likely poorly preserved in the geologic record and/or remained unstudied. An exception could be the Arce tephra erupted from Coatepeque caldera in El Salvador, which produced two large silicic eruptions (~26 and 41 km 3 tephra volume) separated only by a couple of hundreds of years (Kutterolf et al., 2019) with an age of 72 ± 2 ka (Rose et al., 1999) that overlaps that of the YTT and LCY.

Understanding the Climate Impacts of YTT and LCY Supereruptions
Resolving whether the time-space relationship between the YTT and LCY was not purely random but influenced by external factors would critically benefit from refining the absolute dating We propose that these processes may promote new migration pathways or induce instabilities in the magma reservoir and crust ultimately leading to eruption.
for both supereruptions (and other close supereruption pairs), preferentially by applying the same geochronological method, and from high-fidelity stratigraphic records where both eruptions could have deposited tephra. This also holds for assessing the climatic consequences of such paired supereruptions (Brenna et al., 2020;Black et al., 2021;Osipov et al., 2021) connoting that in combination both supereruptions would potentially have more impact on global climate than each eruption on its own (see further discussion in Supplementary Material). We think, however, that ultimately determining the precise time-lapse between the YTT and LCY would come from the identification of volcanic glass shards from both supereruptions within the ice-core layers or high-resolution paleoproxy records that could permit high-precision dating. We deem such an endeavor promising because glass compositions from the YTT and LCY tephra are unambiguously distinct in trace element abundances (Supplementary Figure S9).

FUNDING
AC acknowledges funding support from the Deutsche Forschungsgemeinschaft (DFG) grants SCHM 2521/6-1 and KU 2685/7-1; TM acknowledges funding support from the Crosby Postdoc Fellowship at MIT; SdeS acknowledges support from National Science Foundation grant EAR 1551187 for studies of the Toba caldera.

ACKNOWLEDGMENTS
We thank editor Valerio Acocella for handling this manuscript. The differing opinions and perspectives of the reviewers have been crucial in challenging our thinking and for this, and we are very appreciative of journal reviews RC, TW and John Stix; not all who agree with our speculations herein. Nick Petford and Darren Mark are also thanked for their comments on an earlier version of this work.