Next Article in Journal
Identifying Probable Submarine Hydrothermal Spots in North Santorini Caldera Using the Generalized Moments Method
Next Article in Special Issue
Ionospheric Total Electron Content (TEC) Anomalies as Earthquake Precursors: Unveiling the Geophysical Connection Leading to the 2023 Moroccan 6.8 Mw Earthquake
Previous Article in Journal
Flood-Prone Zones of Meandering Rivers: Machine Learning Approach and Considering the Role of Morphology (Kashkan River, Western Iran)
Previous Article in Special Issue
Improving the Estimation of the Occurrence Time of an Impending Major Earthquake Using the Entropy Change of Seismicity in Natural Time Analysis
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Fractal Patterns in Groundwater Radon Disturbances Prior to the Great 7.9 Mw Wenchuan Earthquake, China

1
Centre for Earthquake Studies, National Centre for Physics, Shahdra Valley Road, P.O. Box 2141, Islamabad 44000, Pakistan
2
Department of Industrial Design and Production Engineering, University of West Attica, Petrou Ralli & Thivon 250, Aigaleo, GR-12244 Athens, Greece
3
Key Laboratory of Geo-Detection, Ministry of Education, School of Geophysics and Information Technology, China University of Geosciences, Beijing 100083, China
*
Author to whom correspondence should be addressed.
Geosciences 2023, 13(9), 268; https://doi.org/10.3390/geosciences13090268
Submission received: 19 June 2023 / Revised: 27 August 2023 / Accepted: 1 September 2023 / Published: 4 September 2023
(This article belongs to the Special Issue Precursory Phenomena Prior to Earthquakes 2023)

Abstract

:
This study reports a fractal analysis of one-year radon in groundwater disturbances from five stations in China amidst the catastrophic Wenchuan ( M w = 7.9) earthquake of 12 May 2008 (day 133). Five techniques are used (DFA, fractal dimensions with Higuchi, Katz, Sevcik methods, power-law analysis) in segmented portions glided throughout each signal. Noteworthy fractal areas are outlined in the KDS, GS, MSS data, whilst the portions were non-significant for PZHS and SPS. Up to day 133, critical epoch DFA-exponents are 1.5 α < 2.0 , with several above 1.8. The fractal dimensions exhibit Katz’s D around 1.0–1.2, Higuchi’s D between 1.5 and 2.0, and Sevcik’s D between 1.0 and 1.5. Several power-law exponents are above 1.7, and numerous are above 2.0. All fractal results of the KDS-GS-MSS are further analysed using a novel computerised methodology that locates the exact out-of-threshold fractal areas and combines the outcomes of different methods per five, four, three, and two (maximum 13 combinations) versus nineteen M w 5.5 earthquakes of the greater area. Most coincidences using different techniques are before the great Wenchuan earthquake and after the earthquake. This is not only with one method but with 13 different methods. Other interpretations are also discussed.

1. Introduction

When it comes to severe natural disasters, earthquakes stand out since they may result in great losses of lives and property. Residents of large cities may experience severe effects from the massive quantity of energy produced during an earthquake, especially if the epicentre of the earthquake is near. Catastrophic earthquakes are unavoidable as natural phenomena but are highly challenging to foresee [1]. Consequently, the search for trustworthy seismic precursors is one of the greatest difficulties of science, and significant efforts have been made for many years in this area [2,3,4,5,6,7,8,9]. The issue of earthquake prediction is still open [10]. Earthquakes are inherently complex; thus, several techniques and multi-level strategies are required for prediction [7]. In regions where severe earthquakes are possible, associated predictions call for a gradual downscaling of time, location, and magnitude [2]. Along with the electromagnetic disturbances in ULF, LF, HF, and VHF frequencies, which can indicate approaching earthquakes [4,9], radon-222 (hereafter, radon) is a long-established precursor of impending seismic activity [2,7,11]. Radon is an inert gas created by the radioactive decay of the 238U series with a half-life of 3.86 days. When it disintegrates, it dissolves in soil pores and liquids before moving on to the surface, subsurface water, and atmosphere. Because radon can travel long distances in water and soil, it can be detected far from its generation location, and due to this, it is very significant in earthquake-related studies [1]. As the above reviews mention [2,7,11], there is a great number of papers reporting pre-seismic changes in radon in soil gas, groundwater, wells, thermal spas, and atmosphere. As a result, there is a substantial body of studies examining the relationship between the emission of radon and seismic activity [10].
The motivation of the present research is to investigate whether the abnormal behaviour of groundwater radon may be due to the catastrophic Wenchuan ( M w = 7.9) shallow (depth = 19 km) earthquake, which occurred on 12 May 2008 (calendar day; hereafter, d a y 133) along the Longmenshan fault (31.0° N, 103.4° E) in Sichuan Province, China [12,13] (Figure 1). This quake was the most destructive one in China since 1976 and the second most devastating seismic shock of this century after the great Sumatra earthquake of 2004 [14]. The groundwater data consist of one-year measurements received from China between 1 January 2008 and 31 December 2008 by five different stations (Figure 1) with epicentral distances between 105.6 km and 526.0 km (Table 1). The possible connection between the variability in geochemical signals and the seismicity due to the Wenchuan earthquake can help to delineate the underlying geophysical processes since the data include information about the subsurface dynamics. Detrended fluctuation analysis, the use of fractal dimensions (FDs) with different methodologies, and power-law fractal analysis are very reliable fractal analysis methods [1,15,16]. These methods can recognise the underlying scaling and long-range features that are characteristic for the unavoidable occurrence of earthquakes even in noisy and non-stationary time series.

2. Materials and Methods

2.1. Experimental Aspects

2.1.1. Earthquake Activity

The Wenchuan earthquake occurred along the Longmenshan fault (LMSF), which is parallel to eastern Tibet and the Sichuan Basin in the northeast and southwest directions (Figure 1) and is about 500 km long and 40–50 km wide [17]. According to Liu et al. [12], the Wenchuan-Maowen, Yingxiu-Beichuan, and Guanxian-Jiangyou faults, as well as a number of thrust faults produced by the compression of the eastern Tibetan Plateau and the Yangtze Craton, dominate the LMSF zone. Only one significant seismic event ( M w = 6.1) occurred in the LMSF zone before the Wenchuan earthquake, and that was in 1989 [17]. The LMSF zone was dormant until 2008 [12]. The 290 km long segment of the LMSF that was ruptured by the Wenchuan earthquake propagated independently 270 km to the NE and 20 km in the SW direction. The co-seismic surface rupture was 80 km. According to Liu et al. [12], the Wenchuan earthquake was caused by a shift in the LMSF’s dip angle (30°–50° SW to 60°–70° NE) and fault motion (SW thrusting motion to NE strike–slip motion). Within 7 days following the mainshock, a series of aftershocks occurred, increased by 5.3 and 10 times more than the standard and relocated catalogues, respectively, [18]. The Wenchuan earthquake disaster resulted in 69,225 fatalities, 374,640 injuries, 17,939 unaccounted-for persons, and over 5 million displaced people.
It is important to emphasise that earthquakes of this magnitude and intensity trigger very important primary and secondary effects, some of which may contribute to the examination of important warning elements [19,20]. Primary effects, such as the surface fault of the Wenchuan earthquake, are associated with changes in concentrations of He, H2, CO2, CH4, O2, N2, Rn, and Hg in soil gas [21]; total electron content anomalies [22]; and well water levels [17]. Primary and secondary effects comprise the need for the introduction of an alternative ESI-07 intensity scale [19,23] and rupture imaging using teleseismic P waves [14].

2.1.2. Measurement Setup

Over the past 50 years, China has developed and established a structured network to track radon levels in groundwater for seismic studies. Almost all of China’ s provinces are represented in this network, which is made up of stations with a variety of instruments. The network is operated, maintained, and expanded with financial support from the China Earthquake Administration. The supplied radon series may be hourly (high sampling rate—HSR) or daily (low sampling rate—LSR), depending on the sensors deployed at each station. Table 1 lists the LSR stations as Ganzi, Mingshan, Panzhihua, and Sonpan, whereas the only HSR station is Kangding.
To measure the concentration of radon in groundwater, two instruments are used at China’s stations. Groundwater is sampled and degassed via bubbling at LSR stations. Then, using a specialised apparatus with the code FD-125, it is directed into an ionisation chamber or ZnS (Ag) detector where the concentration is determined by ionisations or scintillations. For the HSR Kangding station, groundwater is forced via a de-gasser and a gas-collecting device into a ZnS (Ag) detector, where radon is detected by a specialised instrument coded SD-3A with an hourly sample interval [24]. The accuracy for HSR and LSR measurements is 0.1 Bq/L.
As reported elsewhere [21,25,26,27,28], the recorded radon in groundwater concentrations is consistent with the crustal deformative process and/or with the stress diffusion procedure at depths that are able to squeeze deep-seated geofluids towards the surface. Therefore, the related time series include useful information to investigate their potential pre-seismic behaviour.

3. Mathematical Aspects

3.1. Fractal and Long Memory

Numerous natural physical systems can be explained by fractals. Usual fractal behaviour is often seen when the whole system, or a part of it, is translated, rotated, or stretched in space. Depending on the mathematical description of the changes, the system is characterised either as self-affine or self-similar. Self-affine and self-similar natural systems are fractals in the consensus that any component of them is a little or big copy of the total, but at various scales. As a result, fractal systems may be investigated by targeting their scaled parts. Additionally, the scaling properties of a fractal system are closely associated with the long memory [29,30,31,32] and complexity of the system [30,33,34,35,36]. Because fractal behaviour, long memory, and complexity are interconnected concepts, examining one of them will usually result in the other. For example, examining a system’s complexity will help to determine its long memory and vice versa [1,37,38]. All of these characteristics can show whether a system’s past, present, and future are strongly connected.
The direct methods are effective for calculating a system’s fractal features [16,37]. The techniques of Katz, Higuchi, and Sevcik are used in this paper because they offer direct estimates for fractal dimension computations [1,39]. Power-law dependencies are also present in fractal systems with long memory and complexity. DFA and power-law analysis are useful tools for the delineation of the related connections [1,37]. These techniques are used in this paper, as well, because of this. By the use of the associated Hurst exponent, all approaches may be compared. These techniques are going to be thoroughly explained in the sections that follow. The Hurst exponent is introduced first, followed by a robust DFA, techniques for the direct calculation of fractal dimensions, and, finally, power-law analysis.

3.2. Hurst Exponent

A metric known as the Hurst exponent (H) may be used to identify long-lasting connections in both time and space [40,41]. The time evolution of fractal phenomena, as well as the roughness of the related time series, can be identified using the Hurst exponent [1,42,43]. Various research topics have been investigated with the use of Hurst exponent, such as hydrology [40,41], astrophysics and applications [44,45], processes of capital markets [46,47,48,49], noisy observations of traces in traffic [50,51,52], seizures prior to epilepsy [53,54,55], climatic dynamics [56], and precursory time series before impending earthquakes [1,7,8,57].
The Hurst exponent value offers important details about the time series [1,37,39,40,41,58,59,60,61]:
(i)
If 0.5 < H 1 , then the series has a positive long-range autocorrelation. A series’ high value is followed by a series’ low value, and vice versa. High Hurst exponents suggest persistent interactions that are predicted to occur in the series’ far future;
(ii)
If 0 H < 0.5 , then the low values follow high values in the time series, and vice versa. There is an ongoing exchange between low and high values for low H values in the time series’ future (this is known as anti-persistency);
(iii)
If H = 0.5 associated processes are random, then the time series are totally uncorrelated.

3.3. Detrended Fluctuation Analysis (DFA)

Long-range power-law connections as well as erratic oscillations observed in time series appear prior to earthquakes [1,15,62,63]. DFA is a reliable method for spotting long-range power-law relationships in noisy, brief, non-stationary signals [1,64]. DFA has been successfully employed in different scientific domains such as the study of changes in the weather and climate [65,66,67,68], DNA sequences [69,70], heart dynamics [71,72,73,74], urban air pollution [37,39,61], pre-earthquake recordings of radon in soil [25,26], and electromagnetic variations in ULF, kHz, and MHz ranges [1,75,76,77,78,79].
Theoretically, DFA can show if a temporal signal has concealed long-range linkages that result in a self-similar process. Calculating the scaling exponent of the integrated time series allows one to discover these long-term relationships in the initial time series [1,39,57,66,70,71,80,81,82,83].

3.3.1. Application of DFA

The initial time signal is first integrated. The integrated signal’s fluctuations, F n , are then identified within a window of size n. The integrated time series’ scaling exponent (self-similarity parameter), α , is then calculated by fitting the linear l o g F n l o g n transformation via least squares. The l o g F n l o g n line may display one crossover at a scale n, where the slope exhibits an abrupt change, two crossovers at two different scales, n 1 and n 2 [57], or not even show a crossover at all, depending on the system dynamics.
The DFA of a one-dimensional temporal signal, y i ( i = 1 , , N ), can be implemented by the following procedure [1,39,57]:
(i)
The initial time series is, first, integrated:
y k = i = 1 k y i y
The entire average value of the time series is denoted in Equation (1) by the symbol <…>, and k stands for the various time scales.
(ii)
The integrated time series, y k , is then separated into equal, non-overlapping bins of length n.
(iii)
The function y k that represents the trend in the bin is then fitted. Simple linear trends or polynomials of second order or higher order may be used. Here, the linear function is used. This linear function’s y coordinate is denoted by the notation y n ( k ) in each box n.
(iv)
The local linear trend, y n ( k ) , is then subtracted from the integrated time series, y ( k ) , which is detrended in each box of length n. The detrended time series, y d n ( k ) , is determined in this manner and for each bin as follows:
y d n ( k ) = y ( k ) y n ( k )
(v)
The integrated and detrended time series’ fluctuations’ root mean square (rms) is then computed for each bin of size n as follows:
F ( n ) = 1 N k = 1 N y k y d n ( k ) 2
where F ( n ) are the rms fluctuations of the detrended time series, y d n ( k ) .
(vi)
For various sizes ( n ) of the scale boxes, the method steps (i)–(v) are repeated. This reveals the specific sort of connection between F ( n ) and n. If there are long-term relationships in the time series, then F ( n ) and n have an exponential relationship.
F ( n ) n α
The DFA scaling exponent α of Equation (4) assesses the strength of the time series’ long-term relationships.
(vii)
A linear association between l o g F ( n ) and l o g ( n ) is found via the logarithmic transformation of Equation (4). A strong linear connection suggests that the accompanying variations are long-lasting and, consequently, have a long memory. The square of Spearman’s ( r 2 ) is used in this paper to measure the accuracy of the linear fit. Good linear fits are defined as having r 2 0.95 or above [1,15,39,57,84].

3.3.2. Sliding Window DFA

The following six-step procedure was followed in order to implement the sliding window DFA:
(a)
The time series was segmented into windows of 64 samples. This segmentation approximately yields a two-month series’ part for the PZHS, SPS, GS, and MSS LSR stations, which record one measurement per day. The 64-sample window was also employed for fractal analysis of the data from three monitoring stations of urban air pollution with precisely the same measurement recording rate, namely, one measurement per day [39]. In a recent paper for the PZHS, a 256 segmentation DFA was employed [25], whereas, for radon in soil measurements, an approach of 128-sample window was utilised [85]. Nevertheless, since the windows are shifted 1 sample forward (sliding window technique), the whole signal is analysed, except from a 64-sample window, which was the final one. On the other hand, the 64-sample windowing yields a 64 h window for the HSR station of KDS, i.e., an analysis of about 2.5 days. Despite this differentiation, it is noteworthy that for a radon station in Pakistan, with the same recording rate as the one for KDS, a 64-window analysis was also utilised [16]. DFA from the data of KDS was analysed with 64 sample windows for consistency.
(b)
Every window was fitted using the least-squares fit of l o g F ( n ) vs. l o g ( n ) in accordance with Equation (4). The data were fitted to a straight line without seeking cross-overs, as in the related literature [1,25,39], with the restriction that the slope of the fit displays a square of Spearman’s correlation coefficient above or equal to 0.95.
(c)
The window was advanced by one sample, and the steps (a) and (b) were repeated until the signal’s end.
(d)
DFA slopes, α , were plotted against time, and the associated plot data were exported to ASCII output files for further use.

3.4. Fractal Dimension Analysis

3.4.1. Katz’s Method

To determine the fractal dimension, D, according to Katz’s method, the transpose array, [ s 1 , s 2 , , s N ] , of the series, s i , i = 1 , 2 , , N , is first determined, where s i = ( t i , y i ) and y i are the measured series values at the time instances, t i [86,87].
The value pairs ( t i , y i ) and ( t i + 1 , y i + 1 ) correspond to the two following points of the time series ( s i and s i + 1 ), for which the Euclidean distance is:
d i s t ( s i , s i + 1 ) = t i 2 t i + 1 2 + y i 2 y i + 1 2
The distances in Equation (7) add up in a curve, the total length of which is as follows:
L = i = 1 i = N d i s t ( s i , s i + 1 )
This curve will stretch in the planar to d, if it does not cross itself, where d is represented as follows:
d = m a x ( d i s t ( s i , s i + 1 ) ) , i = 2 , 3 , , N
By combining Equations (5)–(7), the Katz fractal dimension, D, becomes the following:
D = l o g ( n ) l o g ( n ) + l o g ( d / L )
where n = L / a ¯ and a ¯ is the average value of the distances of the points.

3.4.2. Higuchi’s Method

To calculate the fractal dimension, D, of a time series
y ( 1 ) , y ( 2 ) , y ( 3 ) , , y ( N )
recorded at intervals i = 1 , 2 N , a new sequence, y m k , is constructed as follows [88,89,90]:
y m k : y ( m ) , y ( m + k ) , y ( m + 2 k ) , , y ( m + N m k k )
The length of the curve associated with the time series is given by [88]:
L m ( k ) = 1 k i = 1 N m k y ( m + i k ) y ( m + ( i 1 ) k ) N 1 N m k k
In both equations, m and k are integers that specify the time interval between the series’ samples and are connected by the formula m = 1 , 2 k , where is the Gauss notation, namely, the bigger integer part of the included value.
By inserting the normalisation factor
N 1 N m k k
the average value, L ( k ) , of the lengths of Equation (13) exhibits a power law of the following form:
L ( k ) k D
The Higuchi’ s fractal dimension, D, is finally calculated by the slope of the linear regression of the logarithmic transformation of L ( k ) versus k, where k = 1 , 2 , , k m a x . It must be noted that the time intervals are k = 1 , . . , k m a x for k m a x 4 , i.e., k = 1 , 2 , 3 , 4 for k m a x = 4 and k = 2 ( j 1 ) / 4 , where j = 11 , 12 , 13 , for k > 4 ( k m a x > 4 ). Again, is the Gauss notation [87].

3.4.3. Sevcik’s Method

The fractal dimension of a time series according to the method of Sevcik [91] is approximated from the Hausdorff dimension, D h , as follows [87]:
D h = lim ϵ 0 l o g ( N ( ϵ ) ) l o g ( ϵ )
where N ( ϵ ) is the number of segments of length ϵ that add up to a curve that is associated with the time series. If the length of the curve is L, then N ( ϵ ) = L / 2 ϵ [87] and D h is the following:
D h = lim ϵ 0 l o g ( L ) l o g ( 2 ϵ ) l o g ( ϵ )
By applying a linear transformation twice, the N points of the curve, L, can correspond to a unit square of N × N cells of the normalised metric space. With this transformation, Equation (15) provides the Sevcik’s fractal dimension [87,91]:
D h = lim N 1 + l o g ( L ) l o g ( 2 ϵ ) l o g ( 2 ( N 1 ) )
The calculation improves as N .

3.4.4. Computational Methodology of Fractal Dimension

The next methodology was followed to calculate the fractal dimensions:
(i)
As in Section 3.3.2, the time series was segmented into windows of 64 samples. As mentioned, this segmentation approximately corresponds to, for the LSR stations (PZHS, SPS, GS, and MSS), a two-month signal. The 64-sample windowing was also employed in the fractal dimension calculation (with the same methods) from the data of the three LSR urban pollution stations with identical rates of measurements, i.e., one measurement per day [39]. As also mentioned in Section 3.3.2, for the HSR station KDS, the 64-sample segmentation corresponds to approximately 2.5 days. In a previous fractal dimension analysis (with the same methods), a 256-window approach was implemented for the PZHS [25]; however, in a very recent fractal dimension analysis with an identical methodology for an HSR radon station in Pakistan with the same rate of measurements as the one for KDS (one measurement per hour), a 64-window approach was utilised [16]. Finally, as in Section 3.3.2 and for consistency with the windowing of the other stations, a 64-sample window was chosen here as well for the KDS station.
(ii)
The fractal dimensions of each method were calculated as follows:
  • For the Katz’s method: The fractal dimension is the D of Equation (8) for n = 64 and a ¯ = 1 collected sample per measurement interval (1 day for PZHS, SPS, GS, and MSS and 1 h for the KDS) since a ¯ corresponds to the distance between the points of the series that constitute the parameter L [1,16,39].
  • Higuchi’s method: Equal to the slope, D, of the first-order least-squares fit of the logarithmic transformation of Equation (13), namely, the relation of l o g ( L ( k ) ) versus l o g ( k ) , for k m a x = 16 . In the aforementioned analysis for the urban air pollution stations [39], k m a x = 4 , whereas in the analysis of radon in Pakistan [16] and of the electromagnetic disturbances of the Ileia station, Greece, the approach k m a x = 16 was used. Based on the two latter papers, k m a x = 16 was also selected here.
  • Sevcik’s method: Equal to the Hausdorff dimension of Equation (16) ( D = D h ) for N = 64, namely, equal to the number of samples in each window, which constitutes parameter L.
(iii)
Each window was forwarded one sample (sliding window technique), and procedures (i)–(ii) were iterated until the end of the time series.
(iv)
Time evolution plots of the fractal dimensions in accordance with the Katz’s, Higuchi’s, and Sevcik’s methods were generated, and the partial data were extracted to ASCII files for further use.

3.5. Power-Law Analysis

The fractal power-law approach is another robust technique to identify the hidden long-lasting trends that are connected with the long-term links between space and time and are detectable before earthquakes [1,15,57,76,84,92,93,94,95,96]. As with all fractal-based techniques, power-law analysis describes the main core of fractality, namely, the existence of a power-law. Another reason is that the earthquake-generating systems progress gradually to self-organised critical (SOC) states exhibiting fractal evolution in space and time [95]. There have been approaches of fractal power-law analysis based on the Fourier transform [95,96]. The advantageous use of wavelets instead of the Fourier transform has been pointed out in several publications (e.g., [75,97,98,99]).
A time series’ power spectral density, S ( f ) , will follow a power law, if the series is a temporal fractal.
S ( f ) = a · f β
In Equation (17), β is the power-law exponent that quantifies the strength of the power-law connection, a is the amplification of the spectral density, and f is the frequency of a transform. According to several publications, this transform was chosen to be the wavelet one based on the Morlet bases and, specifically, f to be the central frequency of the Morlet wavelet [1,15,75,84,92,97,98,99].
The logarithmic transform of Equation (17) gives the following:
log S ( f ) = log a + β · log f
Since Equation (18) is a straight line, the values of β and a can be determined by fitting the corresponding data with the least-squares method. As with Section 3.3.2, the goodness of fit of the least-squares fit is quantified by the square of Spearman’s ( r 2 ) coefficient under the constraint r 2 0.95. The technique has been also described in other publications as spectral fractal analysis or spectral power-law fractal analysis. Hereafter, the phrase “power-law analysis” will be used.

Computational Methodology of Power-Law Analysis

Following the logic of Section 3.3.2 and Section 3.4.4, the next steps were followed to implement the power-law analysis:
(a)
The time series was separated into 128-sample windows. This is a double window size in comparison to the other two methods. This is performed because power-law analysis does not work well with small-sized windows. For the LSR stations (PZHS, SPS, GS, and MSS), this segmentation corresponds to a 4-month signal and, roughly, a 5-day signal for the KDS. In previous publications, a 128-sample window was employed in β parameter estimations [85] for recordings of similar recording rates, whilst in others, a 512-sample window [99] with a recording rate of one measurement every 10 min was employed.
(b)
The power spectrum, S ( f ) , based on the Morlet wavelet, as well as the central Morlet frequency, f, were calculated in each segment.
(c)
The parameters log S ( f ) and log f were fitted via least squares. Exponents, β , and power amplification, α , were computed for every window under the constraint that Spearman’s r 2 0.95 .
(d)
Steps (a) through (c) were iterated to the end of the time series. At each iteration, the window was shifted one sample forwards. As with the other techniques, the whole time series was covered except the last window.
(e)
The β and log a data were tabulated and saved in ASCII format for further use.

3.6. Further Issues

3.6.1. Formation of Analysis Classes

Two classes are formed to further organise the analysis results:
(a)
Class I: This class comprises windows that are associated with DFA least-squares log–log fits with a Spearman’s coefficient of r 2 0.95 and, simultaneously, a scaling exponent between 1 < α < 2 , i.e., modelled by the fBm class [84].
Class I segments:
  • With anti-persistency ( 1.35 < α < 1.5 )–persistency ( 1.5 α < 2 ) changes, some are of precursory worth [15,57,84,92,99].
  • With persistent behaviour ( 1.5 α < 2 ), they are characterised by others as footprints of impending seismic activity (e.g., [7,8] and the references therein).
(b)
Class II: This class includes time series segments with DFA’s r 2 < 0.95 (i.e., they do not adhere to the prominent fBm class) or 0 < α < 1 (i.e., they adhere to the fractional Gaussian noise (fGn) class).
Significantly, Class II segments:
  • They are of low precursory worth and low predictability [1,15,57,84,85,92,99].
  • They are complements of Class I segments.

3.6.2. Comparisons of the Fractal Results

As shown in previous papers [1,57,84,85], the best approach to compare the results of all fractal methods is through the Hurst exponent.
For Class I segments, the Hurst exponent (H) is calculated as follows (e.g., [1,84] and the references therein):
(1)
From DFA’s α exponent as:
H = α 1
(2)
From fractal dimension (D) as:
H = 2 D
(Berry’s equation)
(3)
From power-law β as:
H = 0.5 · ( β 1 )
It should be emphasised that deviations from the straightforward linear connection of Equations (19)–(21) can be seen in the in situ data [1,15,57,84]. The relationship between the fractal analysis parameters remains linear, possibly of a slightly different form, as indicated in the aforementioned works.

3.7. Meta-Analysis

The so-called meta-analysis [1,61] is implemented by combining the outcomes from the ASCII files of all five methods, namely, DFA; Higuchi’s, Katz’s, and Sevcik’s fractal dimensions; and power-law analysis. A two-step process is followed:
(a)
(Step 1): According to user-defined thresholds, each ASCII output results file is computationally scanned for out-of-threshold values. The ASCII files carrying the fractal dimension values are searched for under threshold values, whilst the ASCII files containing DFA’s exponents and the power-law b e t a values are searched for over threshold values. New ASCII step 1 files are generated that contain the out-of-threshold values.
(b)
(Step 2): Under the restriction that each segment’s first sample date is arbitrarily considered as the date of the whole segment, the step 1 ASCII files are computationally filtered to find areas with common dates. The above computational process results in the full coverage of all dates, except the one of the last window. The whole procedure is iterated over the results of all possible combinations of the following:
  • DFA versus fractal analysis or versus at least two fractal dimension calculation techniques (six combinations);
  • Fractal analysis versus at least two fractal dimension calculation techniques (four combinations);
  • One fractal dimension calculation technique versus the other two (three combinations).
Through the above repetitive process, 13 unique combinations of techniques per 5, 4, 3, and 2 are produced. This is very important because it is practically equivalent to the coupling of different mono-fractal methods, and this fact provides a synthetic view of the results of the fractal techniques, increasing the scientific evidence regarding the underlying nature of the identified fractal disturbances. This has been pointed out in recent publications [1,15,37,39,61,84].

4. Results and Discussion

Figure 2, Figure 3, Figure 4, Figure 5 and Figure 6 present the variation in the DFA exponent, α , over time with respect to the evolution of the associated square of the Spearman’s correlation coefficient versus the measured disturbances of radon in groundwater. As can be observed from these figures, the DFA scaling exponent profile is entirely distinct from the one of the time series. This has been noted in earlier works as well [1,15,25,27,37,39,57,84]. The reason relies on the fact that DFA manages to locate effectively hidden forms in time series, even in non-stationary ones [62,64,100]. Numerous α exponents are within the Class I range (Section 3.6.1). This means that the corresponding 64-sample windows are successive fBm ones ( r 2 0.95), and this has been acknowledged as a notable sign of pre-seismic activity [15,57,84,85,92,99].
The specific details of each station show interesting information. At first, Figure 2 is a noteworthy case of completely different DFA exponents and signal profiles for the KDS. By inspecting this figure, a first period can be observed starting from window 1 ( d a y 1) up to window 2100 (approximately d a y 93). During this period, the DFA scaling exponents, α , are enhanced with 1.5 α < 2.0 ( 0.5 H < 1.0 , Equation (19)), whereas, most significantly, several exponents are above 1.8 ( 0.8 H < 1.0 , Equation (19)). All of these are successive fBm segments since the corresponding Spearman’s coefficient in each window is r 2 0.95. These successive fBm segments exhibit persistence since the Hurst exponents are above 0.5, whilst a number of them show great persistent behaviour ( 0.8 H < 1.0 ). As mentioned in Section 3.6.1, the references therein, and those in this section, the successive fBm segments, and especially those with great persistent behaviour, correspond to radon in groundwater areas with a high potential of association with seismic activity. In addition, in Figure 3, a very interesting period can be observed between windows 50 and 100 (approximately between d a y 60 and d a y 121), with DFA exponents 1.5 α < 2.0 ( 0.5 H < 1.0 ) and several above 1.8 ( 0.8 H < 1.0 ) for the GS. There is a similar case for the MSS. A corresponding high DFA period ( 1.5 α < 2.0 , several exponents higher than 1.8) is between windows 70 ( d a y 85) and 115 ( d a y 139). What makes these three cases of great significance is that the above periods with high DFA exponents are rather synchronous. They are parallel in time, and there are similar periods in these three stations. Indeed, the period for KDS is from d a y 1 to d a y 93, the period for GS is between d a y 60 and d a y 121, and the period for MSS is from d a y 85 to d a y 115. Importantly, these periods correspond to strong persistency since several Hurst exponents are in a range of 0.8 H < 1.0 . They also correspond to the Class I category, exhibiting strong fBm behaviour. The significance of identifying strong, persistent fBm behaviour has been emphasised in several publications [1,15,39,84,97,99] as a footprint for ensuing earthquakes. The interpretations of these precursory footprints in these publications are based on an asperity model [101], according to which it is the roughness of fBm profiles, their relative motion, and their association with micro-crack branching and acceleration that explain why especially persistent fBm behaviour is a noteworthy sign of precursory activity during the preparation of earthquakes. These findings are very important and have to be emphasised, especially because the DFA results of the other two stations (PZHS and SPS) have almost all areas with low DFA exponents, mostly in the Class II category or at the lowest part of the Class I category. The latter fact shows, very clearly, that the PZHS and SPS disturbances are of low precursory value, if not totally insignificant. From the above DFA evidence, it can be highlighted, from a time perspective, that the important time periods identified (KDS from d a y 1 to d a y 93, GS from d a y 60 to d a y 121, and MSS frrom d a y 97 to d a y 130) comprise pre-earthquake signs of the great Wenchuan earthquake, which occurred on d a y 133. According to the literature (see reviews [2,8]), the KDS-GS-MSS time window is within the precursory window range of radon precursors. Particularly, when the high magnitude ( M w = 7.9) is accounted for, the above time window extends well before the earthquake’s occurrence (even from d a y 1).
The differentiation between the KDS, GS, and MSS and the PZHS and SPS becomes more important when the differentiations in the distance and the underlying geology are taken into consideration. Under the aspect of distance, a potential claim for the PZHS could be that its long epicentral distance (526.0 km) makes it difficult to show significant DFA disturbances. This claim, however, is unsupported when it is considered that the PZHS has shown previous precursory DFA variations for earthquakes with epicentral distances above 533 km [25]. In support, the GS shows precursory DFA variations despite being 325.5 km from the epicentre, especially when the SPS did not show precursory DFA results even though it was located closer (182.5 km) and at a comparable distance with the other two stations (KDS and MSS), which both showed significant DFA outcomes. This is an important observation that has to be outlined. As has been pointed out in the above references [1,15,39,84,97,99], the preparation area of earthquakes includes special preferable precursory paths. Moreover, the selectivity effect that has been proposed and utilised for precursory activity (e.g., [102,103,104,105]) suggests that during seismic preparation, there are selective paths that are followed by the disturbed activators. The selection of certain paths has been demonstrated for MHz electromagnetic variations [57,106] and radon precursors [57,107] and has also been expressed in reviews on the subject [2,5,7,8,11]. In this sense, the geological path from a station to the earthquake’s epicentre gains special meaning. Hence, it is very important that the GS and KDS are located on the same big fault. More importantly, the MSS is also on the prolongation of this GS-KDS path and is, surely, in line with the Wenchuan’s epicentre. This may be supported by the above information that the high elevation of the underlying geology, in association with the big fault on which the GS-KDS-MSS operate and the proximity of the KDS to the Wenchuan epicentre, may explain why these stations are very sensitive to recording radon disturbances with hidden traces of significant precursory values.
This is reinforced by the fact that the SPS-PHZS are on a completely different fault line, despite being on a common fault with the KDS. The fact (as mentioned above) that KDS is also in line with the other two stations (GS and MSS) and near the earthquake’s epicentre makes the recordings of this station even more important. The huge magnitude ( M w = 7.9) and the low depth (19 km) that resulted in a huge energy release, in association with the geologically sensitive background and the proximity of the KDS to the earthquake’s epicentre, may explain why post-activity was observed only in this station. Indeed, there is a great period from approximately window 3000 (approximately d a y 133) up to window 6100 (approximately d a y 271) in Figure 2 where the DFA scaling exponents, α , are much higher than 1.5, with several above 1.8 ( 0.8 H < 1.0 ). These DFA variations were identified just after the great Wenchuan earthquake ( d a y 133), and this is the first time that such post activity has been found. The reader should note here that there is also an extended period near sample 2000 (approximately d a y 88) and a shorter one near sample 8200 (approximately d a y 363) of completely non-successive ( r 2 < 0.95 ) segments. This latter period provided false DFA exponents above 2, and for this reason, they were cut off from the figure. The small period around 8200 is a Class II one. Except for these two non-precursory cases, there are also scattered Class II exponents of low-precursory value. These are non-successive ( r 2 < 0.95 ) fBm segments or fGn segments. Finally, with regards to the pre-seismic signs, from a geological perspective, these are in line with the literature findings since, according to the reviews for radon and electromagnetic precursors [2,4,5,6,7,8,9], the epicentre’s distance of precursory activity of KDS-GS-MSS is in accordance with the reported ones. In conclusion, based on the visual observations from the DFA results, the reader should refer to the following points: (a) DFA outlines hidden trends not (usually) observed from the signal; (b) the time window of the identified activity is sufficient to accept it as precursory; (c) the distance and the pathway provide some explanations for the identification in three out of five stations; and (d) most importantly, even for the robust DFA method, visual observation is not enough, and a combination of different methods is necessary in order to find the precursory time periods with enhanced importance and combined evidence [1,37,39,61]. All of these will be presented later in the text, together with the combined evidence.
Figure 7, Figure 8, Figure 9, Figure 10 and Figure 11 demonstrate the temporal evolution of the fractal dimensions calculated using Katz’s, Higuchi’s, and Sevcik’s methods. There are noticeable differences found. The three methods’ computed values for the fractal dimension show variations as well. All discrepancies are due to the varied calculation techniques used by the three fractal dimension techniques. This has been acknowledged in recent publications [1,15,39]. As pointed out, the methods of Katz and Higuchi estimate higher fractal dimensions than the ones of Sevcik. As general trends, the Katz’s fractal dimensions are around 1 and 1.1 ( 0.9 < H < 1 , Equation (20)), the ones of Higuchi’s method are roughly between 1.5 and 2 ( 0 < H < 0.5 , Equation (20)), and those of Sevcik’s method are, approximately, between 1 and 1.5 ( 0.5 < H < 1.0 , Equation (20)). However, what is important is not the value range but the abrupt decrease in the calculated fractal dimensions. For this reason, the details of every sub-figure in Figure 7, Figure 8, Figure 9, Figure 10 and Figure 11 have great importance, especially in association with the results of the DFA method.
In reference to the KDS and in the first 200 windows (from d a y 1 to d a y 8) of Figure 7, a noteworthy decrease can be observed in Sevcik’s fractal dimension between 1 and 1.4 ( 0.6 < H < 1.0 ) and in Higuchi’s D between 1 and 2 ( 0 < H < 1.0 ). Katz’s method estimates, unexpectedly, higher fractal dimensions between 1.0 and 1.2( 0 < H < 0.8 ) in this area. Simultaneously, there is an increase in groundwater radon. This is a rare, serendipitous finding that has been reported in other publications [2,8]. Abrupt drops in Higuchi’s fractal dimension are observed, next, between windows 1900 ( d a y 84) and 2200 ( d a y 97). The drops in Sevcik’s fractal dimension are between windows 1800 ( d a y 79) and 2000 ( d a y 88). The Sevcik’s fractal dimensions between windows 2000 and 2200 are above 2, and for this reason, they were cut off, being considered errors in calculations. Once again, the Katz’s fractal dimensions are, peculiarly, higher around window 1800. With regards to the GS, a significant decrease in Higuchi’s fractal dimension can be spotted in Figure 8 between window 40 ( d a y 48) and window 110 ( d a y 133). The decrease in Sevcik’s fractal dimension is roughly synchronous and of the same profile as the one of Higuchi’s, although milder and with a smaller duration, between windows 60 ( d a y 73) and 110 ( d a y 133). The fractal dimensions of the MSS (Figure 9) also exhibit decreased profiles through the calculations with Higuchi’s and Sevcik’s methods. These decreases in D values can be observed between windows 70 ( d a y 85) and 140 ( d a y 170). In summary, here, the key periods from the fractal dimension calculations for the KDS are between d a y s 1 and 8 and between d a y s 84 and 97. For the GS, the period is between d a y 48 and d a y 110, and for the MSS, the period is between d a y 85 and d a y 170. Since the Wenchuan earthquake occurred on d a y 133, the fractal dimension variations of KDS and GS can be considered, most probably, pre-seismic. The same is valid for the d a y 85– d a y 133 variations of the MSS. These findings reinforce the claims expressed in the discussion of DFA above since the fractal dimension calculation techniques are completely different from the ones of DFA. Moreover, two fractal techniques show these tendencies and, significantly, in comparable time intervals. These facts further support the necessity of using different fractal techniques in parallel. This has been emphasised in recent publications [1,15,61,84]. However, as with the DFA outcomes of the KDS, there is a wide period between windows 3000 ( d a y 133) and 6100 ( d a y 271) during which the fractal dimensions from Higuchi’s and Sevcik’s methods exhibit very significant variations. This period is the same as the corresponding one discussed for the DFA results and refers to post-seismic variations. Additional post-seismic variations are addressed here via Higuchi’s and Sevcik’s fractal dimension from the MSS between d a y 133 and d a y 170. The reader should note here that as with the outcomes of DFA, SPSS, as shown in Figure 11, does not show certain D patterns. Despite the fact that DFA did not provide trends in DFA profiles for the PZHS, the corresponding Higuchi’s fractal dimension shows a small decrease between windows 200 and 250.
Figure 12, Figure 13, Figure 14, Figure 15 and Figure 16 show the results from the power-law method. This method is one of the wider used techniques and has been considered as one of the most powerful to identify the hidden patterns in time series [43,75,76,77,92,95,98,99,106,108,109,110,111,112,113,114,115]. Hence, when certain trends are found with the power-law method, there is strong evidence for the underlying long memory of the associated geosystem. At first, as with the outcomes of the DFA and fractal dimension methods, the time evolution of the power-law exponent, β , differs from the one of the time series. However, in order to discuss the interesting results of each figure, the following information should be taken into consideration for the successive fractal segments:
  • If 1.0 < β 3.0 , then the associated time series is a temporal fractal and follows the Class I category:
    • If 1.0 < β < 2.0 , then the time series follows anti-persistent paths;
    • If 2.0 < β < 3.0 , then the time series follows persistent paths.
  • If 1.0 β < 1.0 , then the time series is of low predictability and follows the Class II category:
    • If β = 1.0 , then the fluctuations in the related processes are not growing and, hence, a stationary system describes the series;
    • If β = 2.0 , then the underlying dynamics are random and the related system has no memory.
In Figure 12, for the KDS, as with DFA and the fractal dimension calculation techniques, a first period is observed up to window 2200 (approximately d a y 120) with scattered successive fractal windows with 1.7 < b < 2.2 . These fractal epochs correspond to predictable Class I segments with an interchange between persistency and anti-persistency. According to several publications (e.g., [99]), this is a sign of precursory activity. Interestingly, this epoch is almost identical to those identified with DFA and the three fractal dimension techniques, and this is very important. Figure 13 for the GS shows a first period between windows 1 and 40 ( d a y 1 to d a y 58) and windows 50 and 90 ( d a y 73 to d a y 131). Both these periods have areas with b > 2 (Class I, persistent) and can be considered as precursory of the great Wenchuan earthquake. The last period at the end of the analysed windows shows b < 1.6 and, as with DFA and fractal dimension techniques, has low predictability and is, hence, of low precursory value. This combined finding provides more evidence. In Figure 14, for the MSS, there are two periods, between window 30 ( d a y 43) and window 60 ( d a y 87) and between window 70 ( d a y 102) and window 110 ( d a y 160). Although these periods are anti-persistent, they are within the middle part of the predictable value range of the Class I category. Since the periods match with those of the other techniques, there is a probability that they might be signs (pre and post) of the great Wenchuan earthquake. As with the DFA and fractal dimension techniques, there is a great post-seismic region within the same period. It is very important that the findings of the techniques match even for the two other stations, the PZHS and SPS. Both showed no successive fractal window. This fact reinforces the findings of the DFA and Higuchi’s and Sevcik’s methods for the MSS.
Evaluating the results presented so far, there is a period up to d a y 133 (the day of occurrence of the great Wenchuan earthquake) that can be systematically identified in all fractal analysis results, that is, from all methods for the KDS, GS, and MSS, whereas the PSZH and SPS do not show any such noteworthy period. These fractal epochs are discussed and considered pre-seismic signatures of the great Wenchuan earthquake. The recognition of signs in the KDS, GS, and MSS and the lack of systematic signs in the PSZH and SPS can be attributed to the different geological paths, in association with the theory of asperities and the selectivity effects in earthquake-related research. It is also extensively discussed that finding periods of over- or under-threshold fractal values is significant, but most significant is the synthetic finding of common over- or under-threshold areas with different techniques (meta-analysis, Section 3.7). When, and if, such common locations are discovered, the scientific arguments for the existence of a seismic warnings concealed in the time series increases, making these claims stronger. It is, therefore, very important to apply meta-analysis to the fractal results of the KDS, GS, and MSS so as to enhance the presented evidence. Meta-analysis for the PSZH and SPS results is unnecessary since the corresponding signs are not enough. On the other hand, an earthquake as great as the Wenchuan is expected to be associated with other earthquakes. Moreover, other earthquakes in the overall area might also explain the warnings, and there is a possibility that the presumed post-seismic signs might also be the pre-seismic activity of other earthquakes in the area, and this is the first time such a view is mentioned in this paper. For this purpose, the data on earthquakes of 2008 with M w 5.5 were accessed from USGS over an area greater than the inserted one in Figure 1 [116]. The latter reference presents online the selected area together with data on 19 earthquakes and a map representing it on Chinese terrain. The data on these 19 earthquakes are presented in Table 2, whereas Figure 17 presents the earthquakes in a Google Earth Map after creating the corresponding map kml file from USGS and importing it to Google Earth. In this way, Figure 17 presents the earthquakes alternatively in comparison to the USGS map [116].
In the consensus expressed above, the final step of the related analysis should include the following: (a) the above 19 earthquakes (Wenchuan included); (b) the fractal results and the discussion on threshold setting presented in Figure 2, Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9, Figure 10, Figure 11, Figure 12, Figure 13, Figure 14, Figure 15 and Figure 16; and (c) the logic of connecting the fractal methods expressed in Section 3.7. To achieve this, Figure 18 presents combined plots for each station. This multi-figure conceals all the important findings. Similar representation has been adopted in other publications as well [1,25,37,39]. The subplots are mixed and have several symbols. For this reason it is crucial to delineate the significant information given:
(a)
All over- or under-threshold results of all fractal methods (step 1, meta-analysis) for the KDS, MSS, and GS. The threshold results of each station are combined per 2, 3, 4, and 5 methods (step 2, meta-analysis; a total of 13 combinations) versus all 19 earthquakes of Table 2 and Figure 17.
As mentioned in the previous paragraph, it is important not only to identify footprints using one or more techniques (already conducted here) but more importantly to link the different techniques focusing on similar aspects of the problem at hand. To achieve this:
(1)
The exact over- or under-threshold dates were located computationally from the fractal outputs of each station (step 1, meta-analysis). These dates are year, month, day, and hour for the HSR and KDS, and year, month, and day for the MSS and GS. This is conducted through a serial search.
(2)
The common threshold dates from all different techniques were found through an incremental computational search. The outputs used are from all methods and with up to 13 different combinations of these. All outputs were generated through special software and were stored in a computer for use.
(3)
The earthquake data from USGS [116] were transformed into an adequate ASCII file for the generation of the final plot.
(b)
Wherever the symbols of different methods coincide in time, this means that the signs of seismicity are provided by more than one method. If all 13 methods coincide, this means that the evidence is maximised. The more the techniques point to similar findings, the more rigid the evidence is. It should be emphasised to the reader that this coinciding is conducted on the step 1 results, that is, on the fractal outputs.
In the above sense and starting from the combined findings of lesser importance, it can be observed from sub-figure a of Figure 18 that there are several windows of the combination of DFA versus Higuchi’s and Sevcik’s methods (⊡, the plot needs to be zoomed to show the details) that are in the period of all 19 earthquakes. The fact that there are three techniques that show this behaviour makes these fractal disturbances noteworthy. However, it is not possible to discriminate whether these disturbances are post-seismic activity of the great Wenchuan earthquake or pre-seismic activity of another one in the area, and this is a limitation of the present methodology. On the other hand, observing the data of Table 2 carefully, it can be seen that earthquakes 13–18 are practically the seismic sequence of the Wenchuan earthquake since all occurred on the same day as earthquake 19. Especially 16–18 happened at the same hour. Moreover, earthquakes 11 and 12 are also within the post sequence of the Wenchuan earthquake because they occurred just one (11) and two (12) days after. Moreover, all earthquakes from 8 to 19 occurred in the same month (April 2008). These peculiar facts complicate discriminating the above concurrent findings using the three techniques as post- or pre-seismic. This is not the case for the three co-occurrences of the combinations of the fractal dimension methods ( p e n t a g o n —yellow, namely, Sevcik’s versus Katz’s, and Higuchi’s methods (three techniques); p e n t a g o n —blue, viz. Katz’s versus Higuchi’s and Sevcik’s methods (three techniques); and green ◁, namely, Higuchi’s versus Katz’s and Sevcik’s methods (three techniques)). These occurred prior to events 7, 6, and 5 but could be prior to 4, 3, 2, and 1 as well (decreasing time distance between fractal disturbance and earthquake occurrence) (all these sub-figures need to be zoomed as well).
The most significant finding of this paper is left for the last section. It can be clearly observed in all sub-figures that the vast majority of coincidences are prior to the great Wenchuan earthquake (19) and its synchronous post-earthquakes (10–19). This is the most important observation that the reader should emphasise. It is not only one method but all 13 methods that coincide. Moreover, during this preparation phase of the great Wenchuan earthquake, there are far more coincidences for even three or fewer techniques as well (such as red ▷ in sub-figure b, namely, power-law analysis versus Katz’s and Sevcik’s methods (three techniques), and yellow ⊡, that is, DFA versus Higuchi’s and Sevcik’s methods (three techniques)). In the most emphatic manner, it is here declared that all these fractal disturbances (importantly, with results from the meta-analysis) are, most possibly, due to the great Wenchuan earthquake. This is the most important finding of this paper, which we leave for last to emphasise. It is the great magnitude of the Wenchuan earthquake, its post-seismic activity, and the novel combining–linking of several different fractal techniques that have managed to allow such an important finding to be outlined.
As final remarks of this paper, it should be noted that the overall approach of combining different fractal techniques by employing thresholds on the final ASCII outcomes is a significant novelty of this paper that should be emphasised. Moreover, according to the references of this work, there are only a few papers that report different methods in evidencing possible anomalies associated, most probably, with seismic activity. This is also a significant novelty. In addition, this approach has also worked well with air pollution P M 10 time series [37,39,61]. This latter fact makes the methodology and interpretations useful for other scientific areas as well and is very promising for the continuation of this work.

5. Conclusions

This study investigated the fractal patterns hidden in one-year radon in groundwater disturbances derived from five stations in China before and after the devastating Wenchuan ( M w = 7.9) shallow (depth = 19 km) earthquake that occurred on 12 May 2008 ( d a y 133). The data were analysed with five distinct fractal techniques (DFA; fractal dimensions with Higuchi’s, Katz’s, and Sevcik’s methods; and power-law analysis). Sliding windows of step 1 were utilised in segmented portions glided throughout each signal. Via literature-based thresholds, several notable areas were found in the fractal variations of the KDS, GS, and MSS data, whilst non-significant fractal portions were found in the signals of PZHS and SPS. Up to d a y 133 (12 May 2008), several critical epochs were located in the signals of KDS, GS, and MSS. Specifically, the DFA exponents during these epochs were 1.5 α < 2.0 ( 0.5 H < 1.0 ), whereas several exponents were above 1.8 ( 0.8 H < 1.0 ). The fractal dimension epochs exhibited Katz’s fractal dimensions around 1.0 and 1.2 ( 0.8 < H < 1 ), Higuchi’s dimensions between 1.5 and 2 ( 0 < H < 0.5 ), and Sevcik’s dimensions between 1 and 1.5 ( 0.5 < H < 1.0 ). Several power-law b exponents were above 1.7, and numerous were above 2.0. All these are in agreement with precursory fractional Brownian motion signal parts. The differentiations between the KDS-GS-MSS and the PZHS-SPS ones could be attributed to the geological background with the theories of asperities (fractional Brownian motion profile relative roughness) and selectivity. As a systematic last action, all results of KDS-GS-MSS were analysed using a novel, two-step, fully computerized methodology that located the exact out-of-threshold fractal areas and combined the outcomes of the different methods per 5, 4, 3, and 2 (maximum 13 common combinations) in association with the 19 earthquakes of M w 5.5 in the greater area. The vast majority of the different combinations of techniques showed coincidences prior to the great Wenchuan earthquake and its synchronous post-earthquakes. This important finding was justified not only with one method but also, in many cases, with all 13 different methods. Critical epochs after the Wenchuan earthquake could be attributed to other earthquakes in the area, whereas a post-seismic view can also be accepted. The combined results, in association with the great earthquake’s magnitude and small depth, make the findings of this paper promising for earthquake-related studies.

Author Contributions

Conceptualisation, D.N.; data curation, A.A. and D.N.; formal analysis, A.A., D.N. and N.W.; investigation, A.A., D.N. and N.W.; methodology, A.A. and D.N.; resources, A.A. and N.W.; software, D.N.; supervision, D.N. and N.W.; visualisation, A.A. and D.N.; writing—original draft, D.N.; writing—review and editing, A.A. and N.W. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Data Availability Statement

Not applicable.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Nikolopoulos, D.; Petraki, E.; Yannakopoulos, P.H.; Priniotakis, G.; Voyiatzis, I.; Cantzos, D. Long-lasting patterns in 3 kHz electromagnetic time series after the ML = 6.6 earthquake of 2018-10-25 near Zakynthos, Greece. Geosciences 2020, 10, 235. [Google Scholar] [CrossRef]
  2. Cicerone, R.; Ebel, J.; Britton, J. A systematic compilation of earthquake precursors. Tectonophysics 2009, 476, 371–396. [Google Scholar] [CrossRef]
  3. Hough, S. The Great Quake Debate: The Crusader, the Skeptic, and the Rise of Modern Seismology; University of Washington Press: Washigton, DC, USA, 2020. [Google Scholar]
  4. Hayakawa, M.; Hobara, Y. Current status of seismo-electromagnetics for short-term earthquake prediction. Geomat. Nat. Hazards Risk 2010, 1, 115–155. [Google Scholar] [CrossRef]
  5. Molchanov, O.A.; Hayakawa, M. Seismo-Electromagnetics and Related Phenomena: History and Latest Results; Number A8; Terrapub: Tokyo, Japan, 2008; p. 189. [Google Scholar]
  6. Ouzounov, D.; Pulinets, S.; Hattori, K.; Taylor, P. Pre-Earthquake Processes: A Multidisciplinary Approach to Earthquake Prediction Studies Posters; Wiley: Hoboken, NJ, USA, 2018; p. 384. [Google Scholar]
  7. Petraki, E.; Nikolopoulos, D.; Nomicos, C.; Stonham, J.; Cantzos, D.; Yannakopoulos, P.; Kottou, S. Electromagnetic Pre-earthquake Precursors: Mechanisms, Data and Models-A Review. J. Earth Sci. Clim. Chang. 2015, 6, 250. [Google Scholar] [CrossRef]
  8. Petraki, E.; Nikolopoulos, D.; Panagiotaras, D.; Cantzos, D.; Yannakopoulos, P.; Nomicos, C.; Stonham, J. Radon-222: A Potential Short-Term Earthquake Precursor. J. Earth Sci. Clim. Chang. 2015, 6, 282. [Google Scholar] [CrossRef]
  9. Uyeda, S.; Nagao, T.; Kamogawa, M. Short-term earthquake prediction: Current status of seismo-electromagnetics. Tectonophysics 2009, 470, 205–213. [Google Scholar] [CrossRef]
  10. Conti, L.; Picozza, P.; Sotgiu, A. A Critical Review of Ground Based Observations of Earthquake Precursors. Front. Earth Sci. 2021, 9, 676766. [Google Scholar] [CrossRef]
  11. Ghosh, D.; Deb, A.; Sengupta, R. Anomalous radon emission as precursor of earthquake. J. Appl. Geophys. 2009, 187, 245–258. [Google Scholar] [CrossRef]
  12. Liu, C.; Dong, P.; Zhu, B.; Shi, Y. Stress Shadow on the Southwest Portion of the Longmen Shan Fault Impacted the 2008 Wenchuan Earthquake Rupture. J. Geophys. Res. Solid Earth 2018, 123, 9963–9981. [Google Scholar] [CrossRef]
  13. Yin, Y.; Wang, F.; Sun, P. Landslide hazards triggered by the 2008 Wenchuan earthquake, Sichuan, China. Landslides 2009, 6, 139–152. [Google Scholar] [CrossRef]
  14. Xu, Y.; Koper, K.D.; Sufri, O.; Zhu, L.; Hutko, A.R. Rupture imaging of the Mw 7.9 12 May 2008 Wenchuan earthquake from back projection of teleseismic P waves. Geochem. Geophys. Geosystems 2009, 10, Q04006. [Google Scholar] [CrossRef]
  15. Nikolopoulos, D.; Matsoukas, C.; Yannakopoulos, P.H.; Petraki, E.; Cantzos, D.; Nomicos, C. Long-Memory and Fractal Trends in Variations of Environmental Radon in Soil: Results from Measurements in Lesvos Island in Greece. J. Earth Sci. Clim. Chang. 2018, 9, 1–11. [Google Scholar] [CrossRef]
  16. Rafique, M.; Iqbal, J.; Shah, S.A.A.; Alam, A.; Lone, K.J.; Barkat, A.; Shah, M.A.; Qureshi, S.A.; Nikolopoulos, D. On fractal dimensions of soil radon gas time series. J. Atmos. Sol.-Terr. Phys. 2022, 227, 105775. [Google Scholar] [CrossRef]
  17. Shi, Z.; Wang, G.; Wang, C.y.; Manga, M.; Liu, C. Comparison of hydrological responses to the Wenchuan and Lushan earthquakes. Earth Planet. Sci. Lett. 2014, 391, 193–200. [Google Scholar] [CrossRef]
  18. Yin, X.Z.; Chen, J.H.; Peng, Z.; Meng, X.; Liu, Q.Y.; Guo, B.; Li, S.C. Evolution and Distribution of the Early Aftershocks Following the 2008 Mw 7.9 Wenchuan Earthquake in Sichuan, China. J. Geophys. Res. Solid Earth 2018, 123, 7775–7790. [Google Scholar] [CrossRef]
  19. Audemard, F.; Azuma, T.; Baiocco, F.; Baize, S.; Blumetti, A.M.; Brustia, E.; Clague, J.; Comerci, V.; Esposito, E.; Guerrieri, L.; et al. Earthquake Environmental Effect for Seismic Hazard Assessment: The ESI Intensity Scale and the EEE Catalogue; ISPRA—Servizio Geologico d’Italia: Roma, Italy, 2015; Volume XCVII, pp. 5–8. [Google Scholar]
  20. Keller, E.A.; Pinter, N. Active Tectonics, Earthquakes, Uplift and Landscape, 2nd ed.; Prentice Hall: Upper Saddle River, NJ, USA, 2002; 362p. [Google Scholar]
  21. Zhou, X.; Du, J.; Chen, Z.; Cheng, J.; Tang, Y.; Yang, L.; Xie, C.; Cui, Y.; Liu, L.; Yi, L.; et al. Geochemistry of soil gas in the seismic fault zone produced by the Wenchuan Ms 8.0 earthquake, southwestern China. Geochem. Trans. 2010, 11, 5. [Google Scholar] [CrossRef] [PubMed]
  22. Liu, J.Y.; Chen, Y.I.; Chen, C.H.; Liu, C.Y.; Chen, C.Y.; Nishihashi, M.; Li, J.Z.; Xia, Y.Q.; Oyama, K.I.; Hattori, K.; et al. Seismoionospheric GPS total electron content anomalies observed before the 12 May 2008 Mw7.9 Wenchuan earthquake. J. Geophys. Res. Space Phys. 2009, 114. [Google Scholar] [CrossRef]
  23. Huayong, N.; Hua, G.; Yanchao, G.; Blumetti, A.M.; Comerci, V.; Di Manna, P.; Guerrieri, L.; Vittori, E. Comparison of Earthquake Environmental Effects and ESI intensities for recent seismic events in different tectonic settings: Sichuan (SW China) and Central Apennines (Italy). Eng. Geol. 2019, 258, 105149. [Google Scholar] [CrossRef]
  24. Ren, H.; Liu, Y.; Yang, D. A preliminary study of post-seismic effects of radon following the Ms 8.0 Wenchuan earthquake. Radiat. Meas. 2012, 47, 82–88. [Google Scholar] [CrossRef]
  25. Alam, A.; Wang, N.; Zhao, G.; Mehmood, T.; Nikolopoulos, D. Long-lasting patterns of radon in groundwater at Panzhihua, China: Results from DFA, fractal dimensions and residual radon concentration. Geochem. J. 2019, 53, 341–358. [Google Scholar] [CrossRef]
  26. Alam, A.; Wang, N.; Zhao, G.; Barkat, A. Implication of radon monitoring for earthquake surveillance using statistical techniques: A case study of Wenchuan earthquake. Geofluids 2020, 2020, 2429165. [Google Scholar] [CrossRef]
  27. Alam, A.; Wang, N.; Petraki, E.; Barkat, A.; Huang, F.; Shah, M.A.; Cantzos, D.; Priniotakis, G.; Yannakopoulos, P.H.; Papoutsidakis, M.; et al. Fluctuation Dynamics of Radon in Groundwater Prior to the Gansu Earthquake, China (22 July 2013: M s = 6.6): Investigation with DFA and MFDFA Methods. Pure Appl. Geophys. 2021, 178, 3375–3395. [Google Scholar] [CrossRef]
  28. Ma, T.; Wu, Z. Precursor-Like Anomalies prior to the 2008 Wenchuan Earthquake: A Critical-but-Constructive Review. Int. J. Geophys. 2012, 2012, 583097. [Google Scholar] [CrossRef]
  29. Mandelbrot, B.B.; Ness, J.W.V. Fractional Brownian motions, fractional noises and applications. J. Soc. Ind. Appl. Math 1968, 10, 422–437. [Google Scholar] [CrossRef]
  30. Morales, I.O.; Landa, O.; Fossion, R.; Frank, A. Scale invariance, self-similarity and critical behaviour in classical and quantum system. J. Phys. Conf. Ser. 2012, 380, 012020. [Google Scholar] [CrossRef]
  31. Musa, M.; Ibrahim, K. Existence of long memory in ozone time series. Sains Malays. 2012, 41, 1367–1376. [Google Scholar]
  32. Vadrevu, K.P. Fractal analysis revealed persistent correlations in long-term vegetation fire data in most South and Southeast Asian countries. Environ. Res. Commun. 2023, 5, 011001. [Google Scholar] [CrossRef]
  33. May, R.M. Simple mathematical models with very complicated dynamics. Nature 1976, 261, 459–467. [Google Scholar] [CrossRef]
  34. Sugihara, G.; May, R. Nonlinear forecasting as a way of distinguishing chaos from measurement error in time series. Nature 1990, 344, 734–741. [Google Scholar] [CrossRef]
  35. Liu, C.; Liang, J.; Li, Y.; Shi, K. Fractal analysis of impact of PM2.5 on surface O3 sensitivity regime based on field observations. Sci. Total Environ. 2023, 858, 160136. [Google Scholar] [CrossRef]
  36. Pastén, D.; Pavez-Orrego, C. Multifractal time evolution for intraplate earthquakes recorded in southern Norway during 1980–2021. Chaos Solitons Fractals 2023, 167, 113000. [Google Scholar] [CrossRef]
  37. Nikolopoulos, D.; Moustris, K.; Petraki, E.; Cantzos, D. Long-memory traces in PM _ 10 time series in Athens, Greece: Investigation through DFA and R/S analysis. Meteorol. Atmos. Phys. 2021, 133, 261–279. [Google Scholar] [CrossRef]
  38. Chelidze, T.; Matcharashvili, T.; Mepharidze, E.; Dovgal, N. Complexity in Geophysical Time Series of Strain/Fracture at Laboratory and Large Dam Scales: Review. Entropy 2023, 25, 467. [Google Scholar] [CrossRef]
  39. Nikolopoulos, D.; Moustris, K.; Petraki, E.; Koulougliotis, D.; Cantzos, D. Fractal and long-memory traces in PM10 time series in Athens, Greece. Environmnets 2019, 6, 29. [Google Scholar] [CrossRef]
  40. Hurst, H. Long term storage capacity of reservoirs. Trans. Am. Soc. Civ. Eng. 1951, 116, 770–808. [Google Scholar] [CrossRef]
  41. Hurst, H.; Black, R.; Simaiki, Y. Long-Term Storage: An Experimental Study; Constable: London, UK, 1965. [Google Scholar]
  42. Lopez, T.; Martınez-Gonzalez, C.; Manjarrez, J.; Plascencia, N.; Balankin, A. Fractal Analysis of EEG Signals in the Brain of Epileptic Rats, with and without Biocompatible Implanted Neuroreservoirs. AMM 2009, 15, 127–136. [Google Scholar] [CrossRef]
  43. Eftaxias, K.; Balasis, G.; Contoyiannis, Y.; Papadimitriou, C.; Kalimeri, M. Unfolding the procedure of characterizing recorded ultra low frequency, kHZ and MHz electromagnetic anomalies prior to the L’Aquila earthquake as pre-seismic ones—Part 2. Nat. Hazards Earth Syst. Sci. 2010, 10, 275–294. [Google Scholar] [CrossRef]
  44. Kilcik, A.; Anderson, C.; Rozelot, J.; Ye, H.; Sugihara, G.; Ozguc, A. Nonlinear Prediction of Solar Cycle 24. Astrophys. J. 2009, 693, 1173–1177. [Google Scholar] [CrossRef]
  45. Chattopadhyay, A.; Khondekar, M. An investigation of the relationship between the CME and the Geomagnetic Storm. Astron. Comput. 2023, 43, 100695. [Google Scholar] [CrossRef]
  46. Granero, M.S.; Segovia, J.T.; Perez, J.G. Some comments on Hurst exponent and the long memory processes on capital Markets. Phys. A Stat. Mech. Its Appl. 2008, 387, 5543–5551. [Google Scholar] [CrossRef]
  47. Musaev, A.; Makshanov, A.; Grigoriev, D. The Genesis of Uncertainty: Structural Analysis of Stochastic Chaos in Finance Markets. Complexity 2023, 2023, 1302220. [Google Scholar] [CrossRef]
  48. Pérez-Sienes, L.; Grande, M.; Losada, J.C.; Borondo, J. The Hurst Exponent as an Indicator to Anticipate Agricultural Commodity Prices. Entropy 2023, 25, 579. [Google Scholar] [CrossRef] [PubMed]
  49. Vogl, M. Hurst exponent dynamics of S&P 500 returns: Implications for market efficiency, long memory, multifractality and financial crises predictability by application of a nonlinear dynamics analysis framework. Chaos Solitons Fractals 2023, 166, 112884. [Google Scholar] [CrossRef]
  50. Dattatreya, G. Hurst Parameter Estimation from Noisy Observations of Data Traffic Traces. In Proceedings of the 4th WSEAS International Conference on Electronics, Control and Signal Processing, Miami, FL, USA, 17–19 November 2005; pp. 193–198. [Google Scholar]
  51. Wang, F.; Wang, H.; Zhou, X.; Fu, R. A Driving Fatigue Feature Detection Method Based on Multifractal Theory. IEEE Sens. J. 2022, 22, 19046–19059. [Google Scholar] [CrossRef]
  52. Zhou, H.; Chang, F. The long-memory temporal dependence of traffic crash fatality for different types of road users. Phys. A Stat. Mech. Appl. 2022, 607, 128210. [Google Scholar] [CrossRef]
  53. Li, X.; Polygiannakis, J.; Kapiris, P.; Peratzakis, A.; Eftaxias, K.; Yao, X. Fractal spectral analysis of pre-epileptic seizures in terms of criticality. J. Neural Eng. 2005, 2, 11–16. [Google Scholar] [CrossRef] [PubMed]
  54. Escobar-Ipuz, F.; Torres, A.; García-Jiménez, M.; Basar, C.; Cascón, J.; Mateo, J. Prediction of patients with idiopathic generalized epilepsy from healthy controls using machine learning from scalp EEG recordings. Brain Res. 2023, 1798, 148131. [Google Scholar] [CrossRef] [PubMed]
  55. Wijayanto, I.; Humairani, A.; Hadiyoso, S.; Rizal, A.; Prasanna, D.L.; Tripathi, S.L. Epileptic seizure detection on a compressed EEG signal using energy measurement. Biomed. Signal Process. Control 2023, 85, 104872. [Google Scholar] [CrossRef]
  56. Rehman, S.; Siddiqi, A. Wavelet based Hurst exponent and fractal dimensional analysis of Saudi climatic dynamics. Chaos Solitons Fractals 2009, 39, 1081–1090. [Google Scholar] [CrossRef]
  57. Petraki, E. Electromagnetic Radiation and Radon-222 Gas Emissions as Precursors of Seismic Activity. Ph.D. Thesis, Department of Electronic and Computer Engineering, Brunel University London, London, UK, 2016. [Google Scholar]
  58. Fujinawa, Y.; Takahashi, K. Electromagnetic radiations associated with major earthquakes. Phys. Earth Planet. Inter. 1998, 105, 249–259. [Google Scholar] [CrossRef]
  59. Hayakawa, M.; Ida, Y.; Gotoh, K. Multifractal analysis for the ULF geomagnetic data during the Guam earthquake. In Proceedings of the IEEE 6th International Symposium on Electromagnetic Compatibility and Electromagnetic Ecology, Saint Petersburg, Russia, 21–24 June 2005; pp. 21–24, 239–243. [Google Scholar]
  60. Hayakawa, M. VLF/LF radio sounding of ionospheric perturbations associated with earthquakes. Sensors 2007, 7, 1141–1158. [Google Scholar] [CrossRef]
  61. Nikolopoulos, D.; Alam, A.; Petraki, E.; Papoutsidakis, M.; Yannakopoulos, P.; Moustris, K.P. Stochastic and self-organisation patterns in a 17-year PM10 time series in Athens, Greece. Entropy 2021, 23, 307. [Google Scholar] [CrossRef] [PubMed]
  62. Skordas, E.S. On the increase of the “non-uniform” scaling of the magnetic field variations before the M(w)9.0 earthquake in Japan in 2011. Chaos 2014, 24, 023131. [Google Scholar] [CrossRef] [PubMed]
  63. Stanley, H.E. Powerlaws and universality. Nature 1995, 378, 597–600. [Google Scholar] [CrossRef]
  64. Sarlis, N.; Skordas, E.; Varotsos, P.; Nagao, T.; Kamogawa, M.; Tanaka, H.; Uyeda, S. Minimum of the order parameter fluctuations of seismicity before major earthquakes in Japan. Proc. Natl. Acad. Sci. USA 2013, 110, 13734–13738. [Google Scholar] [CrossRef]
  65. Becker, M.; Karpytchev, M.; Hu, A. Increased exposure of coastal cities to sea-level rise due to internal climate variability. Nat. Clim. Chang. 2023, 13, 367–374. [Google Scholar] [CrossRef]
  66. Ivanova, K.; Ausloos, M. Application of the detrended fluctuation analysis (DFA) method for describing cloud breaking. Phys. A Stat. Mech. Appl. 1999, 274, 349–354. [Google Scholar] [CrossRef]
  67. Koscielny-Bunde, E.; Bunde, A.; Havlin, S.; Roman, H.E.; Goldreich, Y.; Schellnhuberet, H. Indication of a Universal Persistence Law Governing Atmospheric Variability. Phys. Rev. Lett. 1998, 81, 729–732. [Google Scholar] [CrossRef]
  68. Rahmani, F.; Fattahi, M.H. Climate change-induced influences on the nonlinear dynamic patterns of precipitation and temperatures (case study: Central England). Theor. Appl. Climatol. 2023, 152, 1147–1158. [Google Scholar] [CrossRef]
  69. Linhares, R.R. Fractional poisson process: Long-range dependence in DNA sequences. Model Assist. Stat. Appl. 2023, 18, 33–43. [Google Scholar] [CrossRef]
  70. Peng, C.; Mietus, J.; Havlin, S.; Stanley, H.; Goldberger, A. Long-range anti-correlations and non-Gaussian behavior of the heartbeat. Phys. Rev. Lett. 1993, 70, 1343–1346. [Google Scholar] [CrossRef] [PubMed]
  71. Buldyrev, S.V.; Goldberger, A.L.; Havlin, S.; Mantegna, R.N.; Matsa, M.E.; Peng, C.K.; Simons, M.; Stanley, H.E. Long-range correlation properties of coding and noncoding DNA sequences: GenBank analysis. Phys. Rev. E Stat. Phys. Plasmas Fluids Relat. Interdiscip. Top. 1995, 51, 5084–5091. [Google Scholar] [CrossRef] [PubMed]
  72. Ivanov, P.C.; Rosenblum, M.G.; Peng, C.K.; Mietus, J.E.; Havlin, S.; Stanley, H.E.; Goldberger, A.L. Multifractality in human heartbeat dynamics. Nature 1999, 399, 461–465. [Google Scholar] [CrossRef]
  73. Mateo-March, M.; Moya-Ramón, M.; Javaloyes, A.; Sánchez-Muñoz, C.; Clemente-Suárez, V.J. Validity of detrended fluctuation analysis of heart rate variability to determine intensity thresholds in elite cyclists. Eur. J. Sport Sci. 2023, 23, 580–587. [Google Scholar] [CrossRef] [PubMed]
  74. Rogers, B.; Schaffarczyk, M.; Gronwald, T. Improved Estimation of Exercise Intensity Thresholds by Combining Dual Non-Invasive Biomarker Concepts: Correlation Properties of Heart Rate Variability and Respiratory Frequency. Sensors 2023, 23, 1973. [Google Scholar] [CrossRef] [PubMed]
  75. Eftaxias, K.; Balasis, G.; Contoyiannis, Y.; Papadimitriou, C.; Kalimeri, M.; Athanasopoulou, L.; Nikolopoulos, S.; Kopanas, J.; Antonopoulos, G.; Nomicos, C. Unfolding the procedure of characterizing recorded ultra low frequency, kHZ and MHz electromagnetic anomalies prior to the L’Aquila earthquake as pre-seismic ones-Part 1. Nat. Hazards Earth Syst. Sci. 2009, 9, 1953–1971. [Google Scholar] [CrossRef]
  76. Gotoh, K.; Hayakawa, M.; Smirnova, N.; Hattori, K. Fractal analysis of seismogenic ULF emissions. Phys. Chem. Earth 2004, 29, 419–424. [Google Scholar] [CrossRef]
  77. Hayakawa, M.; Ida, Y.; Gotoh, K. Fractal (mono- and multi-) analysis for the ULF data during the 1993 Guam earthquake for the study of prefracture criticality. Curr. Dev. Theory Appl. Wavelets 2008, 2, 159–174. [Google Scholar]
  78. Varotsos, P.; Sarlis, N.; Skordas, E. Natural Time Analysis: The New View of Time. Precursory Seismic Electric Signals, Earthquakes and Other Complex Time-Series; Springer: Berlin/Heidelberg, Germany, 2011. [Google Scholar]
  79. Varotsos, P.A.; Sarlis, N.V.; Skordas, E.S. Order Parameter and Entropy of Seismicity in Natural Time before Major Earthquakes: Recent Results. Geosciences 2022, 12, 225. [Google Scholar] [CrossRef]
  80. Hu, K.; Ivanov, P.; Chen, Z.; Carpena, P.; Stanley, H. Effect of trends on Detrended Fluctuation Analysis. Phys. Rev. E 2001, 64, 1–19. [Google Scholar] [CrossRef]
  81. Peng, C.; Buldyrev, S.; Simons, M.; Havlin, S.; Stanley, H.; Goldberger, A. On the mosaic organization of DNA sequences. Phys. Rev. E 1994, 49, 1685–1689. [Google Scholar] [CrossRef] [PubMed]
  82. Peng, C.; Havlin, S.; Stanley, H.; Goldberger, A. Quantification of scaling exponents and crossover phenomena in nonstationary heartbeat time series. Chaos 1995, 5, 82–87. [Google Scholar] [CrossRef]
  83. Peng, C.; Hausdor, J.; Havlin, S.; Mietus, J.; Stanley, H.; Goldberger, A. Multiple-time scales analysis of physiological time series under neural control. Phys. A Stat. Mech. Appl. 1998, 249, 491–500. [Google Scholar] [CrossRef] [PubMed]
  84. Nikolopoulos, D.; Yannakopoulos, P.H.; Petraki, E.; Cantzos, D.; Nomicos, C. Long-Memory and Fractal Traces in kHz-MHz Electromagnetic Time Series Prior to the ML = 6.1, 12/6/2007 Lesvos, Greece Earthquake: Investigation through DFA and Time-Evolving Spectral Fractals. J. Earth Sci. Clim. Chang. 2018, 9, 1–15. [Google Scholar]
  85. Nikolopoulos, D.; Petraki, E.; Nomicos, C.; Koulouras, G.; Kottou, S.; Yannakopoulos, P.H. Long-Memory Trends in Disturbances of Radon in Soil Prior ML=5.1 Earthquakes of 17 November 2014 Greece. J. Earth Sci. Clim. Chang. 2015, 6, 244. [Google Scholar]
  86. Katz, M. Fractals and the analysis of waveforms. Comput. Biol. Med. 1988, 18, 145–156. [Google Scholar] [CrossRef] [PubMed]
  87. Raghavendra, B.; Dutt, D.N. Computing Fractal Dimension of Signals using Multiresolution Box-counting Method. Int. J. Electr. Comput. Energetic Electron. Commun. Eng. 2010, 4, 183–198. [Google Scholar]
  88. Higuchi, T. Approach to an irregular time series on basis of the fractal theory. Phys. D Nonlinear Phenom. 1988, 31, 277–283. [Google Scholar] [CrossRef]
  89. de la Torre, F.C.; Ramirez-Rojas, A.; Pavia-Miller, C.; Angulo-Brown, F.; Yepez, E.; Peralta, J. A comparison between spectral and fractal methods in electrotelluric time series. Rev. Mex. Fis. 1999, 45, 298–302. [Google Scholar]
  90. de la Torre, F.C.; Gonzaalez-Trejo, J.; Real-Ramírez, C.; Hoyos-Reyes, L. Fractal dimension algorithms and their application to time series associated with natural phenomena. J. Phys. Conf. Ser. 2013, 475, 1–10. [Google Scholar] [CrossRef]
  91. Sevcik, C. On fractal dimension of waveforms. Chaos Solit. Fract. 2006, 27, 579–580. [Google Scholar] [CrossRef]
  92. Cantzos, D.; Nikolopoulos, D.; Petraki, E.; Yannakopoulos, P.H.; Nomicos, C. Earthquake precursory signatures in electromagnetic radiation measurements in terms of day-to-day fractal spectral exponent variation: Analysis of the eastern Aegean 13/04/2017–20/07/2017 seismic activity. J. Seismol. 2018, 22, 1499–1513. [Google Scholar] [CrossRef]
  93. Ida, Y.; Hayakawa, M. Fractal analysis for the ULF data during the 1993 Guam earthquake to study prefracture criticality. Nonlin. Process. Geophys. 2012, 13, 409–412. [Google Scholar] [CrossRef]
  94. Ida, Y.; Yang, D.; Li, Q.; Sun, H.; Hayakawa, M. Fractal analysis of ULF electromagnetic emissions in possible association with earthquakes in China. Nonlin. Process. Geophys. 2012, 19, 577–583. [Google Scholar] [CrossRef]
  95. Smirnova, N.; Hayakawa, M. Fractal characteristics of the ground-observed ULF emissions in relation to geomagnetic and seismic activities. J. Atmos. Sol. Ter. Phy. 2007, 69, 1833–1841. [Google Scholar] [CrossRef]
  96. Yonaiguchi, N.; Ida, Y.; Hayakawa, M.; Masuda, S. Fractal analysis for VHF electromagnetic noises and the identification of preseismic signature of an earthquake. J. Atmos. Sol. Ter. Phy. 2007, 69, 1825–1832. [Google Scholar] [CrossRef]
  97. Eftaxias, K. Footprints of non-extensive Tsallis statistics, self-affinity and universality in the preparation of the L’Aquila earthquake hidden in a pre-seismic EM emission. Phys. A Stat. Mech. Appl. 2010, 389, 133–140. [Google Scholar] [CrossRef]
  98. Kapiris, P.; Peratzakis, J.P.A.; Nomikos, K.; Eftaxias, K. VHF-electromagnetic evidence of the underlying pre-seismic critical stage. Earth Planets Space 2002, 54, 1237–1246. [Google Scholar] [CrossRef]
  99. Nikolopoulos, D.; Petraki, E.; Marousaki, A.; Potirakis, S.; Koulouras, G.; Nomicos, C.; Panagiotaras, D.; Stonhamb, J.; Louizi, A. Environmental monitoring of radon in soil during a very seismically active period occurred in South West Greece. J. Environ. Monit. 2012, 14, 564–578. [Google Scholar] [CrossRef]
  100. Telesca, L.; Lasaponara, R. Vegetational patterns in burned and unburned areas investigated by using the detrended fluctuation analysis. Phys. A Stat. Mech. Appl. 2006, 368, 531–535. [Google Scholar] [CrossRef]
  101. Eftaxias, K.; Contoyiannis, Y.; Balasis, G.; Karamanos, K.; Kopanas, J.; Antonopoulos, G.; Koulouras, G.; Nomicos, C. Evidence of fractional-Brownian-motion-type asperity model for earthquake generation in candidate pre-seismic electromagnetic emissions. Nat. Hazard Earth Syst. 2008, 8, 657–669. [Google Scholar] [CrossRef]
  102. Varotsos, P.; Alexopoulos, K. Physical properties of the variations of the electric field of the earth preceding earthquakes, I. Tectonophysics 1984, 110, 73–98. [Google Scholar] [CrossRef]
  103. Varotsos, P.; Alexopoulos, K. Physical properties of the variations of the electric field of the earth preceding earthquakes, II. Tectonophysics 1984, 110, 99–125. [Google Scholar] [CrossRef]
  104. Varotsos, P.; Sarlis, N.; Lazaridou, M.B.N. Statistical evaluation of earthquake prediction results. Comments Success Rate Alarm. Rate Acta Geophys. Pol. 1996, 44, 329–347. [Google Scholar]
  105. Varotsos, P.; Sarlis, N.; Skordas, E. Magnetic field variations associated with SES. The instrumentation used for investigating their detectability. Proc. Jpn. Acad. Ser. B 2001, 77, 87–92. [Google Scholar] [CrossRef]
  106. Petraki, E.; Nikolopoulos, D.; Chaldeos, Y.; Koulouras, G.; Nomicos, C.; Yannakopoulos, P.H.; Kottou, S.; Stonham, J. Fractal evolution of MHz electromagnetic signals prior to earthquakes: Results collected in Greece during 2009. Geomat. Nat. Hazards Risk 2016, 7, 550–564. [Google Scholar] [CrossRef]
  107. Pinault, J.; Baubron, J. Signal processing of soil gas radon, atmospheric pressure, moisture, and soil temperature data: A new approach for radon concentration modeling. J. Geoph. Res. Sol. EA 1996, 101, 3157–3171. [Google Scholar] [CrossRef]
  108. Eftaxias, K.; Panin, V.; Deryugin, Y. Evolution-EM signals before earthquakes in terms of mesomechanics and complexity. Phys. Chem. Earth 2007, 29, 445–451. [Google Scholar] [CrossRef]
  109. Eftaxias, K.; Sgrigna, V.; Chelidze, T. Mechanical and electromagnetic phenomena accompanying preseismic deformation: From laboratory to geophysical scale. Tectonophysics 2007, 341, 1–5. [Google Scholar] [CrossRef]
  110. Gotoh, K.; Hayakawa, M.; Smirnova, N. Fractal analysis of the ULF geomagnetic data obtained at Izu Peninsula, Japan in relation to the nearby earthquake swarm of June–August 2000. Nat. Haz. Earth Syst. 2003, 3, 229–234. [Google Scholar] [CrossRef]
  111. Hayakawa, M.; Kawate, R.; Molchanov, O.; Yumoto, K. Results of ultra-low-frequency magnetic field measurements during the Guam earthquake of 8 August 1993. Geophys. Res. Lett. 1996, 23, 241–244. [Google Scholar] [CrossRef]
  112. Hayakawa, M.; Itoh, T.; Hattori, K.; Yumoto, K. ULF electromagnetic precursors for an earthquake at Biak, Indonesia on 17 February 1996. Geophys. Res. Lett. 2000, 27, 1531–1534. [Google Scholar] [CrossRef]
  113. Kapiris, P.; Eftaxias, K.; Nomikos, K.; Polygiannakis, J.; Dologlou, E.; Balasis, G.; Bogris, N.; Peratzakis, A.; Hadjicontis, V. Evolving towards a critical point: A possible electromagnetic way in which the critical regime is reached as the rupture approaches. Nonlinear Proc. Geoph. 2003, 10, 1–14. [Google Scholar] [CrossRef]
  114. Smirnova, N.; Hayakawa, M.; Gotoh, K. Precursory behavior of fractal characteristics of the ULF electromagnetic fields in seismic active zones before strong earthquakes. Phys. Chem. Earth 2004, 29, 445–451. [Google Scholar] [CrossRef]
  115. Smirnova, N.A.; Kiyashchenko, D.A.; Troyan, V.N.; Hayakawa, M. Multifractal Approach to Study the Earthquake Precursory Signatures Using the Ground-Based Observations. Rev. Appl. Phys. 2013, 2, 3. [Google Scholar]
  116. USGS. Available online: https://earthquake.usgs.gov/earthquakes/map/?extent=22.41103,91.51611&extent=35.56798,113.57666&range=search&baseLayer=terrain&timeZone=utc&search=%7B%22name%22:%22Search%20Results%22,%22params%22:%7B%22starttime%22:%222008-01-01%2000:00:00%22,%22endtime%22:%222009-01-01%2000:00:00%22,%22maxlatitude%22:36.844,%22minlatitude%22:22.999,%22maxlongitude%22:116.323,%22minlongitude%22:96.504,%22minmagnitude%22:5.5,%22orderby%22:%22time%22%7D%7D (accessed on 26 May 2023).
Figure 1. Location of the epicentre of the Wenchuan earthquake, showing the radon in groundwater monitoring stations together with elevation data and significant geological faults. The location of the presented area within China is shown in the inserted sub-figure at the bottom right. The codes of stations are given in Table 1. Coloured gradient is the altitude elevation in m.
Figure 1. Location of the epicentre of the Wenchuan earthquake, showing the radon in groundwater monitoring stations together with elevation data and significant geological faults. The location of the presented area within China is shown in the inserted sub-figure at the bottom right. The codes of stations are given in Table 1. Coloured gradient is the altitude elevation in m.
Geosciences 13 00268 g001
Figure 2. Results of DFA. KDS (ID = 3). Window of 64 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of F ( n ) versus n in every 64-sample window; (top) the scaling exponent, α (DFA slope). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Figure 2. Results of DFA. KDS (ID = 3). Window of 64 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of F ( n ) versus n in every 64-sample window; (top) the scaling exponent, α (DFA slope). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g002
Figure 3. Results of DFA. GS (ID = 82). Window of 64 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of F ( n ) versus n in every 64-sample window; (top) the scaling exponent, α (DFA slope). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Figure 3. Results of DFA. GS (ID = 82). Window of 64 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of F ( n ) versus n in every 64-sample window; (top) the scaling exponent, α (DFA slope). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g003
Figure 4. Results of DFA. MSS (ID = 83). Window of 64 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of F ( n ) versus n in every 64-sample window; (top) the scaling exponent, α (DFA slope). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Figure 4. Results of DFA. MSS (ID = 83). Window of 64 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of F ( n ) versus n in every 64-sample window; (top) the scaling exponent, α (DFA slope). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g004
Figure 5. Results of DFA. PZHS (ID = 143). Window of 64 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of F ( n ) versus n in every 64-sample window; (top) the scaling exponent, α (DFA slope). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Figure 5. Results of DFA. PZHS (ID = 143). Window of 64 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of F ( n ) versus n in every 64-sample window; (top) the scaling exponent, α (DFA slope). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g005
Figure 6. Results of DFA. SPS (ID = 149). Window of 64 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of F ( n ) versus n in every 64-sample window; (top) the scaling exponent, α (DFA slope). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Figure 6. Results of DFA. SPS (ID = 149). Window of 64 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of F ( n ) versus n in every 64-sample window; (top) the scaling exponent, α (DFA slope). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g006
Figure 7. Results from fractal dimension analysis. KDS (ID = 3). Window of 64 samples, 16 sub-categories of Higuchi’s method, and step of 1 sample. From bottom to top: (a) the radon in groundwater time series and the fractal dimensions according to the algorithms of (b) Katz (KFD), (c) Higuchi (HFD), and (d) Sevcik (SFD). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Figure 7. Results from fractal dimension analysis. KDS (ID = 3). Window of 64 samples, 16 sub-categories of Higuchi’s method, and step of 1 sample. From bottom to top: (a) the radon in groundwater time series and the fractal dimensions according to the algorithms of (b) Katz (KFD), (c) Higuchi (HFD), and (d) Sevcik (SFD). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g007
Figure 8. Results from fractal dimension analysis. GS (ID = 82). Window of 64 samples, 16 sub-categories of Higuchi’s method, and step of 1 sample. From bottom to top: (a) the radon in groundwater time series and the fractal dimensions according to the algorithms of (b) Katz (KFD), (c) Higuchi (HFD), and (d) Sevcik (SFD). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Figure 8. Results from fractal dimension analysis. GS (ID = 82). Window of 64 samples, 16 sub-categories of Higuchi’s method, and step of 1 sample. From bottom to top: (a) the radon in groundwater time series and the fractal dimensions according to the algorithms of (b) Katz (KFD), (c) Higuchi (HFD), and (d) Sevcik (SFD). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g008
Figure 9. Results from fractal dimension analysis. MSS (ID = 83). Window of 64 samples, 16 sub-categories of Higuchi’s method, and step of 1 sample. From bottom to top: (a) the radon in groundwater time series and the fractal dimensions according to the algorithms of (b) Katz (KFD), (c) Higuchi (HFD), and (d) Sevcik (SFD). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Figure 9. Results from fractal dimension analysis. MSS (ID = 83). Window of 64 samples, 16 sub-categories of Higuchi’s method, and step of 1 sample. From bottom to top: (a) the radon in groundwater time series and the fractal dimensions according to the algorithms of (b) Katz (KFD), (c) Higuchi (HFD), and (d) Sevcik (SFD). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g009
Figure 10. Results from fractal dimension analysis. PZHS (ID = 143). Window of 64 samples, 16 sub-categories of Higuchi’s method, and step of 1 sample. From bottom to top: (a) the radon in groundwater time series and the fractal dimensions according to the algorithms of (b) Katz (KFD), (c) Higuchi (HFD), and (d) Sevcik (SFD). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Figure 10. Results from fractal dimension analysis. PZHS (ID = 143). Window of 64 samples, 16 sub-categories of Higuchi’s method, and step of 1 sample. From bottom to top: (a) the radon in groundwater time series and the fractal dimensions according to the algorithms of (b) Katz (KFD), (c) Higuchi (HFD), and (d) Sevcik (SFD). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g010
Figure 11. Results from fractal dimension analysis. SPS (ID = 149). Window of 64 samples, 16 sub-categories of Higuchi’s method, and step of 1 sample. From bottom to top: (a) the radon in groundwater time series and the fractal dimensions according to the algorithms of (b) Katz (KFD), (c) Higuchi (HFD), and (d) Sevcik (SFD). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Figure 11. Results from fractal dimension analysis. SPS (ID = 149). Window of 64 samples, 16 sub-categories of Higuchi’s method, and step of 1 sample. From bottom to top: (a) the radon in groundwater time series and the fractal dimensions according to the algorithms of (b) Katz (KFD), (c) Higuchi (HFD), and (d) Sevcik (SFD). The horizontal axis is from the beginning (1 January 2008) to the end (31 December 2008) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g011
Figure 12. Results from fractal analysis: KDS (ID = 3). Window of 128 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of Equation (18); (top) time evolution of power-law β exponent. Blue areas represent the successive ( r 2 0.95 ) fractal windows. Red areas are non-successive windows. The horizontal axis is from the beginning (1 January 2018) to the end (31 December 2018) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Figure 12. Results from fractal analysis: KDS (ID = 3). Window of 128 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of Equation (18); (top) time evolution of power-law β exponent. Blue areas represent the successive ( r 2 0.95 ) fractal windows. Red areas are non-successive windows. The horizontal axis is from the beginning (1 January 2018) to the end (31 December 2018) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g012
Figure 13. Results from fractal analysis: GS (ID = 82). Window of 128 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of Equation (18); (top) time evolution of power-law β exponent. Blue areas represent the successive ( r 2 0.95 ) fractal windows. Red areas are non-successive windows. The horizontal axis is from the beginning (1 January 2018) to the end (31 December 2018) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Figure 13. Results from fractal analysis: GS (ID = 82). Window of 128 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of Equation (18); (top) time evolution of power-law β exponent. Blue areas represent the successive ( r 2 0.95 ) fractal windows. Red areas are non-successive windows. The horizontal axis is from the beginning (1 January 2018) to the end (31 December 2018) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g013
Figure 14. Results from fractal analysis: MSS (ID = 83). Window of 128 samples and step of 1 sample. From bottom to top: (bottom) Radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of Equation (18); (top) Time evolution of power-law β exponent. Blue areas represent the successive ( r 2 0.95 ) fractal windows. Red areas are non-successive windows. The horizontal axis is from the beginning (1 January 2018) to the end (31 December 2018) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Figure 14. Results from fractal analysis: MSS (ID = 83). Window of 128 samples and step of 1 sample. From bottom to top: (bottom) Radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of Equation (18); (top) Time evolution of power-law β exponent. Blue areas represent the successive ( r 2 0.95 ) fractal windows. Red areas are non-successive windows. The horizontal axis is from the beginning (1 January 2018) to the end (31 December 2018) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g014
Figure 15. Results from fractal analysis: PZHS (ID = 143). Window of 128 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of Equation (18); (top) time evolution of power-law β exponent. Blue areas represent the successive ( r 2 0.95 ) fractal windows. Red areas are non-successive windows. The horizontal axis is from the beginning (1 January 2018) to the end (31 December 2018) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Figure 15. Results from fractal analysis: PZHS (ID = 143). Window of 128 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of Equation (18); (top) time evolution of power-law β exponent. Blue areas represent the successive ( r 2 0.95 ) fractal windows. Red areas are non-successive windows. The horizontal axis is from the beginning (1 January 2018) to the end (31 December 2018) of measurements. The measurement sampling rate is 1 day−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g015
Figure 16. Results from fractal analysis: SPS (ID = 149). Window of 128 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of Equation (18); (top) time evolution of power-law β exponent. Blue areas represent the successive ( r 2 0.95 ) fractal windows. Red areas are non-successive windows. The horizontal axis is from the beginning (1 January 2018) to the end (31 December 2018) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Figure 16. Results from fractal analysis: SPS (ID = 149). Window of 128 samples and step of 1 sample. From bottom to top: (bottom) radon in groundwater time series; (middle) Spearman’s correlation coefficient of the goodness of the linear fit of Equation (18); (top) time evolution of power-law β exponent. Blue areas represent the successive ( r 2 0.95 ) fractal windows. Red areas are non-successive windows. The horizontal axis is from the beginning (1 January 2018) to the end (31 December 2018) of measurements. The measurement sampling rate is 1 h−1. For the window segmentation, please refer to the text.
Geosciences 13 00268 g016
Figure 17. Location of the earthquakes of 2008 with M w 5.5 over an area greater than the one inserted in Figure 1. This figure was created with Google Earth using a kml file from USGS. Greater circles are earthquakes with greater magnitude M w . The biggest circle near Sichuan is the great Wenchuan earthquake. The 19 presented earthquakes are shown in Table 2.
Figure 17. Location of the earthquakes of 2008 with M w 5.5 over an area greater than the one inserted in Figure 1. This figure was created with Google Earth using a kml file from USGS. Greater circles are earthquakes with greater magnitude M w . The biggest circle near Sichuan is the great Wenchuan earthquake. The 19 presented earthquakes are shown in Table 2.
Geosciences 13 00268 g017
Figure 18. Overview of the full computational meta-analysis results by all thirteen selected combinations of fractal methods per five, four, three, and two methods. Data from (a) KDS, (b) GS, and (c) MSS. Symbols: “+” (black): DFA versus all methods (5 techniques—DFA, Power law analysis, Katz’s, Higuchi’s, and Sevcik’s methods); “⊙” (red): DFA versus all fractal dimension techniques (4 techniques); “∗” (green): power-law analysis versus all fractal dimension techniques (4 techniques); “.” (blue): DFA versus Higuchi’s and Katz’s methods (3 techniques); “⊡” (yellow): DFA versus Higuchi’s and Sevcik’s methods (3 techniques); “” (magenta): DFA versus Katz’s and Sevcik’s methods (3 techniques); “▽” (cyan): power-law analysis versus Higuchi’s and Katz’s methods (3 techniques); “▽” (black): power-law analysis versus Higuchi’s and Sevcik’s methods (3 techniques); “▷” (red): power-law analysis versus Katz’s and Sevcik’s methods (3 techniques); “◁” (green): Higuchi’s versus Katz’s and Sevcik’s methods (3 techniques); “⬠” (yellow): Sevcik’s versus Katz’s and Higuchi’s methods (3 techniques); “⬠” (blue): Katz’s versus Higuchi’s and Sevcik’s methods (3 techniques); “+ ’ (magenta): DFA versus power-law analysis (2 techniques). Horizontal axis is in actual dates.
Figure 18. Overview of the full computational meta-analysis results by all thirteen selected combinations of fractal methods per five, four, three, and two methods. Data from (a) KDS, (b) GS, and (c) MSS. Symbols: “+” (black): DFA versus all methods (5 techniques—DFA, Power law analysis, Katz’s, Higuchi’s, and Sevcik’s methods); “⊙” (red): DFA versus all fractal dimension techniques (4 techniques); “∗” (green): power-law analysis versus all fractal dimension techniques (4 techniques); “.” (blue): DFA versus Higuchi’s and Katz’s methods (3 techniques); “⊡” (yellow): DFA versus Higuchi’s and Sevcik’s methods (3 techniques); “” (magenta): DFA versus Katz’s and Sevcik’s methods (3 techniques); “▽” (cyan): power-law analysis versus Higuchi’s and Katz’s methods (3 techniques); “▽” (black): power-law analysis versus Higuchi’s and Sevcik’s methods (3 techniques); “▷” (red): power-law analysis versus Katz’s and Sevcik’s methods (3 techniques); “◁” (green): Higuchi’s versus Katz’s and Sevcik’s methods (3 techniques); “⬠” (yellow): Sevcik’s versus Katz’s and Higuchi’s methods (3 techniques); “⬠” (blue): Katz’s versus Higuchi’s and Sevcik’s methods (3 techniques); “+ ’ (magenta): DFA versus power-law analysis (2 techniques). Horizontal axis is in actual dates.
Geosciences 13 00268 g018
Table 1. Radon in groundwater station data. The numbers and coding of the stations are those globally adopted in China, and the names refer to the actual locations. Distance is the epicentral distance.
Table 1. Radon in groundwater station data. The numbers and coding of the stations are those globally adopted in China, and the names refer to the actual locations. Distance is the epicentral distance.
Station CodeStation NameLatitudeLongitudeDistance (km)
KDSKangding station30.12102.17152.2
GSGanzi station31.61100.01325.5
MSSMingshan station30.1103.1105.6
PZHSPanzhihua station26.51101.74526.0
SPSSonpan station32.65103.6182.5
Table 2. Earthquakes of 2008 with M w 5.5 in China for the area presented in Figure 17. The last event (i/i 19) is the great Wenchuan earthquake.
Table 2. Earthquakes of 2008 with M w 5.5 in China for the area presented in Figure 17. The last event (i/i 19) is the great Wenchuan earthquake.
i/iYearMonthDayHourMinuteSecond M w LatitudeLongitudeDepth (km)
12008831831105.626.232101.9710
22008830830536.026.241101.88911
320088211224306.025.03997.69710
4200885949176.032.756105.4946
5200881832435.732.033104.72211
6200872493095.732.747105.54210
720087231954445.532.752105.4984
82008527837515.732.71105.5410
92008525821496.132.56105.42318
102008517825485.832.24104.9829
112008516525475.631.355103.3513
1220085137785.830.89103.1949
132008512208505.631.413103.88921.7
142008512111126.131.214103.61810
152008512942245.531.527104.09210
162008512643145.831.211103.71510
17200851264285.731.342104.68210
182008512661565.731.586104.03210
19200851262817.931.002103.32219
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Alam, A.; Nikolopoulos, D.; Wang, N. Fractal Patterns in Groundwater Radon Disturbances Prior to the Great 7.9 Mw Wenchuan Earthquake, China. Geosciences 2023, 13, 268. https://doi.org/10.3390/geosciences13090268

AMA Style

Alam A, Nikolopoulos D, Wang N. Fractal Patterns in Groundwater Radon Disturbances Prior to the Great 7.9 Mw Wenchuan Earthquake, China. Geosciences. 2023; 13(9):268. https://doi.org/10.3390/geosciences13090268

Chicago/Turabian Style

Alam, Aftab, Dimitrios Nikolopoulos, and Nanping Wang. 2023. "Fractal Patterns in Groundwater Radon Disturbances Prior to the Great 7.9 Mw Wenchuan Earthquake, China" Geosciences 13, no. 9: 268. https://doi.org/10.3390/geosciences13090268

APA Style

Alam, A., Nikolopoulos, D., & Wang, N. (2023). Fractal Patterns in Groundwater Radon Disturbances Prior to the Great 7.9 Mw Wenchuan Earthquake, China. Geosciences, 13(9), 268. https://doi.org/10.3390/geosciences13090268

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop