1 A MAGNITUDE INDEPENDENT SPACE-TIME EARTHQUAKE CLUSTERING ALGORITHM (MISTIC)Leptokaropoulos M. Kostantinos Seismology and Physics of the Earth's Interior, Institute of Geophysics, Polish Academy of Sciences, Warsaw, Poland Gkarlaouni G. Charikleia Aristotle University of Thessaloniki, Department of Geophysics, School of Geology, Thessaloniki, Greece This presentation introduces of a magnitude independent algorithm for identifying space-time earthquake clusters. 14th International Conference of the Geological Society of Greece 25-27 May, 2016 Thessaloniki

2 Outline of the presentationSamos-Karaburun Efpalio Sequence Outline of the presentation Existing Declustering Techniques and Motivation Introduction to the algorithm Testing and Evaluation Summary In this presentation there is a brief reference to the theoretical background and the existing declustering techniques and a discussion on the concept that motivated the study. There is an introduction of the code’s process and two applications are demonstrated in the case of Samos- Karaburun seismicity and in the case of Efpalio sequence in 2010 in the western Corinth Gulf. Finally results are summarized. EGE 2016

3 Declustering algorithmsThe methods for discriminating the independent from the dependent fraction of seismic activity are mainly divided into two approaches (van Stiphout et al., 2012 and references within): Conventional Stochastic Window-based and Link-based methods Utsu (1969) Gardner and Knopoff (1974) Reasenberg (1985); Frohlich and Davis (1990) A branching point process Kagan (1991); Ogata (1998); Console and Murru (2001); Zhuang et al. (2002, 2004) There is a number of methodologies that they are used for distinguishing the independent from the dependent fraction of seismicity. Independent seismicity, is assumed to be constant over large time periods, in a given constant fault stressing rate. However, the dependent seismicity is the triggered one, such as aftershock sequences and oftentimes they have to be removed in certain types of analysis. Generally they are developed following either a conventional or a stochastic approach. The conventional methodologies are window-based or linked-based. The window-based methods remove the smaller magnitude earthquakes around a larger event in a space-time window (main shock (e.g. Utsu, 1969; Gardner and Knopoff, The larger the main shock, the bigger the window size is considered). The link-based methods remove events which are found within a compromised space-time interval compared to an earlier event. On the other hand, stochastic declustering methods mostly model space-time-magnitude occurrences, in the form of a branching point process. Most recently there are some more ideas based on short time temporal fluctuations, seismicity rates and multiparametric equivalent dimensions. Temporal fluctuations are noticeable even in short time scales (Hainzl and Ogata (2005) Seismicity rate (CURATE, Jacobs et al., 2013) Multi parameter space clustering after data transformation to equivalent dimensions (Lasocki, 2014). EGE 2016

4 Magnitude Independent Space -TIme Clustering algorithmaim The aim of this study is not to discriminate and remove aftershocks from background seismicity, but to identify clusters in space and time without any magnitude constraints, seeking for seismic activity which exhibits enhanced rates with no obvious relation to a main shock. " MIstic" Algorithm Magnitude Independent Space -TIme Clustering algorithm identifies clusters in the sense that clustered seismicity is more intense in both space and time domains than the average (or background) activity The aim of this algorithm is not to discriminate and remove aftershocks from background seismicity, but to identify earthquake clusters in space and time without any magnitude constraints (except the completeness level of the dataset) seeking for seismic activity which exhibits enhanced rates with no obvious relation to a main shock. There is strong evidence supporting the fact that the magnitude of each individual aftershock is independent of the mainshock magnitude So, a magnitude independent algorithm was implemented named MISTIC from the initials of the words. In this way, swarms and relatively smaller magnitude seismic sequences can still be identified and thereinafter be connected with physical processes. It is based upon the well-established observation when seismicity rate is more than 10 times higher than the average, it is unlikely that the cluster’s gravity center is significantly shifted from the actual location This algorithm is magnitude independent, since high rated seismicity concentrated in a narrow region can be observed in the absence of a characteristic earthquake magnitude, something which is verified in the demonstrated application. EGE 2016

5 "MIstic" Algorithm Tmax - Maximum Inter-event time between subsequent events Xmax - Maximum Distance between the epicentres and the gravity center of the earthquake cluster Nmin - Minimum Number of events for each cluster The cluster analysis procedure depends on 3 Parameters: The Source Code is analyzed into 3 separate Steps: MYSTIC Spatial Test 1 Test 2 Temporal The Algorithm Preliminary Criteria Secondary Tests Final Criterion Implemented in Matlab The cluster analysis depends on three parameters, that they are user defined: The maximum inter-event time between subsequent events, the maximum distance between earthquakes epicentres and the potential clusters’ s gravity center and the minimum number of events per cluster. The source code is implemented in Matlab into 3 steps i) In the beginning there is a preliminary temporal constraint, ii) Secondly, a spatial constraint which is performed either manually or with the use of 2 tests and in the end there is a final temporal criterion for the verification of the results. EGE 2016

6 The preliminary temporal criterionMYSTIC Spatial Test 1 Test 2 Temporal The Algorithm Preliminary criteria Second Tests Final Criteria The preliminary temporal criterion Inter-event time between events is the foremost criterion, due to the dense occurrence of the seismic burst A primary spatial constraint might be misleading 3 seismic events of a sequence (ni-1, ni, ni+1) described by t, r, Δt and Δx Δt1, Δt2 and Δt3 are similar Δx1>>>Δx3 and Δx2>>>Δx3 Clusters with N>Nmin are further investigated In the case of a seismic excitation, the clustering is dense the time differences between the successive events are small. The code puts some primary constraints in interevent time. In any other case, that an inter-event distance criterion is adopted first, there is a likelihood for a misleading result. For example, we assume that there is a sequence of three seismic events, ni-1, ni ,ni+1 whose epicentres are shown in this figure. These events are determined by an origin time (tj) a position vector (rj) 2D or 3D. There are events close in time. However, ni-event although it has occurred close in time with the former (ni-1) and the following (ni+1) events (Δt1, Δt2 and Δt3 are similar), it is located far away from them, such that Δx1>>>Δx3 and Δx2>>>Δx3. An inter-event distance criterion could classify both ni, ni+1 events out of the cluster although the distance between ni+1 and ni-1 is relatively much shorter. On the contrary, a predefined temporal criterion ensures that such events are not removed from the cluster. In this way, many non-clustered events in the datasets are created, however none of the clustered event is subtracted, as an outlier After the first criterion is applied, the compiled catalogs are tested for the number of events they contain and only those for N>Nmin, are further investigated through the second criterion. EGE 2016 5/17

7 MYSTIC Spatial Test 1 Test 2 Temporal The Algorithm Preliminary criteria Second Tests Final Criteria The spatial criterion The spatial constraints are successively performed in 2 individual steps: The events whose epicentres are far from the clusters gravity center in a distance Xi>Xmax are identified and removed. If the spatial constraints are of minor importance, the algorithm identifies and removes the outliers by providing the option for the application of two different Tests: Test 1 and Test 2 Proceeding to the second step we set some spatial criteria in order to define clustering in the space domain. This is accommodated in 2 steps Firstly, The gravity center for the clusters defied in the first procedure is calculated. Earthquakes that are found in a distance X-i>Xmax which is specified by the user are identified and removed. This criterion is adopted in order to manually select a radius that is in agreement with the scope and the needs of the analysis. Secondly, in the case that the spatial constraints are of minor importance or not strict enough, the algorithm identifies and removes the outliers by providing two alternative option These are called Test 1 and Test 2 EGE 2016 6/17

8 The spatial criterion Test 1 Test 2MYSTIC Spatial Test 1 Test 2 Temporal The Algorithm Preliminary criteria Second Tests Final Criteria The spatial criterion Test 1 Identifies events in a distance equal to or larger than the average distance of all events from the cluster center, plus k times the standard deviation of these distances (σ) Test 2 Defines the minimum distance from the center of the cluster according to the formula: The first option is to identify events which lie in a distance equal or greater than the average distance ( X) between the events and the center of the cluster, plus k times the standard deviation (σ) of these distances from the cluster’s gravity center. The alternative option assumes that the minimum distance between two events which belong to the same cluster can be defined from the center of the cluster. The first term in the right part of equation 2 balances the effect of the outliers, since Xmsx is the average distance of the 5% of the most distant events from the center. The second term is proportional to the dispersion of the data Ωισθαλιζατιον EGE 2016 7/17

9 The spatial criterion Example:MYSTIC Spatial Test 1 Test 2 Temporal The Algorithm Preliminary criteria Second Tests Final Criteria The spatial criterion Example: Cluster identification for the same seismic dataset for different spatial and temporal criteria (Test 1 - blue circles) (Test 2 - red circle) 2 the spatial distribution of earthquake epicentres which form a cluster as they have been extracted from the "MISTIC" code is shown. In each plot, the radii of the inner and the outer blue circles are equal to the mean distance of the epicentres from the cluster center plus 2σ and 3σ, (i.e. k=2 and k=3) respectively When Xmax takes relatively low values, Test 2 provides almost the same constraint with the case of 2σ although it is always stricter. As Xmax increases and the outliers lie in greater distances, Test 2 becomes sufficiently stricter than 2σ+ and prevents events that are found in intermediate distances from intruding into the cluster. An alternative way to avoid such a situation is to select a reasonably low value of Xmax (and Tmax) from the beginning of the clustering procedure so that the outliers can be initially filtered. These values can be adjusted by the user after repeating the process. In particular, in Fig. 2b, there are three events inside the red circle which have been excluded from the cluster according to the second temporal criterion described in the following section (Step 3). EGE 2016 8/17

10 The FINAL TEMPORAL criterionMYSTIC Spatial Test 1 Test 2 Temporal The Algorithm Preliminary criteria Second Tests Final Criteria The FINAL TEMPORAL criterion Τhere is a further filter with a final temporal constraint, because after the execution of the Spatial Criterion the removal of some remote events resulted to the increase of the inter-event times between subsequent earthquakes in the specified clusters, so they are compared with Tmax Output Catalogs 1) The first catalog extracted after the application of the Preliminary Temporal Criterion (Step 1) 2) The second catalog derived after the application of the Spatial Criteria (Step 2) 3) The final catalog (after employing Steps 1, 2 and 3) The datasets compiled after the application of the temporal and the spatial criteria, are further filtered with a final temporal constraint. The reason for adding a supplementary test appears because after the execution of the Spatial Criterion the removal of some remote events resulted to the increase of the inter-event times between subsequent earthquakes in the specified clusters. At the end of the process, the final catalogs with the extracted clusters are available for visualization and any other processes for the inspection of the method’s efficiency or for the final decision making. This is a snapshot of the visualization of the clusters EGE 2016 9/17

11 Testing and evaluation 1.samos-Karaburun areaData: events - Local network ( ) Tan et al. (2014) along with data provided from HUSN. The first application is performed on the seismicity of the area Samos – Karaburun area in the eastern Aegean. The seismicity data are obtained from a locally established network along with the recordings of the Hellenic Unified seismological network. The duration of the observations is 6 years, during a seismic excitation that took place. Και επειδη εχει πολύ μικρο σφαλμα και χαμηλη πληρότητα- The completeness magnitude is calculated with the algorithm of Leptokaropoulos et al by applying a maximum likelihood goodness of fit test, is low enough and enables us to seek for minor magnitude seismic clusters. The b-value of Gutenberg-Richter relation was found equal to 0.93, a value close to 1.00 which is typica The calculation of the catalog completeness magnitude, MC, was accomplished by applying a goodness of fit test (Leptokaropoulos et al., 2013) which is a modified version of Wiemer and Wyss (2000) methodology. EGE 2016 10/17

12 Testing and evaluation 1.Samos-Karaburun area clusters propertiesQuantitative properties of the clusters: (Tmax=0.5 days, Xmax=50km and Nmin= 30) Cluster ID Number of Events Duration (Days) M1 Mmax ΔΜ between the 2 strongest events Rank in Sequence b-value (95% confidence bounds) C1 141 6.12 2.1 4.1 0.3 103/141 0.76±0.13 C2 37 2.07 3.1 4.3 0.9 2/37 0.61±0.20 C3 50 3.17 3.3 3.6 0.1 43/50 0.75±0.21 C4 30 1.48 2.5 15/30 0.87±0.31 C5 47 1.83 5.1 1.9 1/47 0.88±0.25 C6 69 4.16 4.8 0.6 1/69 0.70±0.17 C7 94 4.51 1.7 3.8 19/94 0.80±0.16 C8 232 25.79 1.8 5.0 0.8 19/232 0.89±0.11 C9 105 4.2 2.7 0.2 65/105 1.08±0.21 C10 97 5.21 1/97 0.67±0.13 C11 136 31.8 2.6 3.7 0.5 88/136 0.73±0.13 (Mc=1.8) Green cells accommodate main-shock aftershock clusters, whereas Red cells show swarm-like clusters Here you can see the clusters that were identified from the algorithm for the specific parameters. 11 clusters fulfill these criteria. Some of these clusters exhibit main shock-aftershocks characteristic (shaded with blue), whereas some others (C1, C3, C4, C7, C9 and C11) are rather swarm-like sequences. This discrimination is based on the magnitude difference between the two strongest events (<0.5 units for swarms) and the ranking of the occurrence of the strongest shock at the earliest stages of the seismic burst. EGE 2016 11/17

13 Leptokaropoulos et al. (2016)Testing and evaluation 2.Efpalio doublet Clustering and coseismic stress changes Data: August 2008– December 2012 Data: April 2010 – December 2012 Gray triangles: Non- clustered activity White circles: Spatio-temporal clusters Almost all clustered epicenters are located in increased stress areas In the second example the code is applied at the aftershock seismic activity of the Efpalio seismic doublet. The coseismic stress changes due to the occurrence of the strongest earthquakes are shown in the color scale. The clusters with circles are the ones defined with the use of the code. Results show that all of them are located in positive stress zones, where seismic activity is increased. Although several events occurred inside negative lobes, only spatiotemporal clusters characterized by increased seismicity rates were located in stress-enhanced areas. EGE 2016 Leptokaropoulos et al. (2016) 12/17

14 Summary "MISTIC" is a magnitude independent space-time clustering algorithm, which identifies seismic clusters with significantly increased seismicity rates and regardless of the existence of a characteristic magnitude earthquake Extracted clusters are well determined and tightly shaped "MISTIC" provides: Fast and flexible calculations even for relatively large catalogs Supplementary information is provided about statistics of the original dataset and the derived clusters (total number of events, time span of the catalog, mean - median inter-event time and area dimensions, starting and ending date/time, number of events, maximum magnitude difference) Visualization of clusters and lists of events potentially belonging to each cluster User friendly (direct testing and controlling by the user, anytime) When is "MISTIC" suitable? Large areas with multiple fault zones interacting on regional scales (more than several hundreds of km) or global seismicity datasets are not suitable Summing up, we introduced MISTIC, a magnitude independent space –time clustering source code! Its characteristic is that it defines earthquake clusters regardless of the existence of a characteristic magnitude, but only when seismicity rates are significant. Testing and up to now performance has demonstrated that the extracted clusters are well determined and tightly shaped! Regarding the code running: The calculations are performed fast even for big data sets, [the computation time is less than 10sec, even for a dataset of approximately events] So it grants the user with the flexibility to test many parameters or value combinations. Additionally there is supplementary information, automatically calculated, on the statistical properties both of the original dataset and the extracted clusters. The process is visualized, figures of clusters and lists of events potentially belonging to each cluster are provided. Finally it is considered user friendly, since the user can manually interfere anytime in the process Its worth to say that MISTIC is rather not accurate for large regions or for global seismicity catalogs Of course, further testing will be beneficial for the verification of the process and will certainly contribute to the optimization of the program. EGE 2016 13/17

15 EGE 2016 references Aktar, M., H. Karabulut, S. Özalaybey, and D. Childs (2007), A conjugate strike-slip fault system within the extensional tectonics of Western Turkey, Geophys. J. Int., 171, Benetatos, C., A. Kiratzi, A. Ganas, M. Ziazia, A. Plessa, and G. Drakatos (2006), Strike slip motions in the Gulf of Siğaçik (western Turkey): Properties of the 17 October 2005 earthquake seismic sequence, Tectonophysics, 426, 263 – 279. Bozkurt, E. (2001), Neotectonics of Turkey - a synthesis, Geodin. Acta, 14, 3-30. Bozkurt, E. (2003), Origin of NE-trending basins in western Turkey, Geodin. Acta, 16, Console, R., and M. Murru (2001), A simple and testable model for earthquake clustering, J. Geophys. Res., 106, 8699– 8711. Çiftçi, N.B., and E. Bozkurt (2009), Pattern of normal faulting in the Gediz Graben, SW Turkey, Tectonophysics, 473, , doi: / j.tecto Felzer, K. R., R. E. Abercrombie and G. Ekström (2004), A common origin for aftershocks, foreshocks and multiplets, Bull. Seismol. Soc. Am., 94, Frohlich, C., and S. D. Davis (1990), Single link cluster analysis as a method to evaluate spatial and temporal properties of earthquake catalogues, Geophys. J. Int., 100, Gardner, J., and L. Knopoff (1974), Is the sequence of earthquakes in Southern California, with aftershock Removed, Poissonian?, Bull. Seismol. Soc. Am., 64, Habermann, R. E., and M. Wyss (1984), Background seismicity rates and precursory seismic quiescence: Imperial Valley, California, Bull. Seismol. Soc. Am., 74, 1743–1755. Hainzl, S., and Y. Ogata (2005), Detecting fluid signals in seismicity data through statistical earthquake modeling, J. Geophys. Res., 110, doi: /2004JB Inan, S., S. Ergintav, R. Saatçilar, B. Tüzel, and Y. İravul (2007), Turkey makes major investment in earthquake research, EOS Trans. Am. Geophys. Union, 88, Kagan, Y. Y. (1991), Likelihood analysis of earthquake catalogues, Geophys. J. Int., 106, 135–148. Kagan, Y. Y., and D. D. Jackson (1991), Long-term earthquake clustering, Geophys. J. Int., 104, Lasocki, S. (2014), Transformation to equivalent dimensions – a new methodology to study earth-quake clustering, Geophys. J. Int., 197, Leptokaropoulos, K. M., V. G. Karakostas, E. E. Papadimitriou, A. K. Adamaki, O. Tan, and S. İnan (2013), A homogeneous earthquake catalogue compilation for western turkey and magnitude of completeness determination, Bull. Seismol. Soc. Am., 103, 5, Leptokaropoulos, K. M., E. E. Papadimitriou, B. Orlecka-Sikora, and V. G. Karakostas (2016), An evaluation of Coulomb stress changes from earthquake productivity variations in the Western Gulf of Corinth, Greece, Pure Appl. Geophys., doi: /s Michael, A. J. and L. M. Jones (1998), Seismicity alert probabilities at Parkfield, California, revisited, Bull. Seismol. Soc. Am., 88, Mountrakis, D., A. Kilias, E. Vavliakis, A. Psilovikos, and E. Thomaidou (2003), Neotectonic map of Samos Island (Aegean Sea, Greece): Implication of Geographical Information Systems in the geological mapping. In: 4th European Congress on Regional Geoscientific Cartography and information Systems, Bologna, Italy, Ogata, Y. (1998), Space-time point-process models for earthquake occurrences, Ann. Inst. Stat. Math., 50, 379– 402. Omori, F. (1894) On aftershocks of earthquakes, J. College Sci. Imp. Univ., Tokyo 7, 111–200. Reasenberg, P. (1985), Second order moment of central California Seismicity, , J. Geophys. Res., 90, 5479–5495. Tan, O. (2013), The dense micro-earthquake activity at the boundary between the Anatolian and South Aegean microplates, J. Geodyn., 65, Tan, O., E. E. Papadimitriou, Z. Pabucçu, V. Karakostas, A. Yörük, and K. Leptokaropoulos (2014), A detailed analysis of microseismicity is Samos and Kusadasi (Eastern Aegean Sea) areas, Acta Geophys., doi: /s Utsu, T. (1969), Aftershock and Earthquake Statistics (I): Some Parameters Which Characterize an Aftershock Sequence and Their Interrelations, J. Fac. Sci, Hokkaido Univ. Ser. VIII, Geophys., 3, van Stiphout, T., J. Zhuang, and D. Marsan (2012), Seismicity declustering, Community Online Resource for Statistical Seismicity Analysis, doi: /corssa Wessel, P., and W. H. F. Smith (1998). New, improved version of the Generic Mapping Tools Released, EOS Trans. AGU, 79, 579. Wiemer, S., and M. Wyss (2000), Minimum magnitude of completeness in earthquake catalogues: Examples from Alaska, the Western United States, and Japan, Bull. Seismol. Soc. Am., 90, 4, Zhuang, J., Y. Ogata, and D. Vere-Jones (2002), Stochastic declustering of space–time earthquake occurrences, J. Am. Stat. Assoc., 97, 369–380. Zhuang, J., Y. Ogata, and D. Vere-Jones (2004), Analyzing earthquake clustering features by using stochastic reconstruction, J. Geophys. Res., 109, B05301, doi: /2003JB 16/17

16 Thank you acknowledgementsThe authors would like to thank the editor and two anonymous reviewers for their valuable comments and suggestions. Maps were made with the QGIS Geographic Information System (http://qgis.osgeo.org) and GMT The magnitude distribution parameters were estimated by the "McCalc" Program, which can be downloaded from the IS-EPOS web platform (https://tcs.ah-epos.eu/) 14th International Conference of the Geological Society of Greece 25-27 May, 2016 Thessaloniki