Natural Hazards and Earth System Sciences Investigating the earthquake catalog of the National Observatory of Athens

The earthquake catalog of the National Observatory of Athens (NOA) since the beginning of the Greek National Seismological Network development in 1964, is compiled and analyzed in this study. The b-value and the spatial and temporal variability of the magnitude of completeness of the catalog is determined together with the times of significant seismicity rate changes. It is well known that man made inhomogeneities and artifacts exist in earthquake catalogs that are produced by changing seismological networks and in this study the chronological order of periods of network expansion, instrumental upgrades and practice and procedures changes at NOA are reported. The earthquake catalog of NOA is the most detailed data set available for the Greek area and the results of this study may be employed for the selection of trustworthy parts of the data in earthquake prediction research.


Introduction
Earthquake catalogs are a valuable result of fundamental seismological practice and they form the basis for seismicity, seismotectonic, seismic risk and hazard investigations.
Before one proceeds in such investigations it is essential to examine and report on the spatial and temporal homogeneity and completeness of the catalog.This is because earthquake catalogs are produced by the recording of seismic waves in seismological networks that change in time and space with varying operational practices and procedures.
Greece has the highest seismicity in Europe and the first "Agamemnon" type seismograph was installed in the National Observatory of Athens (NOA) in 1898.In spite of this a Greek seismographic network did not begin its opera-Correspondence to: G. Chouliaras (g.choul@gein.noa.gr)tion until 1964 and at this time also begins the seismological bulletin production and the network earthquake catalog.This catalog spans over four decades of uninterrupted seismological network data and contains more than 75 000 events, so it is arguably the most detailed recent instrumental earthquake catalog of Greece.
A number of scientific papers have used the entire or parts of the NOA catalog, however little work has been reported concerning artifacts in the homogeneity and completeness of the data set.Significant network expansions, infrastructure upgrades and staff changes have taken place in the last four decades of catalog production and the knowledge of the chronological order of the steps towards the improvement in the detectability of NOA's network and the way that earthquake parameters have been reported is of great importance in order to detect more accurately possible seismicity anomalies related to earthquake preparatory processes.
It is the purpose of this investigation to report on such inhomogeneities, artifacts and biases which may distort the earthquake catalog of NOA and this may provide further insight when examining in detail the seismicity patterns of Greece.

History of network development at NOA
The historical development of the permanent seismological station installations by NOA in Greece are reviewed by Bath (1983) and more recently by Papazachos and Papazachou (2003).The first traditional WWSSN station in Greece began operating in Athens at the Institute of Geodynamics (NOA) in 1962 and a year later the second station was installed in Patras, to be followed by the installation of four more stations on the Greek islands: Crete, Kefalonia, Lesvos and Rodos.In 1964 a Wood Anderson seismograph was Published by Copernicus Publications on behalf of the European Geosciences Union.installed in Athens and the respective local magnitude has been used ever since.The seismological network expanded to 13 stations by 1973 and the determination of the source parameters was carried out manually by using appropriate travel time curves.
Advances in the telecommunication infrastructure technology in Greece allowed for data transmission from the remote stations to Athens via leased telephone lines in realtime by 1981.Around this time the rapid increase in seismic data flow, required the change of the operating practices at NOA and a main frame computer with the Hypo 71 computer program (Lee and Lahr, 1975) was employed for the daily analysis and bulletin production procedure.
In 1988 begins the second period of station expansion with the addition of 14 seismological stations by 1990.At this time portable computers were introduced at NOA and the analysis and bulletin production and archiving gradually changed from analog to digital.More time was spent on signal detection and less on routine analysis and this improved the detectability and led to an increase in reporting.The third station expansion occurred around the end of 1994, by the installation of the first digital seismographic network (Chouliaras and Stavrakakis, 1997).A three component, shortperiod, seismographic network comprised of nine remote stations was operated by a dial -up telemetry server at NOA by 1995.
Simultaneous to the analog to digital instrumentation transition period at NOA in 1995, digital signal analysis practices and procedures were introduced.During the next few years, an increase in the research staff occurred and upgrades in the station infrastructure of the NOA network continued with 12 new, three component digital seismic stations that began their operation around 1998.After the catastrophic earthquake on 7 September 1999, and in view of the preparations for the 2004, Olympic Games, the NOA digital seismographic network rapidly expanded around Athens.At the same time NOA participated in various projects coordinated by ORFEUS and EMSC, that  concerned the real time exchange of broad -band waveform data and parametric results for the rapid determination of epicenters (Stavrakakis et al., 1999;Van Eck et al., 2002).The transition from short period to broad band sensors and related software changes lasted until 2004 and this was accompanied by the largest staff increase in NOA's history in 2005 in order to handle the large influx of seismological data.

Seismological data
The monthly bulletins of NOA from 1964 to 2009 have been used to compile the NOA network earthquake catalog for the region 34 • -42 • N and 19 • -29 • E, as described recently by Chouliaras (2009a, b).The seismicity analysis is performed using the ZMAP software which includes a set of tools written in Matlab®, a Mathworks commercial software language (http://www.mathworks.com), with an open code and driven by a graphical user interface (GUI), (Weimer, 2001).
The epicentral distribution of the earthquakes in the NOA network catalog delineate three areas of dense seismicity as seen in Fig. 1: the Hellenic Arc subduction zone to the south and its northwestern extension towards the Ionian islands, the North Anatolian Fault extension into the northeastern Aegean till Central Greece and the Gulf of Corinth.
Figures 2 through 5 show the respective time, hour, magnitude and depth histograms of the data set.From these results the following are pointed out : (a) the significant increases in the registered earthquakes around 1968, 1981, 1995 and 2004, (b) the 20-25 % decrease in the registered earthquakes during the daylight hours due to the increased noise at the   recording stations, (c) the asymmetrical peak around M=3.1 due to the insufficient detectability for small events and d) the depth of the vast majority of events in the catalog is below 50 km, with the most frequent depth range that of 5 to 15 km.
The time variation of the reported magnitudes in the catalog as seen in Fig. 6 shows periodic "bursts" of large earthquakes, around the same periods mentioned earlier for Fig. 2 and the same time periods show up with sudden increases in the cumulative seismicity curve of the entire catalog shown in Fig. 7 (blue line), from a compilation of 75 449 seismic Cumulative seismicity curves.The blue curve is the NOA earthquake catalog and the red curve is for its de-clustered equivalent.events.Large crustal earthquakes usually produce aftershock sequences that contain hundreds to thousands of events thus introducing bias in rate change investigations.Following Chouliaras (2009a), the removal of aftershock clusters from NOA's catalog using the declustering method of Reasenberg (1985) eliminates these sudden increases as shown by the declustered curve (red line) in Fig. 7 which will be used further on in this study.In order to assess the quality of the earthquake catalog the criterion of the magnitude of completeness (Mc) is employed.Mc is defined as the lowest magnitude of the catalog at which all of the events are detected in space and in time (Rydelek and Sacks, 1989;Taylor et al, 1990;Wiemer and Wyss, 2000).Different approaches and methodologies as described by Schorlemmer and Woessner (2008) exist in determining Mc and the assumption of self similarity for the earthquake process, implying that for a given volume a simple power law can approximate the frequency magnitude distribution (FMD), will be adopted here.
The FMD describes the relationship between the frequency of occurrence and the magnitude of earthquakes (Ishimoto and Iida, 1939;Gutenberg and Richter, 1944): Where N is the cumulative number of earthquakes having magnitudes larger than M, a and b are constants.
The b-value describes the relative size of the events and it is determined by the linear least squares regression or by the maximum-likelihood technique: (2) <M> is the mean magnitude of the sample and Mbin is the binning width of the catalog (Aki, 1965;Bender, 1983;Utsu, 1999).Two different Wiemer and Wyss (2000) methods will be used to determine Mc in this investigation: 1.The Maximum Curvature (MAXC), which simply computes the maximum value of the first derivative of the FMD curve.
2. The Goodness of Fit Test (GFT), which compares the observed FMD curve with a synthetic one.A model is found at which a predefined percentage 90% or 95% of the observed data is modeled by a straight line.
The results of the MAXC method in Fig. 8 show that the cumulative frequency magnitude distribution and its first derivative indicate a magnitude of completeness at M=3.1, with a b-value of 1.14.In the previous section of this study a historical account of the upgrading of NOA's seismological network as well as the related changes in analysis procedures and practices was given since these may introduce biases and inhomogeneities in NOA's earthquake catalog.Identification of artificial or  man made seismicity anomalies in earthquake catalogs have been discussed in several studies (Habermann 1982(Habermann , 1983;;Habermann and Wyss, 1984;Wyss and Bufford, 1985;Wyss, 1991;Zuniga, 1989;Zuniga and Wyss, 1995;Zuniga et al., 2000Zuniga et al., , 2005)).
The Genas algorithm (Habermann, 1983(Habermann, , 1987) ) is the appropriate tool to investigate such artificial rate changes and this task is performed on declustered catalogs in order to avoid false alarms from rate changes due to aftershock sequences and clusters.This algorithm identifies significant changes in seismicity rate (number of events larger and smaller than a given magnitude with respect to time) by comparing the mean rate before the time under study to that of the period which follows.The procedure is repeated for increased values up to the end of the seismicity record.Every time a significant change is found, the catalog is marked and splited into two segments which are iteratively analyzed in the same fashion.The result provides the times which stand out as the beginning of periods were increases and/or decreases of seismicity are detected as well as the magnitude range affected by these changes.
Figure 11 shows the output of the Genas algorithm as times of significant rate changes, where circles indicate rate increases and crosses rate decreases.Since the network catalog was in the built up phase around 1966 and 1967 this period may be regarded as inferior to the rest and ignored however significant rate changes around 1973, 1990, 1995 and 2000-2004 are observed.These periods were also mentioned previously as periods of seismological network expansion and upgrading and for this reason we see their imprints in the seismicity catalog as rate changes.These rate changes may or may not be accompanied with changes in the reporting of magnitudes namely "magnitude shifts" (Zuniga and www.nat-hazards-earth-syst-sci.net/9/905/2009/Nat.Hazards Earth Syst.Sci., 9, 905-912, 2009   The temporal variation of the magnitude of completeness for NOA's earthquake catalog. Fig. 13.The temporal variation of the magnitude of completeness for NOA's earthquake catalog.Wyss, 1995) which are crucial in determining the homogeneity and completeness of the catalog.The rate change from the "clean periods" 1973.5 to 1977 is compared to that from 1983 to 1987 in Fig. 12.The cumulative and noncumulative FMD's show a small (11.54%) increase in the seismicity rate and no magnitude shift in the data set, indicating a homogeneous magnitude reporting for the entire period.In a similar fashion the periods before and after the significant increases in the seismicity rate indicated by Genas were successively compared and even though significant rate increases were observed, no apparent magnitude shifts existed.
The following results are found: 1 . 1973.5   The seismological network of NOA has gone through periods of station expansion, instrumentation and staff changes as well as changes in the procedures and practices in data analysis.As described earlier on in this study, these factors usually cause inhomogeneities, artifacts or biases in earthquake catalogs and may be misinterpreted as natural processes.The periods of rate increases around 1973, 1983, 1990, 1995 and 2000-2004 as seen by Genas coincide with periods of significant station increases in NOA's seismological network.These rate increases are also attributed to the instrumentation and software upgrading as well as staff increased of the NOA network.The results show that such changes affect the seismicity rate of the data base equally in all magnitude ranges.
In a recent investigation Uyeda et al. (2009) used the Japan Meteorological Agency earthquake catalog to investigate the idea of "natural time" as first presented by Varotsos et al. (2005a, b) in order to forecast the occurrence time of an impending earthquake.The high density of seismological network stations in Japan allowed that investigation of a relatively low threshold magnitude of M=2 in their catalog analysis.The results obtained after considering the seismicity subsequent to a Seismic Electric Signals activity they recorded (which was similar to the ones recorded in Greece, e.g., see Varotsos et al., 2003), succeeded in predicting the occurrence of an M=6.0 event within a time window of a few days.
The latter is a characteristic example showing that the detectability of the seismological network is an important factor in searching for critical stage seismicity.In view of this fact and since NOA's earthquake catalog for Greece is also used extensively in on going earthquake prediction research, the properties of this catalog were investigated.Among others, the chronological order of changes in the network station infrastructure was studied and it is found that they influence the results of the analysis procedures.Thus the present study enables the selection of appropriate parts of the data set for trustworthy analysis.
Edited by: M. E. Contadakis Reviewed by: two anonymous referees

Figure 1
Figure 1 Earthquake epicenters (red) between 1964 and 2009 from the earthquake catalog of NOA.

Figure 2
Figure 2 Time histogram of NOA's earthquake catalog.Arrows indicate times of increased seismic activity.

Fig. 2 .
Fig. 2. Time histogram of NOA's earthquake catalog.Arrows indicate times of increased seismic activity.

Figure 3
Figure 3 Hourly variation of the seismicity in NOA's earthquake catalogue.Arrow indicates maximum day/night variation.

Fig. 3 .
Fig. 3. Hourly variation of the seismicity in NOA's earthquake catalog.Arrow indicates maximum day/night variation.

Figure 5 Depth
Figure 5 Depth Histogram of NOA's earthquake catalog.

Fig. 7 .
Fig. 7. Cumulative seismicity curves.The blue curve is the NOA earthquake catalog and the red curve is for its de-clustered equivalent.

Figure 8 Cumulative
Figure 8 Cumulative Frequency -Magnitude Distribution and it's first derivative (green line) for NOA's earthquake catalog.Crosshair is at a value of Mc=3.1.The b-Value of 1.14 and its error (+/-0.03)are determined with the weighted least squares method.

Fig. 8 .
Fig. 8. Cumulative Frequency -Magnitude Distribution and its first derivative (green line) for NOA's earthquake catalog.Crosshair is at a value of Mc=3.1.The b-value of 1.14 and its error (±0.03) are determined with the weighted least squares method.

Figure 9
Figure 9 Results for the Goodness of Fit Test (GFT) for NOA's earthquake catalog.The arrow indicates a 90% confidence for a magnitude of completeness (Mc) of 3.1.

Fig. 9 .
Fig. 9. Results for the Goodness of Fit Test (GFT) for NOA's earthquake catalog.The arrow indicates a 90% confidence for a magnitude of completeness (Mc) of 3.1.

Figure 10 Map
Figure 10 Map of the spatial distribution of the Magnitude of Completeness (Mc) of NOA's earthquake catalog.

Fig. 10 .
Fig. 10.Map of the spatial distribution of the Magnitude of Completeness (Mc) of NOA's earthquake catalog.

b
= log 10 (e)/[< M > −(Mc − Mbin/2)] Figure 9 shows the GFT results which indicates that Mc has a 90% and 95% confidence for the values of 3.1 and 3.5, respectively.These two results in combination show that a value of Mc=3.1 may be chosen as representative of the investigated data set, however we must consider the study of Woessner and Wiemer (2005) which indicates that two methods used in this study may underestimate Mc by about 0.2.The differences in Mc as a function of space is influenced by the seismological network configuration as shown by Chouliaras (2009a, b) and the spatial mapping of Mc may identify regions in outer margins of the network that give radically different reporting and should not be used in quantitative studies.This methodology is applied to the NOA network catalog data and Fig. 10 shows the spatial variability of Mc in Greece.The lowest Mc region is that around Athens with a Mc around 2 and this value increases outwards, with values Mc values around 3 for Central Greece, Peloponesus and Northern Greece and values around 4 or more in the bordering regions all around Greece where NOA's seismic network is sparse (http://www.gein.noa.gr/services/net figure.gif).

Figure 11
Figure 11 Genas algorithm result for NOA's de-clustered earthquake catalog.Circles indicate rate increases and Crosses rate decreases.

Figure 12 A
Figure12A comparison of the cumulative and non-cumulative Frequency -Magnitude distributions from two successive periods of NOA's earthquake catalog.The period of1973.5 -1977  is compared to 1983 -1987.

Fig. 12 .
Fig. 12.A comparison of the cumulative and non-cumulative Frequency -Magnitude distributions from two successive periods of NOA's earthquake catalog.The period of1973.5-1977 is compared  to 1983-1987.
4. 1973.5 to 1977 is compared to 2005 to 2009 with a rate increase of 550% The temporal variation of the magnitude of completeness Mc, for the data set of the NOA network is seen in Fig. 13.It is expected that Mc changes with time in data sets from expanding networks with inhomogeneous reporting practices and procedures and this is also observed here where Mc decreases from 3.8 to 3.0 around 1970 and also from 3.2 to 3.0 around 1982.This result is also confirmed by the Genas results of Fig. 11 and is attributed to significant increases (almost doubling) in the number of seismological stations comprising the NOA network.The unexpected increase in the Mc value from 3.0 to 3.2 around the year 2000, in a period of station expansion and software upgrades will be explained in a separate study regarding the detectability of NOA's network.3 Discussion and conclusion The earthquake catalog of NOA for the Greek area since the beginning of its seismological network installation in 1964 is analyzed in this study.The uninterrupted operation and fundamental seismological practice at NOA during the last four decades provided for a data base of more than 75 000 seismic events until 2009.Statistical analysis of this catalog using two different methods indicates a magnitude of completeness Mc=3.1 as an indicator of the detectability of network and a frequency-magnitude relation with a b-Value of 1.14.These values are shown to vary spatially in Greece and are strongly influenced by the seismological network configuration as well as local seismicity.