Evidence of megathrust earthquakes and seismic supercycles in subtropical Japan from millennia-old coral microatolls

Domain: neuroscience
Published: 2026-02-10
Editorial outcome: major_revision

Abstract
Megathrust earthquakes in subduction zones often go unreported because they are rare and the historical record is short. On the Ryukyu subduction zone of southwestern Japan, unlike neighboring Nankai Trough, the history and future potential of great interplate earthquakes are not well known. While the geodetic measurements on the islands suggest that the plate coupling is very weak, recent observations of slow seismic events as well as offshore geodetic measurements imply the presence of coupled patches along the megathrust. Furthermore, the historical and geological studies indicate evidence of great tsunamis. Here, we use fossil microatolls in Ishigaki island to reconstruct the relative sea level in the Holocene. The coral record reveals several relative emergence episodes clustering between 5-4 and 3-2 thousand years ago (ka). Elastic modeling shows that the observed motions can correspond to coseismic uplift associated with megathrust earthquakes. The clusters of megathrust events suggest possible supercycles of earthquakes with a recurrence interval of more than 2 ka. Such results imply a strong seismic hazard for the upcoming centuries. The devastating 1771 Meiwa earthquake and associated tsunami may mark the onset of the most recent seismic supercycle.

Introduction
The seismic behavior of the Ryukyu subduction zone in Japan (Fig. 1) is a matter of debate. Geodetic data in the last decades tend to indicate that the plate interface is weak (coupling rate lower than 5%), and that the convergence is accommodated by aseismic creep1,2 or slow slip events3,4, large megathrust earthquakes being unlikely5,6. However, two major earthquakes occurred on April 24, 1771 (Meiwa) in the southern Ryukyus and on June 15th, 1911 offshore Kikai island (Fig. 1), although their source remains debated. Both events generated large tsunamis with wave run-ups reaching tens of meters for the Meiwa tsunami7,8 and up to 5 m for the Kikai event9. In addition, geological records provided evidence for prehistorical large tsunamis in the southern Ryukyus7,8,10–14 and for a sudden coseismic coastal uplift at ~0 B.C. in Ishigaki island15,16. Moreover, interseismic strain accumulation on patches of the subduction interface has been inferred for the last one to two centuries from modern microatoll coral analysis17,18, and over the past decades from gravimetric or offshore geodetic data19,20. Aseismic coupling of 20% has been estimated for those patches12.

In this work, we aim to better understand the seismogenic behavior of the Ryukyu subduction zone and to bridge the gap between short-term geodetic series and paleo‑tsunami record. We focus on an exceptionally long and unique record of centimetric-scale sea-level changes from fossil microatolls in Ishigaki island. We found, measured, sampled and dated several generations of fossil and living microatolls in Nagura site. We show that the development of these corals is controlled by seismic supercycles and permanent deformation. This study provides an example of multiple earthquake cycles over the past five millennia in a subduction zone, based on coral microatolls analysis. It challenges previous inference from inland geodesy that this subduction zone is mainly aseismic.

Results
Nagura microatolls measurement and sampling
The Nagura site extends over a 3.4 km2 area, 130 km far from the Ryukyu trench in the southwestern part of Ishigaki island (Fig. 1). We found multiple generations of in-situ Porites fossil corals within a mangrove on the large alluvial fan of the Nagura river, and on a sandy tidal flat forming a bay open to the ocean (Fig. 2a, b). Microatolls have been previously documented around Ishigaki21 but to our knowledge, Nagura is the only site where fossils specimens have been observed22.

Microatolls are massive corals that grow radially from a single polyp in a tropical intertidal zone. Once the coral reaches its highest level of survival (HLS), its uppermost part dies. The lowermost part keeps growing below the HLS, thus recording annual relative sea-level (RSL) changes with a centimetric precision (see Method).

Based on field observations and using high-resolution drone imagery, a total station and real-time kinematic Global Navigation Satellite System (RTK-GNSS), we measured the position and elevation of 73 fossil microatolls, 1–11 m wide in diameter, up to 40 cm thick, as well as a living coral located 1.2 km offshore (Fig. 2c). Using a hydraulic chainsaw, we sampled seven of the largest and best-preserved fossil microatolls by extracting a ~10-cm-thick slab from which 1–2-cm-thick slices were X-rayed, CT-scanned and interpreted. In addition, we sampled the outer part of 12 specimens from various generations of microatolls over the entire site for U/Th dating23 (see Method and Supplementary Tables S1, S2). Photos of the sites, slabs or samples of each microatoll, X-ray interpretation and HLS curve are presented in Supplementary Information S1, S2.

Relative sea-level changes over the past 5 ka
The upper surfaces of the microatolls, whose altitudes have been measured with a precision of ±5 to ±7 cm, define 23 plateaus that represent topographic features (Fig. 2c and Supplementary Table S3). They are first order-paleo sea levels, with the highest—and generally landward—being the oldest. The overall increase in elevation and age landward supports the interpretation that the microatolls are in-situ. Two deeper fossil microatolls (Yama_1 and Yama_2) measured at −1.80 m and −2.00 m and radiocarbon-dated at A.D. 790±97 and A.D. 1157±92 by Yamaguchi22 are included in our analysis (see Method).

Individual coral slabs record a complex RSL history, alternating between 1) periods of free growth below the HLS, 2) periods of growth close to or at the HLS, and 3) sudden drops or increases of the RSL (see Method and Supplementary Fig. S12). The drops were large enough to partially or entirely kill the corals. The microatolls sliced in the mangrove recorded seven sudden RSL drops with values ranging from 11 to 40 cm over approximately three centuries between 2500 and 2000 B.C. (Fig. 3a). Coeval corals (within age uncertainties) recorded the same HLS variations, indicating that the record is reliable. Microatolls in the bay recorded five sudden RSL drops, with values between 9 and 46 cm over five centuries between 600 and 0 B.C. (Fig. 3b).

These RSL drops are followed by years of HLS stability (with RSL changes at small rates of less than 10 mm/yr, Supplementary Table S4), overall indicating that they are of tectonic origin and due to earthquakes24. Long lasting RSL changes at small rates are promoted by interseismic tectonic motion (long term strain accumulation or transient displacements18,25). They likely also include a climatic component (eustatic sea-level changes or glacial isostatic adjustment, GIA) for which we have no precise information. Here, RSL contributions from GIA was inferred to be 0–0.2 mm/year for the last 4–5 ka26,27, which is one order of magnitude smaller.

Although NAG3_21E is coeval with NAG_B, it is 30–40 cm lower and did not record RSL drops of similar magnitude. Variation of tens of centimeters in altitude among corals of the same generation have been observed elsewhere, reflecting the variability of the coral record28. Ponding of NAG_B coral could explain such differences, although we did not observe any traces of past geomorphic features supporting this hypothesis.

Seismic supercycles
Overall, the RSL changes display an irregular and truncated sawtooth pattern of interseismic subsidence or uplift interrupted by sudden coseismic uplifts (Supplementary Table S5). A long-term centennial-scale RSL decrease rate of about 1.0 mm/yr is superimposed on this sawtooth pattern (light thick orange line in Fig. 3). At the mangrove site (Fig. 3a), NAG6_21C and NAG6_18Q are older than NAG_M and NAG6_18H yet lie below them. This may indicate a sudden sea-level increase of ~35 cm.

Altitudes and dating of the external parts of other microatolls (dated or undated plateau, squares, Fig. 4, Supplementary Table S2 for ages and Supplementary Information S1 for information on samples) complete the sawtooth record. We used an interseismic RSL rate of 2.0±0.6 mm/yr (as determined from the best-preserved signal from NAG3_18H, Supplementary Table S4) to estimate the magnitude of the events that killed the microatolls. We consider our approach is reasonable, although variations in the rate of interseismic sea-level changes were observed at smaller time scale in the RSL records of the slices and documented in the literature29 (Supplementary Table S3). We inferred that the ages of the undated plateaus are bounded by those of the plateaus just higher or lower in elevation (Supplementary Table S3).

In addition to the 12 major earthquakes documented from the slices (Fig. 3), at least nine additional sea-level drops with amplitudes ranging from 15±43 to 262±75 cm (events named with letters A to U, Supplementary Table S5) are required to explain our observations. It is striking that over a time period exceeding 1300 years, between 2000 and 700 B.C., there are almost no indication of the RSL. The RSL record is also missing between 75±6 B.C. (death of NAG3_21E) and A.D. 790±97 (Yama_1). This could be due to coral absence or scarcity. Only few fossil microatolls (most undated) suggest the occurrence of earthquakes during these periods (event L and dotted portion in Fig. 4). It is possible that these long periods of missing corals followed major earthquakes, i.e. with greater uplift and a longer recovery time for the reef flat, as observed in Sumatra30. From our reconstruction, we estimate values of 136±41 and up to 262±75 cm for the main events K and R (Supplementary Table S5).

The most recent part of the sawtooth signal includes an undated microatoll marking the lowest plateau P21 (Supplementary Table S2). As its altitude is lower than NAG3_21E, it is probably more recent and may have developed between 75±6 B.C. and the present day. If it is younger than A.D. 1157±92 (age of Yama_2), this coral may have been killed during the Meiwa event in 1771.

The coseismic motions we infer have an average amplitude of 59 cm. At a similar distance to the trench, this amplitude is comparable to that modeled from past Mw-8 class subduction earthquakes in Sumatra (e.g. between 0 to 1 m31). We verified with elastic models (Supplementary Information S3) that ruptures on either shallow normal faults on Okinawa Trough or on splay faults in the forearc domain cannot account for the observed deformation. The sea-level changes recorded are most probably due to strain accumulation and release during the seismic cycle on the Ryukyu megathrust.

At first order, the sawtooth pattern characterizing the RSL signal appears to define supercycles. It shows periods of failure sequence with frequent and discrete emergence events (C-J and M-Q, Fig. 4) before large events (K, R). Those periods are followed by periods of quiescence (characterized by coral absence or scarcity). Supercycle 1 lasted at least ~900 years, contains at least 11 earthquakes with amplitudes ranging from 7 to 136 cm and a mean recurrence time of ~88 years during its failure sequence. It ended by the large 1970±65 B.C. K event. We cannot confirm that this cycle is complete but we have not found any coral older than 2900 B.C. Supercycle 2 lasted about 2000 years and contains at least 6 seismic events of 15–262 cm of amplitude. It ended with the large 75±6 B.C. R event. During the failure sequence, the mean recurrence time of earthquakes is ~110 years. The few earthquakes that occurred during the more quiescent periods at the start of supercycles 2 and 3 have a recurrence time of ~671 and 615 years respectively. This behavior was also observed in Sumatra where the recurrence time of earthquakes within supercycles is ~150 years, and of several hundred years during periods of quiescence following a major earthquake32,33. In this study, supercycles are characterized based on a single paeloseismic site only. To improve their description, it would be necessary to assess the spatial extent of the seismic ruptures documented here33, which is beyond the scope of this paper.

Source of RSL: elastic models
We used dislocation in a 2D elastic half-space to model the vertical tectonic deformations34–36 (see Method, Fig. 5). We adopted the geometry of the SLAB 2.0 model37, which is consistent with that derived from seismic multichannel profiles offshore southern Ryukyus5,6,38 (Fig. 1b). In Fig. 5, RSL drops could be promoted by coseismic uplifts due to a uniform slip of 8–25 m along faults 1+2, or by 1–9 m between 33 and 60 km on fault 2. By considering an arbitrary yet reasonable and conservative length of ~300 km for the ruptured patch (between Yonaguni and Miyako islands39, Figs. 1 and 5), the formulation from Aki40 and Hanks and Kanamori41 with a rigidity of 32 GPa results in magnitude between Mw 7.7 and 8.3 if the rupture is restricted to the deepest segment (fault 2), and between Mw 8.3 and 9.1 if the earthquake ruptured the entire plate interface (faults 1+2). The hinge line (0 deformation), indicating the downdip limit of the rupture, would be located somewhere between Ishigaki island and the edge of the Okinawa Trough (Fig. 5). We also modeled the ~0 B.C. coseismic uplifts values documented by Kawana15,16 and found that they could be reproduced by 25 m of slip on the entire plate interface between 0 and 50 km depth. A rupture of 300 km long would imply a magnitude Mw 9.0 for this earthquake.

In Fig. 4, some old corals are lower than younger generations. Our models suggest that this could be due to a sudden coseismic subsidence generated by a shallow megathrust rupture or a transient interseismic deformation, possibly SSE18 (Supplementary Information S3). We cannot, however, exclude the possibility of local coseismic compaction of the soil.

Interseismic surface deformation varies over time, as also observed in Sumatra29. This can be explained by the occurrence of transient events and loading of more or less coupled asperities at various depths42. From modeling, we estimated the coupling rate of asperities of ~15–100% (Supplementary Information S3). The lower bound is comparable to the 20% estimated by Ando et al.12. Slow interseismic uplift or subsidence can be explained by alternating loading of the shallowest and deepest segments, respectively (Supplementary Information S3).

Methods
Sampling site study
We measured and sampled numerous sets of well-preserved fossil and living microatolls over an extended area in the Nagura bay (Fig. 2). Fossil microatolls provide highly useful and valuable information to retrieve past RSL variations due to tectonics or climate30,53–55. They were first characterized by Stoddart and Scoffin56, as they present the same morphology as atolls, typically a round shape at the surface with a dead plateau at the center. Microatolls grow at a first-order constant rate underwater and their development is controlled by RSL variations at a scale of centimeters. Their upward growth is limited by their highest level of survival, or “HLS”. The HLS is the theoretical elevation above which the exposed part of the coral dies due to excessive sub-aerial exposure. Once the microatoll records its first HLS, the growth of its surviving lower part is conditioned by the variations of the HLS. We distinguish HLS and HLG57–59 (for “highest level of growth”). “HLG” is used for years during which the coral grew freely upward and the HLS was not reached (green points in Supplementary Fig. S12); the distance to the theoretical coral HLS is unknown60. In the coral stratigraphy, a HLS year preserved from erosion corresponds to an unconformity over which younger growth bands curl over, while the growth band is complete during a HLG year57. The coral upward growth is often disrupted by temporary centimetric sea-level drops, commonly named “diedowns”24,59.

The microatoll morphology evolves accordingly and therefore records the variations of the RSL. When the RSL is stable, the microatoll develops laterally and its surface evolves into a flat plateau. When the RSL decreases, the upper part of the microatoll is emerged and die, and the part which lies below the new HLS continues to grow mostly laterally. The emergence can be sudden or progressive. In the case of sudden emergence, the morphology of the coral will show a clear lower step in its external part. In the case of progressive emergence, the coral height will decrease towards its most external part. In both cases of emergence, the coral has a hat-shape morphology. In the case of submergence, when the RSL increases, the coral grows upward and has a cup shape morphology. As the coral upward growth is limited by the coral growth rate, no distinction can be made between a sudden or a progressive submergence context. Alternation between periods of relative submergence and relative emergence along the coral life implies the formation of concentric rings on its top surface. Such morphological features can often be observed and measured.

The exact relationship between HLS and RSL is not fully understood today. The growth rate is 5–25 mm/yr for Porites species28. Each year, the coral microatoll records a dark and a bright growth band that correspond to variations in aragonite density related to seasonal changes. Each of the growth bands will provide information on the HLG/HLS of the coral for its corresponding year. In our study, we use bright growth bands which correspond to rainy season in the middle of the year. The HLS is often associated with the lowest tide level, or Extreme Lowest Water24,57,61 (ELW). The survival of a microatoll may also be associated with other phenomenon such as the sunlight exposure24. Thus, the HLS may lie tens of centimeters above the corresponding ELW17,18,60. This observation is very specific to each study site as the survival of the coral depends on the species and on the annual tidal characteristics of the sampling site.

The RSL changes deduced from the morphology of a microatoll can be of climatic, oceanographic and/or tectonic origin. During a period of stable absolute sea level, aseismic uplift would induce a RSL decrease. At the contrary, aseismic subsidence would correspond to a RSL increase. An event of tectonic origin is a few to few tens of centimeters of persistent vertical motion of the land in a restricted geographic area. It induces a RSL change and may be coeval with a historical earthquake24,58. Annual or multidecadal sea-level changes induced by climatic and/or oceanographic events, such as the El Niño Southern Oscillation17 are recorded over a wider area. In the coral stratigraphy, they are characterized by an abrupt decrease followed by a rapid increase of the RSL. The modern RSL record of microatolls can be compared with that of the tide gauge, satellite and GNSS data to decipher between the tectonic and climatic and/or oceanographic signals. At a longer time scale (Holocene period, few millennia), the microatolls may record sea-level changes related to the glacial isostatic adjustment (GIA) following the ice sheet melting in northern hemisphere or Antarctica54,55,62,63.

During field missions in July and December 2018, and November 2021, we mapped each set of microatolls by using a total station, a tool combining a theodolite for angle measurements with an integrated laser distance meter, allowing rapid acquisition of relative three‑dimensional coordinates. We also used a Real-Time Kinematic GNSS with WGS 84 coordinates (RTK-GNSS, Fig. 2b, c). We measured in a systematic way the altitude of the external ring (younger part of the microatolls, i.e. most recent HLG), the center of all microatolls and the elevation of the substrate on which they developed, when possible. When visible, we also measured the altitude of other morphological features such as markers of possible older HLS records (older concentric rings) along diameters in several directions. We measured the altitude of the HLS of one living microatoll (Fig. 2c), which was necessary to compare the elevation of the old HLS recorded by the fossil microatolls and the modern HLS.

Altitudes and distances relative to the total station base are measured with millimetric precision. We used reference screws and reference microatolls to adjust total station measurement to GNSS results with a closure error of <1 cm, and with a resulting precision related to the GPS-RTK method of ±5 cm. The vertical variability on a given concentric ring is estimated at ±3 cm in average. This last uncertainty is considered as a statistical type A uncertainty, σmean, and the error on a single concentric ring elevation (above the 0 of the WGS 84 ellipsoid), σElev, is the sum of the two errors added in quadrature, i.e. ±5 cm to ±7 cm.

σElev = sqrt(σmean^2 + 0.05^2) (1)

Other errors related to erosion or compaction of the site are complex and difficult to quantify. As we consider that they mostly affect all corals similarly, they are unlikely to bias the RSL curve. Reconstruction of the site. To include deeper fossils documented by Yamaguchi22 (Yama_1 and Yama_2, Fig. 2), we plotted the microatolls documented in Yamaguchi22 against our measurements. As their elevation was estimated and not measured with RTK-GNSS, we apply an uncertainty of ±50 cm. In Yamaguchi22, Yama_1 and Yama_2 were dated with 14C method, with ages of 1063-1257 cal BP and 700-885 cal BP, respectively. This corresponds to ages of 790±97 and 1157±92 years A.D., respectively. High-resolution of 3.2 cm at 75 m of altitude drone imagery was also acquired in the site of Nagura with a Promark120 drone. A mosaic of images was made using the Photoscan software.

Corals slices analysis
After the topographic mapping, we sampled, among the largest and best preserved, seven microatoll fossils. Before the sampling we verified that they were not tilted nor overturned. We sampled ~10-cm-thick sliced of coral skeleton along the radii or the diameters of the microatolls with a hydraulic chainsaw. Before sampling, we installed screws at ~10 centimeters interval along the slice and we measured their altitude in order to retrieve the growth position of the slice before sampling. A reference screw was also installed and measured near the microatoll to compare the measurements made during different days and for future investigations if needed. Thin cutting of the coral slabs was further performed at the Ogyu Manufacturing Inc., also known in Japan as Milano Seisakusho Co. Ltd (equipped with stable diamond wire) to obtain 1–2-cm-thick slices. The thin slices were then X-rayed at the Veterinary Medical Center of the University of Tokyo to image in thin section the internal architecture of the coral skeleton and to evidence the alternation of white and darker growth bands for counting. This is also important to distinguish and discard from our analysis new overgrowths that commonly develop on top of dead corals. We also used CT scan imagery available at the National Institute of Advanced Industrial Science and Technology in Tsukuba. CT scan advantage is to provide a 3D tomography of the slice in which numerical virtual slices (tomographic cross-sections) as thin as 2 mm can be extracted by using the Horos software (Horos Project; https://www.horosproject.org/) without cutting the sample. The resolution of those virtual sections is lower than X-Ray however. The X-Ray or CT scan images are reoriented within an uncertainty of ±1 cm in its initial growth position by using the altitude of the screws. We then measured and characterized each growth band to reconstruct the RSL change curve.

Absolute dating is needed to estimate the chronology of the sea-level changes. We collected and dated with the U/Th method23,64 several samples both along the coral slabs and on many other microatolls of different generations (Supplementary Table S2, Supplementary Information S1). Along the slabs, we collected, by pairs with a small core drill, several micro-core samples of the coral skeleton on one growth band. We used a drill to ensure that the hole created from the core drill extends through the whole slice and can be located later on X-ray imageries of the slice. On the other fossil corals not sampled with the chainsaw, we collected with a hammer decimetric samples on the youngest external part of the microatolls. Considering an average growth rate of ~1 cm/yr for Porites corals28, we estimate that their age approximates that of the most external growth band, i.e. the year of their death, by several years. Coral boulders samples were cleaned by ultrasonically bath using tap water for 20 min. After that they were cut within several cm thicknesses. The coral slices sampled were also cleaned by ultrasonically bath using tap water for 20 min and dried in the incubator at 40 °C for several days. All samples were divided into two subsamples, one for U/Th dating, the other one, grounded to powder, for X-Ray Diffraction (XRD) measurements. Those last samples were analyzed with a MiniFlex II (RIGAKU corporation) diffractometer. The scan-angle, speed, and step were from 25 to 32°, 1°/minute and 0.01°, respectively. The aragonite and calcite in each sample were detected using PDXL v.1.8.03 Software (RIGAKU corporation).

All samples were measured by XRD to estimate the percentage of aragonite and screen for calcite. Only samples with more than 95% of aragonite were dated at the High-Precision Mass Spectrometry and Environment Change Laboratory (HISPEC), National Taiwan University23. The results are presented in Supplementary Table S2.

Uncertainties on the chronology of the RSL changes may be large when established by using fossil corals31,53,65. It depends both on the uncertainties on U/Th ages and on the growth band counting. We assign an age and uncertainty to each growth band from a combination of U/Th dating results and growth band counting. To estimate uncertainty related to growth band counting, we used a numerical estimation from Weil-Accardo66. We followed the formula detailed by Meltzner and Woodroffe28.

We reconstructed the RSL curve and identified diedowns using the difference in elevation between the highest HLG point (eroded or not) preceding the diedown60, and of the HLS point (eroded or not) marking the diedown itself. We calculated the rate of sea-level rise of decrease from our microatoll record by using different methods provided in Zachariasen58 and Zachariasen et al.57 (Supplementary Table S4).

The methods used to estimate RSL trends are all based on simple linear regression and depends on the number and type of data used: Zachariasen et al.57 used the maximum height of all growth bands except those corresponding to the initial upward growth of the coral (Zach-1 method), Zachariasen58 used only the elevation from the HLS bands (Zach-2 method), and we propose Zach-3 method, derived from Zach-2, that consider only non-eroded HLS bands. Following Meltzner et al.60, we also use when available the highest eroded or preserved HLG point preceding a diedown (Meltzner method). We distinguish the Zach-3 and Meltzner methods that provide well-constrained results. In contrast, the Zach-1 and Zach-2 methods include HLG points and/or eroded HLS/HLG points, which can lead to an underestimation of the RSL and produce results that do not reflect the true RSL changes. We refer to these results as “apparent” RSL changes. However, when the RSL record is robust, all method gives similar results67,68 (Supplementary Table S3). We use Zach-3 method to reconstruct the RSL history in Nagura and Zach-1 method for slice analysis, when Zach-3 method could not be applied. The HLS curve of a single microatoll has a relative uncertainty of ±1 cm, due to the uncertainty in the horizontal reorientation of X-Ray and CT scan images. The amplitude of the RSL drops that were only measured in one slice has therefore an uncertainty of ±1 cm. When a drop was observed in several slices (e.g. NAG6_21C and NAG6_18Q, Fig. 3), we observe an average range of 15 cm between measurements, leading to an estimated uncertainty of ±7 cm. This value is in the same range as the HLS variability observed in a given site from microatoll study in the Ryukyus17,18, therefore we use an HLS variability of ±7 cm. As an indicator of past RSL, each plateau is assigned an uncertainty of ±7 cm. Each plateau elevation is estimated as the average elevation of 1–10 corals or concentric rings in coral surface. For the uncertainty in elevation of each coral σHLS, we added in quadrature the uncertainty of ±7 cm related to HLS variability.

σHLS = sqrt(σElev^2 + 0.07^2) (2)

To estimate the uncertainty of each plateau σP, we calculate the uncertainty on the weighted mean using the following equation:

σP = 1 / sqrt(sum(1/σHLSi^2)) (3)

The amplitude of the RSL drops we infer from a given plateau Pi is therefore controlled by (1) the age uncertainties of the plateau Pi and the following younger one Pi-1, (2) the vertical uncertainty of both plateaus σPi and σPi-1 and by (3) the uncertainty of the interseismic RSL rate σR we use. Considering those latter, we calculated for each documented drop its related uncertainty with the following equation:

σDrop = sqrt(σPi^2 + σPi-1^2 + σinter^2) (4)

Estimation on the uncertainty related to the interseismic rate R+σR and the age of each plateau AgePi +σagePi and AgePi-1 ±σAgePi-1 is given by:

σinter = R * ΔT * sqrt((σR/R)^2 + (σagePi^2 + σagePi-1^2) / ΔT^2) (5)

With ΔT being the time lapse between the age of the two plateaus.

The diversity of our dataset affects the robustness of our interpretation. To account for this variability, we assign a confidence index ranging from 1 to 6 to each interpreted seismic event (Supplementary Table S5). Earthquakes identified within coral slices receive a confidence index of 6 when affected growth bands are preserved from erosion, and of 5 otherwise. Those interpreted based on the death of a coral observed at the same elevation as several other colonies are assigned a confidence index of 4. Events inferred from the death of a single coral colony are given a score of 3. Earthquakes corresponding to the observations reported by Yamaguchi22 receive an index of 2, and the event attributed to the hypothesized Meiwa earthquake has a score of 1. For all events whose estimated vertical motions are close or smaller than their associated uncertainties, the initial confidence index is reduced by one point. The same stands for events identified in one slab that are not observed in coeval slices. We deduced the subsidence event A’ from the difference in elevation of several corals (NAG6_21C and NAG6_18Q, ang NAG_M and NAG6_18H). It was assigned a confidence index of 4. However, we do not interpret it further due to the uncertainty of our 35 cm amplitude estimation.

Modeling
Our elastic model consists of a simple dislocation in a 2D elastic half-space34 with a rigidity of 32 Gpa. The geometry of the plate interface is determined from SLAB 2.0 model from Hayes et al.37 (Fig. 5).

According to SLAB 2.0, for simplicity, we modeled coseismic deformation on the plate interface with two faults: Fault 1 with a dip of 12.5° between 0 to 33 km and Fault 2, with a dip of 31° between 33 and 50 km or between 33 and 60 km. We also attempted to reproduce the coseismic deformation promoted by slip on a normal fault at the southern margin of the Okinawa trough. In Supplementary Fig. S13, we used a 20 km-wide dislocation dipping northward by 45° with a slip of 3 m (implying a magnitude 7.7 for a length of 100 km).

During the interseismic period, the slip deficit accumulation on the megathrust can be modeled by a back-slip on a planar fault35. Such a model has been developed on curved geometry, and when a kink in the slab dip is observed above the downdip limit of the locked zone, surface deformation can be better reproduced by dislocation on a planar fault tangential to this downdip limit36,69. In Supplementary Fig. S14, we therefore used back-slip dislocation on Fault 3 tangential at a depth of 33 km to the locked zone corresponding to Fault 1, and on Fault 4, which is tangential at a depth of 60 km to Fault 2. We used different dislocations on Faults 1–4, corresponding to seismic coupling between 0 and 100%, 100% indicating a fully-locked zone. When the plate interface is fully coupled, the slip deficit correspond to 12 cm/yr, which is the convergence rate between both plates70.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "Debaecker et al. present relative sea-level (RSL) reconstructions from coral microatolls located along the southeastern coast of Ishigaki Island, located north of the Ryukyu subduction zone. Combining the coral slab interpretations with elastic dislocation modelling, the authors suggest several coseismic uplift events superimposed on millennial-scale long-term relative sea-level (RSL) fall that may be attributed to permanent deformation (irrecoverable uplift) over seismic cycles. Considerable effort was taken to provide detailed documentation to support the findings of this study, both in the main text and the supplementary material. We believe that the coral record presented has important implications for advancing the understanding of the seismogenic potential of the southern Ryukyu subduction zone. We recommend publication of the manuscript after addressing the following concerns. We find that there are important details that need to be addressed, particularly with regards to the inferred coseismic uplifts. We believe that greater focus should be placed on the periods of RSL rise recorded internally within the coral slabs, which may indicate interseismic or postseismic subsidence. The uplift estimates are problematic due to erosion of the coral surface. Glacial isostatic adjustment (GIA) drives far-field RSL fall during the mid to late Holocene, which may explain a non-negligible part of the RSL fall that is currently inferred as long-term net land uplift. The inference of RSL changes between different coral generations (vs internally within a slab) is more challenging as it requires 1) clearer demonstration that the corals were in-situ and were not moved by past tsunamis; and 2) consideration of the natural variability in coral growth that is currently underestimated. We provide detailed comments below. Major comments 1. There seems to be a misinterpretation of diedowns, highest level of survival and highest level of growth, and missed recognition of erosion. Lines 406-407: “... a HLS year corresponds to a truncated growth band, while the growth band is complete during a HLG year”. A truncated growth band indicates erosion, not the HLS. The HLS is the level the coral dies down to during a diedown event. In cross section, the HLS can be recognised as an unconformity, where the annual bands are seen to curl over at the diedown (resulting from outward growth at the HLS after the diedown). In Text S2 of the supplemental information, years that were described as the ‘first HLS’ being reached were marked by the first (oldest) annual bands that were truncated (supplemental lines 143-144; 146). The first HLS should be where the first diedown is recorded by the coral, not where the first annual bands get truncated. In lines 406-409, HLS was described first, followed by a separate explanation of diedowns as ‘centimetric sea-level drops’, suggesting the authors do not recognise that HLS are only formed during diedowns. As a result of the misunderstanding of HLS, several consecutive years were misinterpreted as HLS years and therefore stable RSL (Figure S2b) (supplemental material lines 119; 150; 188; 199), when in fact we cannot know because of erosion of the top of the coral. Additionally, in several instances erosion was misinterpreted as relative emergence/ periods of RSL fall (supplemental material lines 123; 200; 213; 237). In these cases, often the slab profiles indicate the truncation of near-vertical growth bands, which are a clear sign of erosion. Near the base of a coral, the growth direction is near horizontal, so the annual bands (which are perpendicular to the growth direction) are near vertical. If there were no erosion, the surface of each concentric ring of the microatoll should show annual bands that are near-horizontal, as growth direction is near vertical at the ring crest. Figure S2b: • many ‘eroded HLS’ points are actually ‘eroded HLG’ • Some of the ‘HLG’ points are eroded, but there is no distinction between eroded HLG and preserved HLG in this figure. HLG is preserved when the ring crest of a particular annual band is preserved; if the annual band is truncated, then the HLG in that year is eroded. 2. The authors do not account for natural variability in the coral HLG/ HLS. For the RSL reconstruction, only two living microatolls were surveyed (line 343). This likely does not capture the natural variability in microatoll HLG. The authors report levelling uncertainties (line 348) and a vertical uncertainty of +/- 3 cm on a single concentric ring (lines 348-349), but it is unclear what the latter refers to. Published studies indicate that HLG of Porites sp. microatolls may vary by up to about +/- 20 cm across a study site. When reconstructing centimetric-scale RSL changes between different corals this uncertainty in the natural variability of coral HLG becomes important. Related to the above, the authors also do not mention the genus of corals being studied. Different genus can survive at different elevations with respect to the tides. Can this explain the vertical offset between the corals NAG_B and NAG3_18H? Based on the caption of Figure 4.2 in the Yamaguchi et al. 2016 chapter, it seems that there were fossil Porites sp. and Favia microatolls at this site. 3. Possibility of reworking This is an area with historical tsunamis (Figure 4b). Is there a possibility that the corals could have been transported during a tsunami? NAG6_21C looks like a fragment of a coral that was broken off and was found within the main river channel. Is it possible it was transported during the Yaeyama tsunami event? Were there dated fossil corals at similar elevations to those slabbed, which would support the in-situ nature of those slabbed corals? Currently, it is difficult to discern based on the plateaus shown in Figure 3. The plateaus are presumably combined from different dated (but not slabbed) microatolls? It may be more useful to show individual sea-level index points for each dated microatoll. As coseismic uplift events were inferred from the relative elevations of dated coral blocks, we believe it is also important to consider if these dated coral blocks were in-situ. 4. Coseismic uplift events The authors reconstructed coseismic uplift events by assuming an “arbitrary interseismic RSL rate of 1.4 +/- 0.2 mm/yr” (line 139), supplemented by coral analyses (Table S2b). The application of this arbitrary interseismic RSL increase (i.e., interseismic subsidence) forces coseismic events even when the coral elevations and ages themselves may not require so. Events L, M and N were recorded by dated plateaus at rather similar elevations. However, these were inferred as coseismic uplifts because of the assumed interseismic RSL rate. It is possible that this period was marked by gradual RSL fall, rather than coseismic uplift followed by interseismic subsidence. Of the 22 uplift events reported, only three were supported by RSL fall internally within the slab profiles (events E, G, and O). A third was not inferred but is evident in the NAG6_18H slab (the annual band for year 2089 that was sampled for the U/Th date A2); event J was inferred a few years later than year 2089 but the annual bands do not suggest a diedown in the year of event J. For events E, G and O, the slab profiles suggest evidence of a RSL fall, as indicated by a step-down in elevation from the pre-diedown HLG (be it eroded or preserved) to the HLS at the diedown. However, the HLS in events G and at the 2089 band in NAG6_18H were not preserved, so it is also possible that the step down in elevations are from erosion rather than a true RSL fall. The pre-diedown HLG in all of events E, G and O are not preserved, so the inferred uplifts can only be minimum uplift estimates. The diedowns for events O and G also appear from the slab annotations to be a few years earlier than is currently inferred. Event D is obscure because of the “non-uniform growth” that was described in lines 184-185 of the supplementary material. If the diedown was indeed at the top part of the discontinuous growth, then the diedown is small and may be due to other oceanographic/climatic drivers rather than coseismic uplift. However, the annual bands in the lower portion of this discontinuous growth appear to be curling in a different direction than the upper part. Is it possible that the upper part is out-of-sequence overgrowth? For event P, the pre-diedown HLG is also eroded, therefore the uplift for this event should also only be treated as a minimum uplift estimate. There appears to have been continuous upward growth between events O and P. It is possible that postseismic or interseismic subsidence after the uplift in event O outpaced the coral growth rate, and the coral only caught up to sea level during event P. Event P may also be due to other oceanographic/climatic drivers rather than coseismic uplift, as the magnitude of the diedown appears small. For event R, the HLG/HLS is not preserved, so it is not clear there was a diedown (or RSL fall), and it should not be used to quantify or characterise coseismic uplift. Based on the slab, there seems to be two to three earlier diedowns before event R, and one after event R – all of which seem to be separated by a similar number of years. This suggests that if event R were a diedown, it is likely to be caused by periodic oceanographic drivers, rather than coseismic uplift. Six uplift events were inferred from the thickness of the outer margin of the coral, which was assumed to have resulted from an uplift event that killed the coral entirely (events F, H, K, Q, R for NAG_B, S). It is common that coral microatolls produce RSL records that are temporally discontinuous, and which cluster in certain time periods. Therefore, we are not certain that the death of a microatoll necessarily implies that it was killed by an uplift event. Several uplift events (Table S2b) were inferred from the relative elevations of dated coral blocks (i.e. sudden ‘drops’ in the elevations of younger corals): events A, B, I, J, K, L, R, T, U. We suggest further investigation into the possibility that the lower corals were eroded. 5. Attribution to RSL drivers Overall, the elevations and ages of the corals in figure 4 appears to indicate a gradual RSL fall between ~3000 B.C. and present, apart from the two low coral blocks marking events T and U. Can the magnitude of the long-term RSL fall be explained by glacial isostatic adjustment (GIA)? Given that the study site is in the far field, a mid- to late Holocene RSL fall from a highstand can be expected. We suggest showing a RSL predictions from GIA models to illustrate the magnitude of RSL fall that can be expected from GIA. Within this record, the slab profiles of NAG6_18Q and NAG3_18H indicate clear and abrupt RSL fall followed by gradual RSL rise; and that of NAG_B indicates gradual RSL rise. We believe that these periods of postseismic/interseismic subsidence (RSL rise) are more robust than the coseismic uplifts, whose magnitudes are not well-constrained due to erosion of the pre-diedown HLG. Minor comments Line 28: Remove hyphen from “relative sea-levels”. As a noun, sea level should not be hyphenated. Lines 30-31: “Elastic modelling shows that the observed motions are coseismic uplift associated with megathrust earthquakes.” Rephrase this to suggest that there are there are other possible interpretations. In this study, a constant interseismic rate was assumed between majority of uplift events, but this is only one possible mechanism. Line 46: “A seismic coupling of 20% was inferred in this area.” Clarify what “this area” refers to. Line 52: “This study provides the first example of multiple earthquake cycles over five millennia in a subduction zone.” Is this referring to just the Ryukyu subduction zone? Or anywhere in the world? Line 71: “… its upper part dies and it keeps growing laterally” Rephrase this. After the coral dies during a diedown, the coral polyps grows out radially (both upwards and outwards), not just laterally. Line 80: “Photos of the sites, slabs, or samples of each microatolls.” Correct “each microatolls” to “each microatoll”. Line 93: “… measured with an accuracy of +/- 6 cm” Is this referring to the precision of the instrument, rather than accuracy? Line 125-126: “…, it grew 30 cm below and was too deep to record RSL drops of the same magnitude.” Even if so, NAG_B should have died down to a similar elevation (within the natural variability of HLS) as NAG3_21E, just that the magnitude of the coeval diedown in NAG_B is larger. In figure 3, it appears that NAG_B died down to an elevation (-0.3m) that was 40 cm higher than the elevation that NAG3_21E died down to (-0.7m) at ~256 B.C. Line 137: “… the external parts other microatolls”. Add “of” between “external” and “parts”. Line 139: “ We used an arbitrary interseismic RSL rate of 1.4 +/- 0.2 mm/yr (as determind from the best-preserved signal from NAG_B, Table S2a)…” In figure 3b, the interseismic rate labelled for NAG_B is 1.4 +/- 0.4 mm/yr, not 1.4 +/- 0.2 mm/yr as is described in lines 139 and 421; please make consistent. Table S2 a does not indicate the rate 1.4 +/- 0.2 mm/yr, but it is described in the footnote of Table S2 that “Underlined values from NAG_B coral are used to reconstruct the RSL at a larger time scale in Fig. 4”. Explain how the 1.4 +/- 0.2 mm/yr was derived from the underlined values (1.6 +/- 0.2 mm/yr and 1.2 +/- 0.2 mm/yr), and how those specific underlined values were chosen? Line 145: “ages of the undated plateaus are bounded by those of the plateau just above or below” Is there a stratigraphic relationship between the corals (i.e., does “above” and “below” refer literally to corals stacked on top of one another?) The extrapolation of ages of “undated plateaus” would make sense only if the stratigraphy supports it: a coral that grew physically on top of another coral must be younger than the coral beneath it. In Figure S1j, it appears the corals were all in separate locations. If this is the case, the relative ages of undated corals to the dated corals is not clear. Line 207: Typo mistake here “(fault 2) 2” Line 228: As per line 145, the sentence here “Some old corals are below younger generations” is vague. Does “below” mean that some corals are physically stacked on top of other older corals below? Line 231: Given that the site is located at the mouth of a river, and some corals are located amongst mangroves, can there also be sediment compaction? Lines 254-255: “Although we cannot decipher between shallow and deep dislocation models, the absence of tsunami traces could be in favour of deep ruptures for the smaller earthquakes.” The “absence of tsunami traces” seems contradictory given that the paragraph talks about documented tsunamis. Please clarify. Lines 269-270: “… and not of the same category as the largest earthquakes interpreted in the area.” The phrase is vague, it is not clear what “the area” and the “largest earthquakes” refers to. Lines 270-271, 277-278: Add line spacing between paragraphs. Lines 281-282: If interseismic subsidence is currently occurring, is this evident in the geomorphology of the site? Are the living microatolls recording continuous relative sea-level rise? Are there other independent drowned landscapes? Line 297: “Fossil microatolls permit to have extremely useful…” A subject is missing between “permit” and “to”? Line 299: “Stoddard” is spelt wrong Line 300: “Microatolls grow at a constant rate underwater” Corals can grow at variable rates throughout their lifetimes, causing the annual band thicknesses to vary across the microatoll. Line 301: “their development is controlled by sea-level variations…” Please change “sea-level” to RSL, to clarify that corals microatolls also record land-level changes Line 303: “excessive subaerial and sunlight exposure” It is not clear if excessive sunlight is needed for a coral microatoll to experience a diedown Lines 303-304: “Once the HLS is reached, the growth is restricted to the horizontal direction” This statement implies that the coral simply grows upwards until a limit (the HLS), then proceeds to grow horizontally. In reality, the HLS is the elevation a coral dies down to (from its pre-diedown highest level of growth), rather than the elevation it grows up to. Additionally, following a diedown, the coral polyps continue to grow radially outwards (both horizontally and vertically). Line 318: “the coral microatoll record a dark …” Add an ‘s’ to ‘record’ Line 320: “will provide information on the HLS of the coral …” change “HLS” to “HLG/HLS”. Depending on whether there was a diedown in a given year, the annual band could be representing a HLS (during a diedown), or HLG (between diedowns). Lines 332: “The most recent RSL record” The term “most recent” is vague. Modern? Lines 339-340: “the altitude of the external ring (younger part of the microatolls, i.e., most recent HLS)…” The outermost ring crest is not the most recent HLS, it is the most recent HLG. For fossil corals, the outermost ring crest is the most recent pre-diedown HLG (HLG just before the most recent diedown). Line 344: “… one being ball-shaped” What does a ball-shaped microatoll mean? Does it have concentric rings? If not, is it just a hemispherical coral head, and not a microatoll? Line 348-349: “The average vertical uncertainty on a single concentric ring measurement is of +/-3 cm”. This statement is vague. Is this referring to the precision of the levelling equipment? Or the natural variability in HLG across a given concentric ring, derived by multiple elevation measurements of the same ring? Lines 351-352: “Other errors related to erosion or compaction of the site are complex but small …” The coral slabs suggest significant erosion, as evident from the truncation of near-vertical growth bands. Lines 353-354: Were the Yamaguchi et al dates calculated using the same half-lives as the dates in this study? Line 363: Some coral slabs were sampled along the radius of the coral, not diameter. Lines 371-372: “… to distinguish and discard from our analysis new colonies that commonly develop on top of dead corals.” Colonies B and C of NAG_B look like they belong to an independent coral colony thats distinct from colony A. Line 389: Is the uncertainty in the sampling of coral chunks (described here to approximate the most external growth band “by several years”) accounted for in the uncertainty of the coral ages during subsequent interpretation? Line 408: We recommend shifting this part describing the “die downs” together with the part that introduces the HLS, as they are related concepts. Line 417: Clarify that the +/- 1 cm relative uncertainty in the HLS curve of a single microatoll relates to the uncertainty in the reorientation of the X-Ray and CT scan images (described in lines 377-378). What about band-counting uncertainty, associated with uncertain interpretation in the coral growth bands when trying to reconcile the annual bands with the U/Th dates? Lines 418-420: “… we observe an average range of 15 cm between measurements, leading to an estimated uncertainty of +/- 7 cm.” Does this mean that coeval diedowns were observed at different elevations in different corals, therefore deriving an uncertainty of +/- 7cm? If so, is this implying the natural variability in the HLS is +/- 7 cm? Line 426: “… permit to make…” A subject is missing between “permit” and “make” Lines 427-428: It is not clear how the 21 cm is derived. Supplementary line 95: How was mean sea level derived? Supplementary Table S2b: Event P – is there a typo for the age of the event (442 BC not 44 BC)? Supplementary Tabel S2b: what is the difference between the “Age” column and the “All-data age average” column? How were both derived? U/Th ages: • Were the samples pre-screened for secondary calcite? • NAGM_A1 and NAG6_18N have initial d234U values that deviate by more than 5 ‰ from modern seawater (145 ‰), suggesting possible open system behaviour. Were these used in the analyses? • Samples NAG321E_A1, NAG618H_A2, NAG3_18C, NAG3_18D, NAG3_18B, NAG5_21A have anomalously low 230Th/232Th ratio, suggesting presence of non-negligible amounts of initial (non-radiogenic) detrital Th. How sensitive are the ages to the atomic 230Th/232Th ratios used for the age corrections (footnote d)? Figure comments Figure 1 • What do “low frequency” and “very low frequency” earthquakes mean? • Difficult to identify/ distinguish the “Graben” and “Okinawa Trough” colours Figure 2 • Can the coral locations in panel b be drawn shown also in panel a? • It was very difficult to notice Yama_1 and Yama_2 in panel a as they are blue points on a blue background • Is there a clearer image of the reef in panel b? • Panel c: what is the datum for the elevation measurements? Specify in the caption or y axis title. • Panel c: Do the “Fossil coral” symbols refer to the “undated corals”? If so, how were they placed along the x axis? • Panel c: we recommend colouring the error bars according to the symbol colour, so it is easier to read which error bars belong to which samples • Panel c: it is difficult to distinguish the red and orange symbols • Panel c: why is the horizontal error bar for one of the living corals dashed? • Panel c: how was the tidal range measured? It seems odd that the living coral microatolls are out of the tidal range. It may also be useful to show the lowest astronomical tide as an estimate for extreme low water level at this site – which is what controls the microatoll diedowns. Figure 3: • The caption describes this as a relative sea level plot, but the y axes are labelled as “Altitude (m)”. Please clarify if this is indeed showing the absolute elevation of the corals, or relative sea level. • For NAG_B, which radius is shown in panel b, given that the full diameter was slabbed? • Uplift events estimated from the eventual death of coral microatolls should be minimum uplift estimates. We recommend annotating an arrow and indicating “> XXX cm” to illustrate that the uplift (or RSL fall) could be greater than the vertical red line drawn. Figure 4: As with figure 3, please clarify if the y axis is showing the absolute elevations of the corals or relative sea level. Figure 5: The ellipses for the deep and main megathrust in panel a do not match the cross section in panel b. In the bottom panel of b, fault 2 (deep megathrust) ends at the northern end of Ishigaki (middle panel of b), but in panel a the blue ellipse does not extend to the northern end of Ishigaki island. The bottom panel of b also shows the blue ellipse extending further beyond the main megathrust, but the northern limit of both ellipses in panel a are aligned. Figure S1b – S1h: There seems to be a mislabelling of the screw positions. Lines 363-365 describe that elevations were surveyed before the slabs were extracted, so it seems unlikely that the screws would be along the yellow line, and should instead be on the coral surface. Figure S2a: Please indicate scales for each of the coral slabs.", "reviewer_2": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature Communications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career Researchers who co-review manuscripts.", "reviewer_3": "This is an exciting paper that the authors present a comprehensive study of the history of megathrust earthquakes over the past five thousand years in the southern Ryukyu Trench using coral microatolls. The manuscript is well-organized and clearly written. I only have few minor comments and support publication. Lines 193-195: The inferred supercycles in the southern Ryukyu Trench are based on data from a single paleoseismic site. Are there any published studies using sediment cores from the seafloor, marine terraces, or tsunami boulders dated in Okinawa island that support the existence of these supercycles found in southern Ryukyu Trench?\"? Lines 232-238: The Okinawa Trough is rifting southward at a rate of approximately 5 cm/year in this region (e.g. Nishimura et al., PEPI, 2004). As a result, the convergence rate between the Philippine Sea plate and the southern Ryukyu forearc exceeds 8 cm/year (e.g. Chen et. al. https://doi.org/10.1029/2022GL098218). Authors need to refine interseismic modeling accordingly.", "reviewer_4": "Summary\nThis manuscript presents a new data set of relative sea-level over the last 5000 years on Ishigaki\nIsland (tropical Japan). The authors acquire the data from coral microatolls and interpret their\nnoisy record in form of three supercycles of 1’000 years with more earthquakes occurring in the\nlast parts of those supercycles. They use the 1771 Meiwa earthquake to suggest the end of the\nthird supercycle is approaching. I find the data somewhat over-interpreted"}, "author_response": ""}, {"round": 2, "reviewer_comments": {"reviewer_1": "Summary\nThis manuscript presents a new data set of relative sea-level over the last 5000 years on Ishigaki\nIsland (tropical Japan). The authors acquire the data from coral microatolls and interpret their\nnoisy record in form of three supercycles of 1’000 years with more earthquakes occurring in the\nlast parts of those supercycles. They use the 1771 Meiwa earthquake to suggest the end of the\nthird supercycle is approaching. I find the data somewhat over-interpreted", "reviewer_2": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature Communications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career Researchers who co-review manuscripts.", "reviewer_3": "The authors have thoroughly addressed all comments raised by me and the other reviewers. I have no further comments.", "reviewer_5": "The manuscript \"Evidence for megathrust earthquakes (...)\" provides clear evidence for paleoseismic activity along the Ryukyu trench that can reasonably be attributed to megathrust earthquakes. The interpretation of the coral record seems robust. However, the paper is poorly written and the figures are of poor quality. I recommend minor revisions to improve the quality of the manuscript. Line 26: replace \"infer\" with \"imply\". Write \"Furthermore, historical and geological studies show evidence of great tsunamis.\" But perhaps indicate when and where, as the sentence, as is, is vague. Line 27. \"Here, we use\". Use present tense and active voice. Line 28: \"relative sea level\". Singular. \"We use fossil microatolls in Ishigaki Island to reconstruct relative sea level in the Holocene. The coral record reveals several emergence episodes clustering 5-to-4 and 3-to-2 thousand years ago (ka) compatible with rapid uplift from megathrust earthquakes.\" Line 33: remove \"would\". \"The paleoseismic behavior implies a strong seismic hazard for the upcoming centuries. The devastating 1771 Meiwa earthquake and associated tsunami may mark the onset of a seismic supercycle.\" Line 36: remove \",\" in the middle of the sentence. Remove \"strong\". Line 43: \"evidence\". Singular. Remove \",\" before \"and\". Throughout the text: There are problems with virtually every sentence of the manuscript. I will stop correcting, but this is a major issue. They are many useful tools to correct and improve the writing. See DeepSeek, ChatGPT, Mistral, and the like. The language and writing must be seriously improved. Line 49: \"evidence\" -> \"record\". \"To this goal\" -> \"Toward this goal\" Line 97: \"motion\". Singular. Line 98: \"motions\"->\"displacement\". Line 215: \"2 ky apart\" Line 216: \"over the last 6 kr\". Line 238: Permanent or irrecoverable deformation is often attributed to splay faulting and folding in the accretionary prism. Ishagaki is a forearc high built by such deformation. See a recent review: Qiu Q. and S. Barbot, \"Tsunami excitation in the outer wedge of global subduction zones\". Earth-Science Review, 2022. Figure 1: add contours of slab2.0 depth. Is it \"Ruykyu\" or \"Ryukyus\"? Probably singular. Use similar color themes for the main map and the two insets. Fix the width of the black rectangles on the side. Some have a black frame, others don't. Date the aseismic events. Figure 2: can probably be split into two figures. Figure 2c: \"Distance along plot (km)\" What is \"plot\". Is it a plot of land? Or plot as in \"figure\"? Perhaps find a more suitable axis label. All the figures are pixelated. Hopefully, vector-graphic figures will be provided eventually."}, "author_response": "We thank the authors for their thoughtful revisions and detailed responses, which have strengthened the manuscript. The introduction of the confidence index and distinction between eroded and preserved HLG/HLS enhances the clarity and robustness of the record, and addition of event labels in Figure S2a makes the descriptions in the text easier to follow. The study provides important evidence of past tectonic uplift and potential subsidence along the southeastern coast of Ishigaki Island. We believe the manuscript remains a meaningful contribution to the field and offer a few additional suggestions for consideration prior to publication. 1. Erosion and sea-level tendency It is now much clearer with preserved/ eroded HLG/HLS illustrated in Figure S2b. The interseismic rates are also now inferred using the authors’ best estimate of the interseismic RSL trend using the longest, best-preserved record from NAG3_18H. However, eroded HLGs were still used to quantify coseismic uplift estimates, which we draw caution to. While it is true that eroded HLG/HLS have been used to estimate RSL changes (as cited in the examples by Sieh et al. 2008 and Meltzner et al. 2015), the focus for this is to quantify interseismic RSL rates. The premise for doing so is that the coral only experiences a diedown once it grows close enough to sea level that a temporary negative sea-level anomaly/ drop in sea level can cause a diedown. Therefore, the HLG immediately preceding a diedown (i.e., highest point on the coral before a diedown) best approximates HLS, even if the highest points are eroded. Additionally, the pre-diedown HLG is less affected by the magnitude of a given diedown event (which can be influenced by the severity of a sea-level anomaly), compared to the post-diedown HLS. (See Meltzner et al., 2010, section 3.3.2) In the examples regarding Sieh et al. 2008 and Meltzner et al 2015, coseismic uplifts were only estimated where there were preserved HLS points to indicate how far the corals died down to. Without a preserved HLS/ diedown, it is not clear if there truly was a diedown that later got eroded, or that there was no diedown to begin with. We recommend that events G and H should be assigned a lower confidence index due to the lack of preserved HLS and lack of agreement with coeval records. During the timing of events G and H when NAG_M was inferred to have died down, NAG6_18H recorded continued upward growth, even though NAG6_18H was at a higher elevation. The two coeval coral records do not seem to agree, even considering the +/- 7 cm HLS variability, so it seems possible that events G and H in NAG_M reflect erosion horizons, rather than diedowns. The concentric nature of this “step-down” marking event G is also not clearly illustrated in the field photo of NAG_M provided in the rebuttal: Therefore, we recommend that event G is downgraded from confidence index of 5 (the most confident) to a lower confidence, and hedge the wording regarding the “major sea-level decrease” and “second sudden sea-level decrease” in supplementary lines 142-144. This is unlike event D, where coeval diedowns were recorded in both NAG6_18Q and NAG6_21C (despite NAG6_21C not having a preserved diedown); and event E, where NAG6_21C died down completely around the same time when NAG6_18Q recorded a diedown (despite there being no preserved HLS in NAG6_21C). Supplementary text lines 117-119: the two emergence events inferred are from eroded parts of NAG_B – there may not have been a diedown as the HLS is not preserved. Figure S2a: The upward growth of a coral does not necessarily indicate sea-level rise tendency (as shown by black lines in Figure S2a indicating “apparent RSL increase”). The coral could still be catching up to sea level, during which upward growth reflects the coral growth rate, not sea level. Instead, it is the stepping upwards of successive concentric rings/ diedowns that indicate RSL rise, not the increase in elevation of successive growth bands. Similarly, RSL fall should be interpreted as successively lower concentric rings/ HLS, not successively lower growth bands. Supplementary text lines 138-139: the words ‘first HLS’ were removed, and eroded HLG/HLS shown in Figure S2b, but the text and interpretation did not change (i.e., erosion is still interpreted as RSL fall, as indicated by the blue arrow in Figure S2b). Supplementary text lines 180-181: similar comment to the above; RSL decrease of -1.4 mm/yrs is inferred from eroded HLG, but HLS are not preserved Supplementary text line 191: similar comment to the above; RSL decrease is inferred from eroded HLG, but HLS are not preserved 2. Offset of NAG3_21E and NAG_B not explained by HLS variability applied It is unclear from the text how the HLS variability of +/- 7 cm was derived from two corals with coeval diedowns. Which corals / diedowns was this based on? In lines 102-104 of the main text, the current explanation for why NAG3_21E is lower than NAG_B seems to suggest that this offset could be explained by natural variability in the coral HLG/HLS. However, the elevation difference between the preserved HLS of the two corals is nearly 35 cm at 301 yrs BC based on figure 3, which is larger and seems to contradict the the +/- 7 cm applied as the HLS variability. Could NAG3_21E have slumped? Or NAG_B be ponded at a higher elevation? Alternatively, is it possible that NAG3_21E was growing too deep compared to NAG_B, such that the uplift event that killed NAG_B (event Q) was not large enough to kill NAG3_21E? If so, then this would provide an upper limit on the amount of uplift for event Q. This would require that the 301 yrs BC and 329 yrs BC inferred diedowns were not diedowns. From figure S2a, the 342 yrs BC diedown is not shown on NAG3_21E, and the 301 yrs BC diedown in NAG3_21E is also not as clear as that in NAG_B. Minor comments 1. Line 99: Should be “glacial isostatic adjustment” not “glacio-isostatic adjustment” 2. Lines 181-183: Please specify which corals were used to infer possible coseismic subsidence (it’s indicated in figure 4 caption as “potential sudden RSL increase event”) but could be described better in the main text. It will also be helpful to add the inferred subsidence event(s) to table S2b 3. Lines 234-236: Consider adding that different glacial isostatic adjustment models could also produce variable RSL predictions. 4. Supplementary text line 119: 26’ +/- 8 yrs BC seems to be a typo 5. Supplementary lines 190-191: NAG6_21C was described as starting to record HLS variations at 2436 years BC, but there was no preserved diedown/ HLS here. 6. Figure 3: Similar to figure 4, also distinguish the confidence index in the vertical lines that indicate diedown events? 7. Figure 4: Rather than bolding letters with confidence index 5, consider having distinct line styles for the different confidence indices. 8. Please check that table S2b is consistent with figure 3. For example, events N and O are indicated as >15 cm and > 18 cm in figure 3 and in the main text, but reported as ~15 cm and ~18 cm in table S2b. Please check for other events as well. 9. NAG3_18H, events N and O: • It seems odd that the erosion line was not extrapolated to the diedown for event N • Why is the vertical downward dashed arrow for event O (which indicates the tectonic diedown) drawn from the HLG/ring crest after the diedown, not the pre-diedown HLG? • Rather than crudely estimating the amount of erosion based on the dashed blue erosion lines drawn, since the uplift amounts are reported and described as minimum uplift estimates anyway, would it be better to simply take the difference between the highest (eroded) pre-diedown HLG and the preserved HLS elevation as the minimum uplift estimate? Alternatively, perhaps extrapolate the pre-diedown HLG from the RSL trend prior to the diedown, similar to the examples cited in the rebuttal for Sieh et al. 2008 and Meltzner et al. 2015? Redacted 10. The complete die-off of any individual coral microatoll colony may not always require a tectonic uplift event. By “temporally discontinuous” records, we refer to the fact that even in relatively tectonically stable locations, coral microatoll ages often cluster into age clusters. This may be due to biological life spans of the coral microatolls, as according to Smithers (2011): “Large microatolls demonstrate continuity of growth for several centuries, but most microatolls are < 2 m in diameter. This probably reflects the relatively low probability of extended survival in shallow reef environments where corals are most vulnerable to a range of debilitating conditions. Microatolls on the northern GBR surveyed by Scoffin and Stoddart (1978, p.105) averaged 0.5 m in diameter, although other than stating that they “vary in size from a few centimeters to a few meters” no statistical description of the variance was offered.” Redacted References Hallmann, N., Camoin, G., Eisenhauer, A., Botella, A., Milne, G. A., Vella, C., ... & Fietzke, J. (2018). Ice volume and climate changes from a 6000 year sea-level record in French Polynesia. Nature communications, 9(1), 285. Meltzner, A. J., Sieh, K., Chiang, H. W., Shen, C. C., Suwargadi, B. W., Natawidjaja, D. H., ... & Galetzka, J. (2010). Coral evidence for earthquake recurrence and an AD 1390–1455 cluster at the south end of the 2004 Aceh–Andaman rupture. Journal of Geophysical Research: Solid Earth, 115(B10). Meltzner, A. J., Sieh, K., Chiang, H. W., Wu, C. C., Tsang, L. L., Shen, C. C., ... & Briggs, R. W. (2015). Time-varying interseismic strain rates and similar seismic ruptures on the Nias–Simeulue patch of the Sunda megathrust. Quaternary Science Reviews, 122, 258-281. Sieh, K., Natawidjaja, D. H., Meltzner, A. J., Shen, C. C., Cheng, H., Li, K. S., ... & Edwards, R. L. (2008). Earthquake supercycles inferred from sea-level changes recorded in the corals of west Sumatra. Science, 322(5908), 1674-1678. Smithers, S. (2011). Microatoll. In D. Hopley (Ed.), Encyclopedia of Modern Coral Reefs: Structure, Form and Process (pp. 691–696). Springer Netherlands. https://doi.org/10.1007/978-90-481-2639-2_111 Tan, F., Horton, B. P., Ke, L., Li, T., Quye-Sawyer, J., Lim, J. T., ... & Meltzner, A. J. (2024). Late Holocene relative sea-level records from coral microatolls in Singapore. Scientific Reports, 14(1), 13458."}, {"round": 3, "reviewer_comments": {"reviewer_1": "No further comments", "reviewer_2": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature Communications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career Researchers who co-review manuscripts."}, "author_response": "The revised manuscript shows significant improvement over the previous version and has satisfactorily addressed most of the comments raised by me and the other reviewers. I have no further comments. We thank the authors for taking the time and effort to address the comments raised. We find that the manuscript is now much clearer in its assumptions and inferences, providing important data to enhance understanding of seismogenic nature of the Ryuku subduction zone. It is especially interesting that the clusters of high-confidence abrupt uplifts coincided with the ages of tsunami deposits, with lack of both coral data and tsunami data between the inferred supercycles. We only have a few remaining minor comments and recommend publication thereafter. Main text 1. Line 101: please specify over what time period the quoted GIA rate is referring to. 2. There are places where there seems to be missing line spacing between paragraphs. Please check to ensure that the spaces are added accordingly (e.g., between lines 102–103, 135–136…) 3. Line 405: Here it says you use Zach-3 method, but in table S2a the footnote mentions you use the Zach-1 method. Please check to make sure it is consistent. 4. Line 436: Earthquakes inferred from within the coral slices are assigned with confidence rating of 5. While we appreciate that the authors downgraded the confidence ranking of event G to grade 4 to address our earlier comment, the confidence rating now is inconsistent with the description in the text (because event G is inferred from the slab but is assigned grade 4 and not 5). We suggest adding a new confidence rank to differentiate between diedowns inferred from preserved HLG/HLS (ranking of 6) in slabs versus those inferred from eroded HLG/HLS in slabs (ranking of 5). 5. Figure 1: panel labels are missing (e.g., 1a, 1b) 6. Figure 4: where are the dotted lines? 7. Figure 4: we appreciate the attempt to address our earlier comment by introducing different line thicknesses for different confidence rankings. However, the different line thicknesses are more difficult to read than the bolded letters initially. Could you use different line symbols (e.g., dotted, dashed, solid, twodash) to represent different confidence rankings, and add a legend for this? 8. Figure 4: please arrange the red line for event R to be on top of the black vertical line as now it is difficult to see the red line. 9. Figure 5: please indicate in the caption which corals/ uplift events were used to constrain the elastic dislocation models. Maybe I missed this but was this mentioned explicitly in the text? If not, please add a line to clarify this. Supporting text 1. Line 209: change “at the sooner” to “at the soonest” 2. Table S2b: to make it easier for the reader, please add to the description in the table header that a confidence ranking of 5 is the best, so that the reader does not have to look for this in the text."}]

Decision letter