Trend Prediction of Seismicity in Lower Reaches of Yangtze River – Yellow Sea (Y-Y) Seismic Belt in 2011-2020

In this study, seismic activity trend of lower reaches of Yangtze River – Yellow Sea (Y-Y) seismic belt during 2011-2020 is predicted based on statistic analysis of historical earthquake data and seismic activity analysis of recent mid-strong earthquakes occurred in this area. Three potential seismic risk zones are indicated where the destructive earthquakes might occur which will cause damages in that area. The three zones and the estimated magnitudes of potential earthquakes are: 1) southern sags of Southern Yellow sea and coastal regions, Ms ≥ 6; 2) Yangtze estuary and its eastern sea area, Ms ≥ 5; 3) areas along Changshu city to Taicang city, Ms ≥ 5. Seismic parametric study of that area is performed to support our conclusions.


Introduction
Major tasks of the analysis of seismicity trend include determine the current phase of the seismicity in studied seismic zones and belts, roughly estimate future seismicity level in those areas, and predict the locations in which strong seismic activities may occur. In general earthquake prediction, the objective of seismicity analysis trend is to correctly describe the time, location, and intensity of the future earthquakes. The prediction opinion will be made based on the possibility that most likely to occur.
This study investigated the comparison of historical seismic activity staging, multiple mathematical statistics models, and recent seismicity in Shanghai and its vicinity, etc. The final predictions are made based on these investigation as well as the analysis of medium's Q value and other seismic source dynamics parameters The seismic zone from lower reaches of Yangtze River to southern Yellow Sea is chosen as the study area and our analysis results will be compared to the previous data.

Analysis Method
Several methods are being used for seismicity trend analysis. One method is to compare characteristics of seismic active period, include the lasting time duration, seismicity level (number of strong earthquakes, total released energy, etc.), and evolution characteristics of the seismicity process, etc; and predict the characteristics of future seismicity in the current seismic active period based on the comparison results. However, due to the limit of historic seismic data and the complexity of seismic process, features of different active period may show some repeatability but usually render distinct differences.
Another method is to apply seismic related and triggering factors in different seismic zones and belts. For example, the changing of the Earth's rotation speed and sunspot activity level may induce the seismic activity. Fig. 1 clearly shows that the seismicity in China synchronizes with the global seismicity. Therefore the seismicity in mainland of China can be inferred from the global seismicity. Fig. 2 explains that when the Earth rotates faster, the activity of shallow earthquakes in North China-Northeast Plain usually decreases, while the seismicity in the regions of Liupan mountain-Qilian mountain-Altyn mountain increases, and vice versa. Thus, the seismicity level in those areas can be predicted from the variation of the Earth's rotation speed.
The third method is based on characteristic analysis of the curves of released strain energy in seismic areas. The released strain energy curves can be classified to three types. (1) Linear type with stable slope, the corresponding seismic areas have comparatively high seismicity level, where strong earthquakes occur frequently with regular time interval. The (Y-Y) Seismic Belt in 2011-2020 maximum earthquake magnitude can be extrapolated from the curve (as shown in Fig. 3). (2) The curve with increasing slope, which indicates a high level of accumulated strain energy in corresponding seismic areas thereby indicating strong seismic activities in the future (Fig. 4). (3) The curve with decreasing slope, from which it can be deduced that the future seismicity level will not be too high in the corresponding seismic areas (Fig. 5).   Another type of method is mathematical statistics, which selects a character of the seismic activity (such as the magnitude, energy, time interval, etc.) and assumes that its data has some type of distribution (normal distribution, Gaussian distribution, Weibull distribution, etc.). The selected data is arranged into a series according to a certain order (magnitude order or chronological order, etc.), and fitted by a curve; the seismic trend then can be determined from interpolation or extrapolation of that curve. There are a number of methods that belong to the group of mathematical statistics method, such as maximum prediction, regression analysis, periodogram analysis, statistical prediction of stationary independent processes, wavelet analysis, entropy spectrum analysis, and so on. In applying those methods, there are two things need to be noticed: (1) the prediction results obtained from different methods may not be identical and even contradict each other, therefore we have to use different methods to obtain the results and find the average value as the final decision; (2) interpolation is superior to extrapolation because the latter can only predict one point at each time.
Among the presented methods, method 1 and 3 are main methods which are indispensible for seismicity trend analysis while method 2 and 4 are auxiliary methods. However, in current seismic safety evaluation, the seismicity trend is predicted based on mathematical statistics (method 4). Even worse, in using mathematical statistics to predict the future seismicity trend, only periodogram analysis is being used where five extra points are extrapolated from current periodogram to predict the seismicity trend in next 50 years. The prediction made using this method obviously underestimated the potential seismic risk, sometimes is inconsistent with the real situations of historic and contemporary seismic activity, and may be in disagreement with the zonation of the potential seismic source area. Fortunately, the seismicity trend analysis results did not affect the results of seismic hazard assessment because the zonation of the potential seismic source area is appropriate and various factors that affect the seismic hazard assessment have been fully considered. In deciding the seismicity parameters, the α value of the central area was increased to 11.35, thus, the ground motion parameters and conclusions made in current seismic evaluation reports of Shanghai are reliable and safe enough.

Working Procedure
(1) Select appropriate area for study. In this step, we need to select an entire seismic active area, neither choose part of that area nor mix one seismic active area with another. If the selected area is not appropriate for seismic study, it will be difficult to obtain regular historic statistical law from the study, therefore fail to provide meaningful historic data for seismicity trend analysis.
(2) Investigate the seismicity in different seismic periods by making diagrams of M-t, E t -å , and N-t. Plot the distribution of epicenter in different seismic periods to study the spatial evolving process of the seismic activity. Important factors that affect the spatial evolving process include: (a) if those seismic periods appeared regularly, the duration of each seismic period, their seismicity level, and similarities and differences among the seismic activities in different periods; (b) characteristics of the spatial distribution of different seismic periods; (c) correlations (especially the sequential relationship) among seismic activities in different seismic zones and belts, and epicenter migration; (d) effects of some environmental factors (rotation of the Earth, sunspot, tides, etc.).
(3) Determine the current stage of the seismic activities and predict the future seismicity level. For example, the seismicity in a seismic zone will be in active period in the near future if one of following conditions is satisfied. (a) The seismicity is in an active period and its lasting time is shorter than the shortest active period in history; (b) it is in a quiet period but its lasting time is close to the longest quiet period in history; or (c) the strain release curve accelerates.
(4) Determine and predict the intensity of seismic activity using extrapolation of the strain release curve or the data obtained from mathematical statistics.
(5) Roughly estimate the risk zones of strong earthquakes, major seismic active areas in current stage, epicenter migration, immunity zones of strong earthquakes, and intervals of time in which strong earthquakes occur.

Data
The study area of this research, the seismic zone from lower reaches of YangtzeRiver to southern Yellow Sea, includes two areas: Eastern China (29º-37º N, 100º-124º E) and Shanghai and its vicinity (29º-34º N, 118º-124º E). All the historical seismic data of these two areas used for this study come from local historical records (for earthquakes occurred before 20 th century) and seismic data recorded by local observatories (for earthquakes occurred since 20 th century, especially since last 70's) [1][2][3]. Lee et al. [4] studied the conserving probability of seismic data in Eastern China based on completeness of the historical records and obtained the conserving probabilities of Eastern Chinese seismic data in different centuries (Table 1). From Table 1, it can be found that a complete historical catalog (with probability close to 1) of Eastern Chinese earthquakes is only available after 15 th or 16 th century (after 1400 or 1500), which can be used for statistical analysis.
Based on all available historical seismic data, it can be concluded that: (1) historical records for earthquakes occurred in that area before 1491 are seriously missed. (2) From 1491 till 1900, the records for land and coastal earthquakes with Ms ≥ 5.0 are fairly complete, but most sea earthquakes were not recorded except four sea earthquakes with Ms ≥ 6.7. (3) Since 1905, with application of modern monitoring instruments, a complete seismic data of all the earthquakes with Ms ≥ 5.0 in this area becomes available. Hence, in this research we investigate the seismic data of earthquakes with Ms ≥ 5.0 which occurred after 1491 and especially focus on those occurred since 1905.

Comparison of Active Periods
Following figures and tables compare the characteristics of earthquakes occurred in the Y-Y seismic zone during different episodes since 1491.

Statistical Forecast of Seismicity in Subareas
Next, we divided the Y-Y belt into three subzones (north, central, and south) and the intensity of seismicity in the decades to come in those subzones were estimated using statistical forecast.  Table 4. An alternative method to determine the earthquake recurrence period is the maximum entropy spectral analysis, whose maximum entropy power spectral density S(f) is: In above equation, A m,i is the filtration coefficient of predicting errors, P m is the variance of parameter estimation based on the autoregressive model, m is the order of the power spectrum, and f i is the autocorrelation function of the stationary process.The maximum entropy spectral analysis results are displayed in Table 5. The seismicity level in those areas can also be predicted through periodogram analysis of Fourier series. In such analysis, we use every 10 years as an interval, and the logarithm of the maximum released energy during each interval, logE, is used as a statistical data. Table 6 lists the prediction results. Next, five statistical methods -fuzzy modular neural network (method 1), gray modeling method (method 2), extreme value theory (method 3), prediction method based on probability of distribution of time interval of strong earthquakes (method 4), and Weibull distribution of time interval of strong earthquakes (method 5) were applied to separately estimate the maximum magnitude of the earthquakes that will occur in the Y-Y seismic belt during next 50 and 100 years. Table 7 displays the estimated maximum earthquake magnitude.

Studied Area and Data
In order to reveal the relationship among earthquake magnitude, Q value, and stress drop of small earthquakes, the study area (29º~34ºN, 119º~124ºE) is divided into five small regions based on distributions of potential hypocentral regions and earthquakes so that each region possesses a certain amount of records of earthquake data (Fig. 13). Average Q value and stress drop Δσ are calculated for each region. Specifically, since the Ms6.1 occurred in Nov-9-1996 at (32ºN, 123ºE) had a evident influence on Shanghai are, the current Q values and Δσ are calculated as well as the then values before and right after that earthquake. Table 7. Statistical estimation of maximum earthquake magnitude during the next 50 and 100 years in Y-Y seismic belt and its three subzones Seismic data are selected from the records of Shanghai earthquake networks. Those clearly recorded and unclipped seismograms are selected. Magnitudes of selected earthquakes range from M L 2.0 ~ 4.0 and epicenter distances vary within 100 ~ 200km, aftershocks are excluded. Additionally, seven digital seismic records are selected, which were recorded by Shanghai digital earthquake network, which was put into operation in 2001. Thus, overall 92 earthquakes are used for our study.

Methods and Data Processing
Solve medium's average Q values Assume n earthquakes are recorded by an observatory within the same area and spectrum of the i th earthquake is: (3) where A i0 (ω) is the seismic spectrum emitted from the hypocenter, G(R i ) is geometric diffusion, R i is hypocentral distance, I(ω) is frequency characteristic of instruments, exp[-(ωR i )/(2v p Q)] is medium absorption, and v p is the velocity of P-wave which takes the value 5.7km/s. Substitute two fixed frequencies, ω 1 and ω 2 , into above equation separately and make a division, we can have Ln|A i (ω 1 )/A i (ω 2 )| = Ln|A i0 (ω 1 )/A i0 (ω 2 )| -(ω 1ω 2 )R i /(2v p Q) (4) where (ω 1ω 2 )R i /(2v p Q) is a constant. Meanwhile, since our selected range of earthquake magnitude is comparatively small and those earthquakes occurred in the same region, therefore the Ln|A i (ω 1 )/A i (ω 2 )| of different earthquakes can also be regarded as a constant. Those assumptions transform Eqn. (4) to Y i = a -b×R i (5) In above equation, Y i can be calculated, R i is known, and coefficients a and b can be solved using least squares method. The average Q value of this region then can be solved from b as: Q = (ω 1ω 2 )/(2v p b) (6) Solve stress drop Δσ By using dislocation model of circular shear to simulate the medium and small strike-slip earthquakes and taking samples from the seismogram and then performing Fourier analysis, the source spectra can be obtained after taking off instrumental and medium effects. From that spectra the spectral amplitude in lower band u(ω) ω→0 and corner frequency f cα can be obtained, therefore the seismic moment m 0 can be determined from where ρ is the medium's density (we use 2.7 × 10 3 kg/m 3 for the crust), r is the hypocentral distance, the radiation pattern R α = sin(2θ)cosΦ (if θ and Φ are unknowns, we can take the average radiation pattern calculated over the focal sphere, 4/15 for the P wave [8], and v p means P-wave velocity (we assume v p = 5.7km/s within the crust [9]). The radius of rupture circle a can be calculated from the corner frequency f cα as: In above equation, the rupture propagation velocity v f = 0.775v s , where v s is S-wave velocity (we assume v s = 3.38km/s within the crust [10]). Finally, the stress drop Δσ can be calculated using Eqn. (9) based on m 0 and a Δσ = (7m 0 ) / (16a 3 ) (9) Data processing For the analogous records, they are scanned and digitized through the seismogram digitization and database management system (SDDMS) [11]. Through that system the sampling rate can be selected. Considering quality of the analogous records and range of interested corner frequency f cα , the sampling rate is set as 200 dpi. Sampling step size Δt is set as 0.0635 seconds and the window length T is 6 seconds. For small earthquakes around M L 3.0, the variation range of the corner frequency is 2 ~ 4Hz. Therefore the window length of 6 seconds is 12 times as long as the reciprocal of minimum corner frequency (1/(f cα ) min = ? ? [4,5) [5,6) [6,7) [7,8) 0/5seconds), which guarantees an enough frequency spectrum's resolution. Meanwhile, the Δt of 0.0635 seconds is less than 1/(2(f cα ) max = 0.125 seconds, so that the high-frequency aliasing can be prevented.
The digitized records and waveforms are then adjusted by the SDDMS system through slant correction, removing DC component of baseline, and arc correction. Afterwards, the Hanning window and fast Fourier transform (FFT) are applied to obtain the spectra of seismic waves. The average Q value of a selected region then can be solved using the method presented above. Also, based on the obtained spectra of seismic waves, the source spectra are obtained after instrument and medium calibration, and the stress drop Δσ can be solved using the method presented above. The Hanning window W(t) is where we assumed the low cutoff frequency ω is 0.25. For the digital records, their spectra of seismic waves can be directly obtained applying the Hanning window and FFT based on a sampling rate 100 dps. Because the digital seismometer is a broadband seismometer, the instrument calibration is not necessary and the source spectra can be obtained after the medium calibration.

Summary
The solved results are displayed in Fig. 14 and Table 8. The relative errors of spectral curves caused from digitalization are about 1~ 2% [12] and the errors of the Q values are in the same magnitude. Thus, the presented differences among the average Q values of different regions are meaningful. The most important factor that causes errors in stress values is the error of earthquake magnitude, whose maximum acceptable value is 0.3. The maximum relative error of stress τ 0 can reach 70% [12]. Since the average stress value of each region is a result of multiple earthquakes, it can be considered that the average stress drops in regions A, B, and C are comparatively higher, next are the stress drops in regions E and F, and the stress drop in region D is the lowest.
From the calculation results, it is observed that the highest average Q value appears in region A's sea area, which is 412. Other high Q values include 348 which appears in the central Jiangsu Province and coastal area of region B; 340 which appears in region F (close to the epicenter of the Nov-9-1996 earthquake); and 305 which appears in Liyang of region C. The lowest average Q value appear in Shanghai and the sea area to the east of Hangzhou Bay (in region D), which is 241. And the Q values in other regions are 289. Comparing the distribution of Q values to the potential seismic source areas, it shows that the higher Q values appear in the potential source areas of higher-magnitude earthquakes. This is because that only under a high Q value, the medium can accumulate high-level stress therefore gives rise to strong earthquakes.
The Q value in region F dropped from 352 to 342 after the Ms6.1 earthquake occurred in Nov-9-1996 on estuary of Yangtze River. Normally, the Q value will drop after the mid-strong earthquakes because of destruction of local medium. Considering measurement errors, the decrease in Q value after that Ms6.1 earthquake is slight, which indicates that the damage caused by the Ms6.1 earthquake in the local medium was small. The current Q value in that area is 340, which is the same as the one right after the Ms6.1 earthquake.
From distribution of the stress drop Δσ, it is found that the drop is high in regions A, B, C, and the next high Δσ appears in regions E and F. The stress drop in region D is comparatively low. From time evolution, it seems that the period of high stress drop corresponds to the time period of high seismicity level in that region.
Based on above discussion, it should be noted the possibility of about Ms6.0 earthquakes in the Southern Yellow Sea and coastal area (regions A and B), the sea area to the east of Yangtze estuary (region F) and Liyang and vicinity (region C).

Conclusions
The conclusions made based on the earthquake prediction can be listed as: (1) Historical earthquakes occurred in the Y-Y seismic zone can be classified into six episodes and each episode can be divided into the main release period and the residual release period. The seismic activity level during the main periods is much higher than that during the residual periods and all earthquakes with Ms ≥ 6.0 occurred in the main release periods. The earthquake occurrence shows an obvious pseudo-periodicity during the main release periods. The duration of a main period is about twice as long as that of the residual release period.
(2) The criteria that indicates the transition from a main release period to a residual release period is: 1) the time interval between the last Ms ≥ 5.0 earthquake during the main release period and the first Ms ≥ 5.0 earthquake of the following residual release period is much longer than the maximum time interval of the Ms ≥ 5.0 earthquakes occurred during the main period; 2) an Ms ≥ 6.0 earthquake may occur about two years before the last Ms ≥ 5.0 earthquake of the main period.
(3) The Y-Y seismic zone is divided into four areas: the sea area of northern sag to the north of 35ºN, the sea area of southern sag to the south of 34ºN and the offshore area of Yangtze estuary, lower Yangtze block to the east of Tanlu fault (PTLF), and Qinling-Dabie block. The seismic activity level in the sea areas is higher than that on the lands and the highest seismicity level is in the sea areas to the south of 34ºN. Except for the current episode, Ms ≥ 5.0 earthquakes have occurred in all four areas during previous episodes.
(4) The current episode, which starts from 1974, likely entered into its residual release period in 1999. That period will end till 2012 and then enter into a seismically quiet episode. During this period, there may have two or three groups of earthquake with Ms ≥ 5.0, but the maximum magnitude will be lower than Ms6.0. According to the possibility order of the occurrence of future mid-strong earthquakes, the four areas can be arranged as: the sea area of southern sag to the south of 34ºN, the lower Yangtze block, the sea area of northern sag to the north of 35ºN, and the Qinling-Dabie block (from high to low).
Finally, based on all above work, we determine three seismic risk zones where mid-strong earthquakeswhich will cause damages may occur within the next 10 years. Those risk zones are: 1) southern sags of Southern Yellow sea and coastal regions, Ms ≥ 6; 2) Yangtze estuary and sea area to the east of it, Ms ≥ 5; 3) Changshu-Taicangarea, Ms ≥ 5 (Fig.  15).