Radon degassing triggered by tidal loading before an earthquake

The presence of anomalous geochemical changes related to earthquakes has been controversial despite widespread, long time challenges for earthquake prediction. Establishing a quantitative relationship among geochemical changes and geodetical and seismological changes can clarify their hidden connection. Here we determined the response of atmospheric radon (222Rn) to diurnal tidal (K1 constituent) loading in the reported 11-year-long variation in the atmospheric radon concentration, including its anomalous evolution for 2 months before the devastating 1995 Kobe earthquake in Japan. The response to the tidal loading had been identified for 5 years before the occurrence of the earthquake. Comparison between these radon responses relative to crustal strain revealed that the response efficiency for the diurnal K1 tide was larger than that for the earthquake by a factor of 21–33, implying the involvement of crustal fluid movement. The radon responses occurred when compressional crustal stress decreased or changed to extension. These findings suggest that changes in radon exhaled from the ground were induced by ascent flow of soil gas acting as a radon carrier and degassed from mantle-derived crustal fluid upwelling due to modulation of the crustal stress regime.


Results and discussion
Atmospheric radon concentration was monitored hourly between 1984 and 1995 except for a period in 1989 at the Kobe Pharmaceutical University (KPU). The monitoring site was located on the RAFZ striking along the foot of Mount Rokko, approximately 27 km north-east from the epicentre (Fig. 1a). Figure 2a shows an 11-year variation in atmospheric radon concentration measured there. The radon concentration increased from September 1994, and subsequently, exceeded the normal variation range 2 months prior to the earthquake, which could be viewed as a precursory change 26 . The detrended data showed that the radon concentration change reached up to 10 Bq m −3 (Fig. 2b). Meanwhile, the crustal strain also changed from contractional to extensional, measured at the Rokko-Takao (RTK) Observatory situated in Mount Rokko 19 (Fig. 1c). The detrended data showed that the change in the areal extensional strain reached up to 5.2 × 10 -6 . A coincidence of these changes indicates that the increased radon exhalation from the ground was caused by crustal deformation preceding the earthquake 28 . Based on these findings, as a measure of the sensitivity of radon to crustal deformation, the ratio of radon (ΔC e ) to areal strain (Δε e ) change related to the earthquake was determined to be 1.9 × 10 6 Bq m −3 , considering a radon change of 10 Bq m −3 associated with an areal crustal strain change of 5.2 × 10 -6 .
For the two sets of hourly continuous radon data in Fig. 2a, fast Fourier transform (FFT) analysis was applied to investigate the radon response to periodic crustal strain change induced by tidal loading. Dominant tidal loading cycles have diurnal, luni-solar (K1), and lunar (O1) constituents, and semidiurnal, principal solar (S2), and principal lunar (M2) constituents (Table 2). Their relation in amplitude-magnitude is M2 > K1 > O1 > S2 in Kobe. Figure 3a shows the amplitude spectra for the radon datasets in 1984-1988 and 1990-1994. A few sharp peaks were identified clearly at ≤ 24-h cycles in both datasets, most of which were related to meteorological factors 29 .
The amplitude spectra around the targeting tidal cycles are enlarged in Fig. 3b. In the shorter-cycle range, no sharp peaks were identified at the cycle corresponding to the M2 constituent. Although a sharp peak was identified at the 12.000-h cycle, it is attributed to meteorological factors (e.g., air pressure and wind speed) rather than the S2 constituent because, despite the larger M2 tidal loading than the S2 loading in the study area, the 12.421-h cycle was not seen in the radon time series.
In contrast, significantly different spectra were observed in the 1984-1988 and 1990-1994 datasets for the longer-cycle range. Only one sharp peak was observed at the 24.000-h cycle for the 1984-1988 dataset, whereas a weak, but clear peak was also identified at the 23.934-h cycle as well as 24.000-h cycle for the 1990-1994 dataset. The radon change at the 24.000-h cycle was attributed to thermal stabilisation/destabilisation in the atmospheric boundary layer driven by solar radiation 30 rather than the solar S1 tidal loading because of no peak at the 25.819-h cycle induced by the larger O1 constituent. However, the 23.934-h cycle was not induced by any meteorological factors, which is confirmed by no corresponding cycle peak for the 1984-1988 dataset. Accordingly, the 23.934-h cycle of radon change was induced by K1 tidal loading. This is the first evidence of periodic release of geochemical gas from the ground to the atmosphere driven by periodic crustal loading. The amplitude spectrum revealed that the radon concentration change reached up to 0.44 Bq m −3 (95% confidence interval: 0.36-0.59 Bq m −3 ). It should be noted that no information regarding the radon response to the earth tides was given in 1989 when the radon monitoring was unavailable.
The ratio of radon concentration (ΔC t ) to areal strain (Δε t ) change related to the K1 tidal loading was evaluated further; the areal strain induced by solid earth and ocean tide loading was calculated theoretically using the GOTIC2 software (ver. 2004.10.25) developed by the National Astronomical Observatory, Japan 31 . The calculated areal strain by the K1 tidal loading was 0.93 × 10 -8 based on the RTK site locality information. Thus, the ratio (ΔC t /Δε t ) was evaluated as 47 × 10 6 Bq m −3 (95% confidence interval: (39-63) × 10 6 Bq m −3 ) considering radon changes of 0.44 Bq m −3 (95% confidence interval: 0.36-0.59 Bq m −3 ) associated with an areal crustal strain change of 0.93 × 10 -8 by the K1 tide. The values were larger than the ΔC e /Δε e value for the case preceding the 1995 Kobe earthquake by a factor range of 21-33. Two potential mechanisms can be considered to induce the response of atmospheric radon to crustal deformation. One qualifier is the development of crustal rock microfractures due to applied loading. Microfracture development would increase (1) radon production rate due to increased rock grain surface area and (2) radon transport pathways (i.e., permeability) due to connections between isolated microcracks 32,33 , which would consequently elevate radon exhalation. However, microfracture development is an irreversible process; therefore, and the Kobe earthquake (star), measurement sites of atmospheric radon concentration and crustal strain (green squares), those of helium isotopes in mineral and hot springs (triangles), and active faults (red lines). The triangles are coloured depending on the corrected 3 He/ 4 He ratio normalised by the value (R a : 1.39 × 10 -6 ) for the air (see the inset). Underneath the western area depicted by an open ellipse, the subducting Philippine Sea Plate might be locked during 1991-1993 (ref. 40 ). ABU Abuyama, AMA Amagase, DNZ Donzurubou, KPU Kobe Pharmaceutical University, RTK Rokko-Takao, RTR Rokko-Tsurukabuto, ATFZ Arima-Takatsuki fault zone, RAFZ Rokko-Awaji fault zone. Data on seismicity were obtained from the Unified Japan Meteorological Agency Earthquake Catalogue (http://evrrs s.eri.u-tokyo .ac.jp/tseis /jma1/index .html) and those on topography, helium isotopes, and active faults were obtained from refs. 51,17,16 , respectively. (b) Variations in principal crustal strain rates observed at RTR from 1977 to 1995 (Ref. 21 ). (c) Time-series variations in crustal strain observed at RTK from November 1993 to February 1995 (Ref. 19 ). The map was created using the Generic Mapping Tools 52 . www.nature.com/scientificreports/ it cannot reproduce the periodic pattern induced by the tidal loading. In addition, radon that originates deeper cannot reach the ground surface because it decays with a mean lifetime of 5.5 days. Additionally, assuming a soil porosity of 0.4, effective diffusion coefficient of 10 -6 m 2 s −1 , and subsurface ascent air flow velocity of 10 -6 m s −1 , the diffusion and advection lengths of radon are estimated to be 0.7 and 1.2 m in soil, respectively 34,35 . It implies that the exhaled radon is not produced in the basement rocks, but only in the uppermost soil layer. Therefore, microfracture development in crustal rocks is an unlikely primary physical mechanism for current radon exhalation responding to crustal deformation. The other potential mechanism is a change in subsurface vertical flow of air carrying the radon. The air flow modulates the vertical soil radon profile; consequently, radon is exhaled from the ground surface without a change in radon production rate by newly forming microfractures. King 4 reported profile changes in soil radon concentration responding to an earthquake with magnitude 4.0 on the San Andreas Fault system (USA). At the San Jacinto Fault (USA), Birchard and Libby 36 observed increases and decreases in soil radon concentrations in the quadrants of seismic compression and dilatation, respectively, following an earthquake. These results could be qualitatively explained by the induction of subsurface ascending or descending air flow from a change in loading stress. Here, for several months prior to the 1995 Kobe earthquake, the increase in atmospheric radon concentration corresponded to an increase in the groundwater discharge rate at the RTK site and changes in the crustal strain from widespread contraction to extension, namely at RTK, Abuyama (ABU), Amagase (AMG), and Donzurubou (DNZ) sites 19,25 . Unlike the case of Birchard and Libby 36 , which indicated a phenomenon of squeezing out soil gas through soil-pore volume compression, the current case indicated that crustal fluids upwelled due to subsurface pathway expansion, and radon carriers degassed from the fluids induced subsurface ascent air flow. In the study area, non-volcanic mantle-derived fluids upwelled at hot springs and mineral springs 17,18,37 . Umeda et al. 37 noted that the NW-SE oriented extensional stress field generating open fractures favoured upwelling of mantle-derived fluids in the area. The fluids contained abundant quantities of gases such as carbon dioxide that can carry radon toward the earth's surface. Similar association of extensional tectonics elevating carbon dioxide degassing from the deep crust into the atmosphere was recently recognised worldwide 38,39 . Tamburello et al. 39 also noted an important role of normal or strike-slip faults with large dip angles on carbon dioxide degassing in addition to extensional tectonics. Combined with the tectonic settings in the study area, induced subsurface ascent air flow due to mantle-derived fluid upwelling is likely a physical mechanism of radon exhalation responding to crustal deformation.

Scientific
Questions arise as to why, unlike 1984-1988, the radon exhalation responding to the K1 tidal loading occurred in 1990-1994 and why its sensitivity was high compared to the case of earthquake occurrence. We propose that initiation of the tidal radon response can be regarded as a precursory phenomenon and it is caused by upwelling deep crustal fluids. During the latter period, precursory phenomena were determined from the geodetic and seismologic observations (elevation, crustal strain, and seismic activity; Fig. 1b,c, Table 1). These observations indicated that the crustal stress regime change occurred almost synchronously in the broad area surrounding Table 1. Occurrence of geophysical, hydrological, and geochemical precursors before the Kobe earthquake (17 January 1995). a Ide and Tanaka 40 did not report the decrease in background seismicity as a precursor. b Nishinomiya is located 30 km northeast of the mainshock epicentre.   www.nature.com/scientificreports/ Okayama, west of the RAFZ, quiescence of background seismicity occurred and positive correlation between the numbers of background seismicity and deep tremors predicted by the tide disappeared from 1991 to 1993 40 , which implied plate boundary contact beneath the area. However, unlike the stress state around the RAFZ, compressional stress decreased at RTR, and elevation uplift and tide phase and amplitude modulation in crustal strain occurred at RTK of the RAFZ 19,21,24 . The locally different stress fields could be explained by a change in pore pressure related to the upwelling of deep crustal fluids. They reached the RAFZ and increased pore pressure around 1990, leading to RTK site uplift as well as compressional stress field reduction. The stress field changes from fluid upwelling could also initiate the radon response to diurnal tide variation. The upwelled fluids would increase pore pressure and reduce effective stress on the fault plane, which consequently generates the response to the small diurnal tide variation magnitude. The stress-induced change in the tidal crustal gas discharge response was evident from the reports by Weinlich et al. 7,8 demonstrating disturbances in diurnal cycles of radon and carbon dioxide discharges induced by earth tides due to crustal stress redistribution following earthquake swarms in the Eger Rift, Czech Republic. The dynamic stress change induced by the passage of seismic waves originating from the Izmit earthquake (magnitude 7.4) in 1999 caused an apparent modification of the spectrum of confined-well water radon concentration and the disappearance of the tidal (M2) signal in the radon time-series 3 . In connection with the role of fluids, recent laboratory triaxial creep experiments for quartzrich sandstone samples revealed that pore-pressure oscillation reduced effective stress and caused cyclic deformation (acoustic emission events) of the samples 41,42 . In the present study area, the water injection experiment conducted at the Nojima Fault on Awaji Island, the southern part of the RAFZ following the Kobe earthquake showed that water injection into the fault amplified the diurnal crustal strain cycle produced by the earth tide 43,44 . Therefore, the diurnal tide phase and amplitude modulation in crustal strain records at RTK, observed as an intermediate-term precursor before the earthquake (Table 1), were possibly affected by external fluid injection 24 , that is, deep crustal fluid upwelling into the RAFZ. The amount of degassing from the fluids was controlled and amplified by the tidal loading thereby initiating the atmospheric radon response to the diurnal tide. The deep crustal fluid contribution was supported from the present study finding that the radon sensitivity (ΔC t /Δε t ) to   www.nature.com/scientificreports/ tidal high-frequency strain change was larger than that (ΔC e /Δε e ) to the monotonic (i.e., low-frequency) strain change during 4 months before the earthquake, implying involvement of crustal fluid movement. Some previous studies 45,46 did not consider changes in radon concentration related to earth tides. This is likely because a peak corresponding to the strongest M2 tide could not be identified from spectral analysis. However, we determined the weak peak of radon at 23.934 h as the earth-tide (K1) signal, as in previous studies by Crockett et al. 47 , Groves-Kirkby et al. 13 and Weinlich et al. 8 who reported signals induced by weaker earth tides. There are two possible reasons for the radon signal at the M2 tidal period being absent. From a meteorological viewpoint, the atmospheric radon concentration reported herein was primarily controlled by the variation in atmospheric stability due to the 24-h solar cycle. The night-time thermal stabilisation attained in the atmospheric boundary layer leads to the accumulation of radon exhaled from the crust, whereas daytime thermal destabilization disperses accumulated and exhaled radon immediately via convection 29,30,48 . Such a strong meteorological effect would possibly mask a radon signal related to the M2 tide, which is stronger in tidal force but higher in frequency than the K1 tide. From a tectonic viewpoint, modulation of the crustal stress regime preceding the earthquake may cause signal amplification at a specific period. As described above, in concordance with the crustal stress modulation, only diurnal variations in crustal strain were altered at the site along the RAFZ 24 . A similar phenomenon was observed in the water injection experiment conducted after the earthquake 43,44 . In addition, the deformation of quartz-rich sandstone samples is dependent on the frequency of pore-pressure oscillation rather than its amplitude 49 ; samples exposed to lower-frequency pore-pressure oscillation exhibited substantial deformation. Accordingly, it can be considered that modulation of the crustal stress regime during 1990-1994 amplified only the radon signal related to the K1 tide.
We presented the two cases of precursory changes in atmospheric radon concentration associated with the crustal strain changes preceding the earthquake. The crustal strain-radon change relationships apparently had different radon sensitivities to crustal strain for the earth tide and earthquake. The sensitivity of radon to crustal deformation has been only surficially examined, even including the present case, so that further observations are required to generalise their quantitative relationship and clarify the physical mechanism behind it. Multi-station observations are also important to confirm the precursory changes in atmospheric radon concentration. Such observations could not be made in the present study because the work was performed at KPU only. After the earthquake, we established a nation-wide network for the measurement of atmospheric radon concentration by radioisotope facilities throughout Japan 50 , which could be used to confirm the validity of identified precursory changes preceding future earthquakes.
The atmospheric radon variation reported here contained several-years-long and several-months-long precursory changes associated with the 1995 Kobe earthquake governed by the compression to extension stress regime change. Additionally, the stress regime changes should accompany the discharge of radon carriers like carbon dioxide, indicating upwelling of non-volcanic mantle-derived fluids through the strike-slip active faults with large dip angles as identified and characterised in the study area. These two characteristic processes constrained the geogenic radon variation before the earthquake.

Methods
Area description. The study area is in a high-strain-rate field (> 0.1 × 10 -6 year −1 ) corresponding to the western margin of the Niigata-Kobe Tectonic Zone 53 . Right-lateral active faults are distributed striking NE-SW (Rokko-Awaji fault zone; RAFZ) and ENE-WSW (Arima-Takatsuki fault zone; ATFZ) 16 . The Kobe earthquake (on 17 January 1995) initiated at north of Awaji Island and ruptured the RAFZ north-eastward. Mount Rokko, composed of granite and granodiorite, is behind the centre of Kobe City. Kobe Pharmaceutical University (KPU), the observation station for atmospheric radon concentration, was located at the foot of Mount Rokko, whereas the underground Rokko-Takao (RTK) 19 and Rokko-Tsurukabuto (RTR) 21 observatories for crustal strain were located inside Mount Rokko at 240 and 100 m in depth, respectively. In the eastern part of the study area, the underground observatories were also located at Abuyama, Amagase, and Donzurubou (ABU, AMA, and DNZ) 25 .
Several-years-long to several-months-long geophysical and hydrological precursors were determined as well as atmospheric radon concentration before the Kobe earthquake. For several-years-long geophysical precursors, the elevation was uplifted by 1.5 mm inside the RTK observatory for 5 years before the earthquake 19 . At the same site, phase and amplitude changes had been observed in the O1 constituent of a laser strainmeter record since 1992 24 . The contractional strain rate estimated by extensometer records changed from ENE-WSW to NE-SW and decreased in magnitude from 1991 to 1994 at the RTR site 21 (Fig. 1b). Along the RAFZ, seismic quiescence had occurred since 1990 for background seismic activity on the RAFZ 22 . In the east of the study area, N-S components of strain records showed an increase in contraction during 1990 to 1992 at the ABU, AMA, and DNZ sites 25 . At off Okayama, west of the RAFZ, background seismicity appeared to decrease during 1991-1993 although Ide and Tanaka 40 did not note it. In wide area (34-36° N; 134-136° E) including the present study area, the b value 23 increased from 1991 to 1993, and then it decreased in 1994 for seismic activities with magnitude > 2. For severalmonths-long geophysical and hydrological precursors, N-S oriented crustal strain changed to extension after mid 1994 at the ABU, AMA, and DNZ sites 25 . Similarly, at the RTK site, crustal strain measured by borehole strain meters switched from contraction to extension from September 1994 (Fig. 1c), which accompanied a gradual increase in groundwater discharge rate from November 1994 with an intermittent peak triggered by the external force of a distant earthquake 19 . The reported precursory observations, including geochemical precursors 6,20,26 , are briefly summarised in Table 1.
Radon monitoring. Atmospheric radon concentration was measured hourly using a flow-type ionisation chamber (Model NAG513, Fuji Electric Systems Co., Ltd, Japan) with a large detection volume of 1.8 × 10 -2 m 3 . www.nature.com/scientificreports/ Air at 5 m above the ground was sampled at a flow rate of 15.6 m 3 min −1 , and a part of the filtered air was introduced into the ionisation chamber at a flow rate of 6.7 × 10 -3 m 3 min −1 . Continuous monitoring occurred between 1984 and 1995 except for 1989 due to system maintenance. Details can be found in Yasuoka and Shinogi 26 .

Scientific Reports
Periodic analysis. FFT was applied to the two radon datasets (1984-1988 and 1990-1994) to derive the amplitude spectra. Welch's method 54 was adopted to reduce estimation uncertainty (i.e., 95% confidence interval). In the 1984-1988 time-series, the first 32,768 (= 2 15 ) data were analysed by FFT with the Hamming window function to obtain its power spectrum density. Subsequently, another 32,768 data shifted from the first data by 720 data (approximately 1 month) forward were analysed to obtain the second power spectrum density. This procedure was repeated until the analysis of the whole 1984-1988 time-series was completed. Then, based on the obtained power spectrum densities, the average power spectrum density with 95% confidence interval was calculated and it was converted into an amplitude spectrum. The same procedure was applied to the 1990-1994 time-series to obtain its amplitude spectrum. In the periodic analysis, the software Signal Processing Toolbox of MATLAB ver. R2018b (the MathWorks, USA) was used. The frequency resolution of the spectral analysis was 0.000143 h −1 . This is equal to approximately 0.013 h of time resolution at the 24.0-h period. In addition, the spectral leakage, as a ratio of the amplitude at the 23.934-h period to that at the 24.000-h period, was determined to be 0.005. Figure 1a was created using data on seismicity, topography, active fault distribution, and helium isotope ratio. The Unified Japan Meteorological Agency Earthquake Catalogue is available at: http://evrrs s.eri.u-tokyo .ac.jp/tseis / jma1/index .html. The topographic data 51 are available at: https ://www.ngdc.noaa.gov/mgg/globa l/. Active fault distribution 16 is available at: https ://gbank .gsj.jp/activ efaul t/index _e_gmap.html. Data on the helium isotope ratio were obtained from Sano and Wakita 17 . Figure 1b,c are from refs. 21,19 , respectively. The data on atmospheric radon that support the findings of this study are available from the corresponding author upon reasonable request. www.nature.com/scientificreports/ Publisher's note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Data availability
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.