Probabilistic Seismic Hazard Maps for the North Balkan Region

R.M.W. Musson

Abstract - A set of seismic hazard maps, expressed as horizontal peak ground acceleration, have been computed for a large area of Central and Eastern Europe covering the North Balkan area (Former Yugoslavia, Hungary, Romania). These are based on (a) a compound earthquake catalogue for the region; (b) a seismic source model of 50 zones compiled on the basis of tectonic divisions and seismicity, and (c) a probabilistic methodology using stochastic (Monte Carlo) modelling. It is found that the highest hazard in the region comes from intermediate focus earthquakes occurring in the Vrancea seismic zone; here the hazard exceeds 0.4 g at return periods of 475 years. Special account has been taken of the directional nature of attenuation from this source.

Key words: Seismic hazard maps, Pannonian Basin, Monte Carlo simulation, Vrancea, intermediate focus earthquakes, seismic hazard methodology.
 
 


Introduction



A study of probabilistic seismic hazard mapping for the Northern Balkan area was conducted as part of the Copernicus Project "Quantitative Seismic Zoning of the Circum Pannonian Basin" (EC Project CIPA-CT94-0238). In this project, seismic hazard analyses of this area were undertaken using two methodologies, a deterministic approach based on numerical synthesis of ground motion, and a probabilistic approach based on analysis of the regional earthquake catalogue through the medium of a seismotectonic seismic source zone model. This can be considered as a hybrid approach, in which the results for the two techniques can be compared and contrasted. In this present paper, maps are produced using probabilistic methods, the deterministic results being discussed elsewhere (Zivcic et al 1999, Bus et al 1999, Markusic et al 1999, Radulian et al 1999). Although this Copernicus project CIPA-CT94-0238 had no formal links with the Global Seismic Hazard Assessment Programme (GSHAP) project, the results of the probabilistic component of the project were made available to GSHAP to cover what would otherwise have been a gap in the geographical coverage.
 
 


Probabilistic and stochastic methodologies for seismic hazard









The probabilistic method (often referred to as PSHA, standing for probabilistic seismic hazard assessment) follows these steps:

  1. The area to be studied is divided into discrete seismic source zones, each of which is deemed to be uniform in the character of its seismicity. There should be an equal probability that an earthquake of given magnitude could occur it any place within a single source zone. The geometry of the source zones is determined by the analyst from a combination of tectonic, geophysical, geological and seismological data.
  2. The seismicity within each source zone is analysed, using the local earthquake catalogue, in order to determine the magnitude-frequency parameters and the maximum magnitude.
  3. A locally appropriate attenuation relationship is chosen, to relate the expected ground motion at site during an earthquake to the magnitude of the earthquake and its distance from site. The uncertainty or scatter of the ground motion values about the predicted mean is an important variable which must also be used in the analysis.
  4. The hazard analysis is based on the fact that the probability that an earthquake of magnitude m occurs in a source zone within any given distance interval is proportional to the fraction of the area of the zone that occurs within this range of the site. Since each source zone is deemed to be homogenous, the fractional occurrences expected in any small sub-area of the zone can easily be calculated. An analytical integration is performed over all ground motion values, magnitudes, and source zones. From the results it is possible to determine the probability of any acceleration value being exceeded, assuming earthquake occurrences to follow a Poisson distribution.
The method is discussed in more detail elsewhere, for example, Reiter (1990).

One important issue is how to treat uncertainties in the basic parameters of the seismicity distribution. The parameters of the magnitude frequency curve can be estimated from historical catalogues of earthquakes, but with a margin of uncertainty which may be significant. The maximum magnitude parameter is also uncertain, as no value can be proved to be a bounding limit. One means of dealing with uncertainties has been found following the work of Coppersmith and Youngs (1986), which is to use a logic tree in which different values for key parameters are assigned, with different weights attached to each. The problem with this is that the choice of values and weights can be subjective and open to debate.

An alternative is to use a stochastic modelling approach, also known as Monte Carlo simulation. In this approach, the initial three stages of the PSHA method as described above are followed in exactly the same way. The only difference is in the method by which the probabilities are calculated. Instead of an analytical integration, the seismic source zone model is used to generate a large number of synthetic earthquake catalogues, each having the same properties as the historical earthquake catalogue, but with random occurrences of earthquakes following the Poisson model, with epicentral locations randomly determined within each source zone. This process generates a very large number of synthetic observations at site, and from these observations probabilities and return periods can be calculated directly in a simple and straightforward way. This technique is not new, and has been used in a number of seismic hazard studies in different parts of the world, as in studies by Rosenhauer (1983), Shapira (1983), Johnson and Koyanagi (1988), Ahorner and Rosenhauer (1993), etc.

This method has a number of advantages over conventional PSHA techniques. The first is the powerful way in which one can handle uncertainties in numerical parameters (such as activity rates or b values) - these can be modelled as distributions rather than discrete alternatives as in the logic tree approach. Each iteration of the program can sample a new value from the distribution. Secondly, special cases, such as non-Poissonian behaviour or local variations in attenuation can be dealt with very easily, whereas implementing such irregularities in an analytical PSHA can be difficult. Thirdly, on a non-technical note, the method can be very easily explained to, and understood by, persons outside the seismological community such as politicians and planners who may have little mathematical ability. When it comes to the use of seismic hazard results for planning and policy making, it may be an advantage if the persons using the results have understanding of, and therefore increased confidence in, the methods used to obtaining those results.

The disadvantages are twofold: firstly, the most extreme values generated in the simulation process may vary from run to run, giving slightly different results. This is not a significant problem as long as sufficient numbers of synthetic catalogues are used, and it can be demonstrated that with adequate numbers of simulations the results are stable. Secondly, the method is computationally intensive; this is not a problem today but might have been a factor in the lack on interest in the technique in the 1980s when adequate computing power was less easily available.
 
 




The earthquake catalogue











A primary component of any probabilistic seismic hazard study is always the earthquake catalogue underlying the analysis. The preparation of the earthquake catalogue for this study has been documented in Musson (1996). A brief summary will be given here.

In the absence of a reliable, homogenised catalogue for the whole region, based on original data treated in a uniform way, it was necessary to compile a working catalogue from readily available sources. The regional limits of this working catalogue are from 12-30 East and 42-50 North. The parameters listed in the catalogue are as follows: day, month, year, hours, minutes, seconds, latitude, longitude, depth, Ms, mb, ML, agency code; the minimum magnitude represented is 4. The total number of earthquakes in the working file is 3946. The start date for the catalogue is nominally 1000; in fact the first earthquake occurs in 1022. The last event in the catalogue is dated 9 March 1993. From statistical analysis, the catalogue as a whole is believed to be complete since 1885. For events larger than 5.5 Ms it is complete for another hundred years before that, and for the largest events (? 6.0 Ms) it is complete since 1590.

Since the publication of the catalogue minor revisions were undertaken to correct errors that came to light subsequent to the report of Musson (1996). These mostly involved earthquakes on Romanian territory. Some other new regional catalogues were investigated (e.g. Herak et al 1996), but these new catalogues are themselves compilations, and there is a limit to the merits of producing compilations of compilations. For some countries in the region it is understood that new improved national catalogues have been compiled which are being kept confidential. Obviously, there is no possibility to take account of such catalogues. While the catalogue used in this study is not considered perfect as regards historical seismology, it is sufficiently robust, particularly with respect to the period since 1885, which informs the bulk of the seismicity analysis, to be adequate for the seismic hazard task for which it was compiled.
 
 


Identification of foreshocks and aftershocks











Any complete earthquake catalogue is clearly non-Poissonian: earthquakes are not entirely time-independent events because any substantial earthquake is usually followed by a cluster of aftershocks whose occurrence is dependent on the appearance of the main shock. Since probabilistic analysis of seismic hazard relies on the assumption that seismicity follows a Poisson process, it is generally considered essential to remove any non-Poissonian behaviour from earthquake catalogues.. If only main shocks are considered, then it has been found that the earthquake behaviour for reasonably large areas is generally found to be described satisfactorily by the Poisson model (eg Gardner and Knopoff 1974), in which case the use of hazard estimation models that assume a Poisson model is justified. The effect on the hazard estimation caused by the elimination of aftershocks from consideration is generally regarded as unimportant or acceptable on the grounds that aftershocks are an order of magnitude smaller than main shocks.

Since it is not valid to derive recurrence statistics from the complete catalogue and apply this to predicting mainshock occurrence, it is necessary to decluster the catalogue by removing all aftershocks and foreshocks (collectively referred to as accessory shocks). Otherwise, one will obtain incorrect estimates of the probability of large main shocks, since the slope of the magnitude-frequency curve will be affected by the appearance of many small events which are not main shocks (in effect, the removal of aftershocks makes the magnitude-frequency curve less steep).

To do this is not entirely a straightforward procedure. As is remarked by Reasenberg and Jones (1989), "aftershocks can only be identified in a statistical fashion: they bear no known characteristics differentiating themselves from other earthquakes". There are two ways in which this can be done: by manual inspection or algorithmically. To successfully remove accessory shocks by hand requires, ideally, first hand knowledge of the earthquake catalogue on an event-by-event basis. This was practical, for example, in the hazard mapping study by Musson and Winter (1997) for the UK, where one of the investigators was also the author of the earthquake catalogue used (Musson 1994). Also, in the low seismicity context of the UK the total number of events to be considered is not very considerable. In most cases, however, it is preferable to use some computational method for identifying accessory shocks, especially where earthquake catalogues are large.

The method used here is original; it is similar to that of Reasenberg (1985), although derived independently. The significant difference is that whereas Reasenberg's technique always considers the first event of a sequence to be the main shock, and a subsequent, larger earthquake becomes a "larger main shock" (Reasenberg and Jones 1989), the method in this study always considers the largest event in a sequence to be the main shock (if two equal events occur, the first is the main shock), enabling one to discriminate separately between foreshocks and aftershocks, which may be useful in future studies.

The algorithm works in the following manner: at the start, all earthquakes in the catalogue are flagged as "unassigned". A space/time "window" of a certain number of days in extent and a certain number of kilometres in radius is defined. The algorithm then considers each earthquake in descending order of magnitude. If it is unassigned, it must be a main shock (it has not been associated with a larger event). This earthquake is therefore flagged as "main shock", and the catalogue is then worked through backwards in time, looking at each event. If an event falls within the space/time window centred on the mainshock, it is flagged as a foreshock and the start of the time window reset to the time of the newly identified foreshock. Once a period is found, equal in duration to the length of the time window, in which no foreshock can be identified, it is concluded that all foreshocks have been found. The same procedure is then used to look for aftershocks, starting with the time of the main shock and working forwards in time through the catalogue. Through this sieving technique, all events in the catalogue are eventually identified as main shock, foreshock or aftershock. The computation is implemented by a computer program AFTERAN, which has options for either writing out a "pruned" catalogue of main shocks or conducting analyses of the aftershock behaviour in terms the number of aftershocks observed for main shocks of differing magnitudes.

The size of the distance window used by the program was magnitude-dependent (see Table 1). A single value of 100 days was used for the time extent of the window. This value was arrived at by inspection of several clusters of seismicity in the catalogue in order to establish a useful value by direct observation. (Although it has to be admitted that there is still some subjectivity in this choice). The declustered catalogue is shown in Figure 1.
 
 
 
Magnitude Radius
4.0-4.4 30
4.5-4.9 35
5.0-5.4 40
5.5-5.9 47
6.0-6.4 54
6.5-6.9 61
7.0-7.4 70
7.5-7.9 81
8.0-8.4 94

Table 1 - Distance in km from main shock at which minor events are considered to be foreshocks or aftershocks, as a function of magnitude. (After Gardner and Knopoff 1974).
 
 


The seismic source model











The seismic source model produced for this study contains 50 source zones. The principles on which it is constructed are the subdivision of the study area into its component tectonic features, while reflecting the distribution of seismicity. The basis of the source zonation was the work done by other members of the Copernicus project CIPA-CT94-0238, as described in Mândrescu and Radulian (1996), Suhadolc (1996), Suhadolc and Panza (1996), Zivcic et al (1996), Mândrescu et al (1997), Zivcic and Poljak (1997). In addition, notice was taken of the GSHAP zonation of adjoining areas as documented in Grünthal et al (1996).

Some effort was required in order to arrive at a single, usable zone model for the whole region. While it was obviously desirable to keep as much as possible to the ideas expressed in the references above, some alterations were needed for the following reasons:

  1. To reconcile interpretations where differences existed;
  2. To ensure that the zonation adequately reflected the seismicity pattern;
  3. To eliminate zones too small to be analysed;
  4. To obtain a seamless coverage over the whole study area;
  5. To make such simplifications as were necessary for the production of a hazard map.
The zonation is shown in Figure 2, and the seismicity parameters of each zone are summarised in Table 2. Blanks in Table 2 indicate insufficient data. It should be born in mind when examining the depth distributions in Table 2 that some of them are based on very few earthquakes, and also some mean figures may be biased by values of 33 km which really indicate shallow indeterminate depths. The Gutenberg-Richter parameters a and b were calculated by the maximum likelihood method, using only portions of the catalogue considered complete following the analysis in Musson (1996). Parts of the frequency-distribution curve suffering from roll-off at low (or high) magnitudes were not used (they were excluded by eye). Uncertainties in a and b values were also calculated (but not listed here). In cases where the number of events is very low, clearly the values obtained are very approximate, but in such low seismicity cases the precise values are not very influential. The maximum magnitudes shown in Table 2 are the observed values; the values eventually used in the hazard calculations were the largest magnitude in the zone or any zone deemed to be potentially similar, plus a small safety margin.
No CODE Area Start year # events a b Mmax Events/sq km Mean depth Max depth
1
MURZ
10127
1267
50
3.223
-0.909
5.5
0.3815
12
43
2
WPER
6324
1443
13
1.288
-0.619
6.4
0.1026
12
36
3
EPER
29415
1453
20
3.144
-1.033
6.6
0.0349
10
40
4
STEI
27856
1590
17
3.575
-1.096
6.6
0.0557
12
20
5
WEAL
11552
1348
16
2.360
-0.900
5.7
0.0498
11
19
6
BELL
2863
1114
36
2.955
-0.906
6.8
0.7485
15
38
7
FRIU
3268
1278
65
2.717
-0.792
6.1
1.0832
13
35
8
SALP
9594
1077
76
4.249
-1.155
6.3
0.4436
11
35
9
PANN
90000
1092
21
2.706
-0.924
5.6
0.0114
12
19
10
VESZ
1428
1100
12
2.316
-0.900
5.6
0.3641
33
81
11
DUNA
5433
1258
17
2.958
-0.929
6.2
0.3213
11
17
12
EGER
3001
1767
9
2.933
-0.973
5.3
0.3662
9
16
13
EEAL
20015
1443
20
3.244
-1.022
6.5
0.0716
10
19
14
EPDM
8153
1781
19
2.830
-0.950
6.5
0.1314
16
50
15
DRAV
5653
1839
13
3.289
-1.015
5.6
0.2997
13
26
16
BACK
15197
1528
10
3.059
-0.969
5.7
0.1003
24
30
17
BANA
21501
1797
29
3.349
-0.966
5.6
0.1421
12
45
18
ZAGB
4862
1459
59
3.918
-1.045
6.2
1.1251
11
30
19
SAVA
27445
1502
34
3.434
-0.971
6.0
0.1293
13
33
20
NWED
7743
1279
37
2.172
-0.726
5.7
0.2394
12
30
21
RIJE
5641
1505
21
3.085
-0.906
5.8
0.5124
18
33
22
VELE
8620
1280
23
3.048
-0.973
5.7
0.1661
14
30
23
HERV
8259
1853
66
3.243
-0.848
6.1
0.8592
20
63
24
DALM
8300
1480
43
2.990
-0.848
6.1
0.4774
16
56
25
TREB
5662
1866
42
3.239
-0.906
6.1
0.7278
16
57
26
DRIN
12770
1563
63
3.597
-0.956
7.0
0.4643
15
46
27
IDIN
56817
1386
102
4.423
-1.091
6.0
0.2016
18
136
28
ADRI
83691
1123
114
3.738
-0.949
7.2
0.1045
24
72
29
EITA
17774
1268
99
4.036
-1.070
5.7
0.3208
17
48
30
APPN
27558
1243
277
3.888
-0.928
6.1
0.5442
14
59
31
APUS
32485
1614
6
2.735
-1.022
5.0
0.0137
20
33
32
TRAD
27889
1517
10
1.200
-0.680
5.6
0.0108
33
ECAR
35449
1599
12
2.535
-0.876
5.0
0.0303
71
230
34
SECM
11059
1473
15
1.995
-0.777
6.4
0.0697
15
23
35
VRAN
1884
1022
143
4.375
-0.902
7.5
31.0398
127
165
36
TRVR
6897
1906
24
2.198
-0.670
6.4
0.4779
64
150
37
MOEP
66997
1906
1
5.1
38
IMOF
11369
1901
13
1.796
-0.634
7.2
0.1601
29
60
39
SFGF
4647
1901
8
2.930
-0.963
5.2
0.2575
24
46
40
SWCM
15599
1665
15
3.050
-0.981
5.6
0.0857
15
20
41
IBAR
17475
1739
84
3.264
-0.845
6.0
0.4381
14
43
42
BALT
17868
1902
10
1.824
-0.693
4.8
0.0631
11
22
43
SKOP
14091
1641
39
3.402
-0.927
5.8
0.3508
15
33
44
SOFI
7316
1818
20
2.464
-0.794
6.6
0.2653
26
77
45
SSPM
18935
1909
3
2.565
-0.950
4.9
0.0307
21
40
46
MARI
4000
1641
16
3.329
-0.975
6.9
0.6713
15
37
47
PLOV
2084
1750
18
2.506
-0.788
7.0
1.0842
17
49
48
ESPM
7502
1875
7
3.269
-0.958
7.0
0.3646
15
33
49
LUKA
12193
1892
8
2.216
-0.758
5.2
0.1253
16
51
50
SZAG
7970
1907
8
1.379
-0.550
5.9
0.1895
26
40

Table 2 - Summary statistics for the 50 seismic source zones used in the model.

For the purposes of discussion, the zone model can be divided into six main areas, each containing several source zones. These areas are referred to, for convenience, as Alpine, Pannonian, Dinarides, Adriatic, Carpathian and Balkan. The division is not intended to be definitive; some zones are transitional and could be easily placed in more than one grouping.

Alpine Group

This area marks the northern boundaries of the Pannonian Fragment and the Adria Plate against the Eurasian Plate, summarised by Gutdeutsch and Aric (1987). The Mur-Muerz Valley, running SW-NE through Austria, is an important locus of seismicity (zone MURZ). This lineation continues eastwards as the Peripieninian Lineament, which is here divided into two zones on the basis of seismicity (WPER and EPER).

To the south, the active centres of seismicity around Friuli and Belluno are well known and important to the hazard in this area. Each is given its own zone (FRIU and BELL). To the east of Friuli is the Southern Alps zone (SALP), a zone of moderate and shallow seismicity, the tectonics of which are described by Zivcic and Poljak (1997). Between these three zones and the Mur-Muerz Valley is inserted a zone of relatively low seismicity (WEAL).

The seismically active area in Austria north of the Mur-Muerz Valley is treated as a single zone (STEI). The seismicity here is lower than in the adjacent zones described above, although it has a larger magnitude earthquake within it (unless this event is mislocated). The area further to the north, on the Czech border, is of sufficiently low seismicity to be ignored altogether, especially as this area is marginal to this study.

Pannonian Group

This group of zones makes up the Pannonian Basin itself, which can be divided into a number of different features.

The largest zone (PANN) occupies the Great Hungarian Depression. The seismicity here is low to moderate, but clusters of higher activity can be discerned, especially in the area of Veszprem (VESZ), the bend of the Danube above Budapest (DUNA) and the Eger district (EGER). The westernmost extension of the Pannonian Fragment, between the Mur-Muerz Valley and the Peradriatic Lineament, includes the Mura Depression and is modelled as a single zone (EEAL).

The southern edge of this system mostly comprises basins with a WNW-ESE trend, becoming more NW-SE towards the west. These include the Drava, Backa, Sava, Slavonia and Srem Depressions. The first two of these each are given one zone (DRAV, BACK). The remainder are similar in seismicity, and grouped together. To the east, the Banat Depression includes significant faulting with a roughly N-S trend (BANA).

A zone of seismicity with an apparent NE-SW trend is found to the northeast of the Pannonian Basin. This is the East Pannonian Basin Margin (EPBM).

Finally, the Zagreb-Balaton Zone (ZAGB) is characterised by a set of regional faults in a NE-SW direction with sinistral horizontal displacement (Royden and Horvath 1988). This is a marginal zone, part Alpine, part Dinaric.

The Dinarides Group

The main division within this group is between the External and Internal Dinarides. The former set of zones (NWED, RIJE, VELE, HERV, DALM, TREB, DRIN) marks the collision zone of the Adria Plate on its NW margin. This is an area of considerable seismicity, which, however, varies with local variations in the style of faulting, the predominant trend over the whole region being a dense pattern of NW-SE trending faults with also some thrusting from NE to SW. The division into individual source zones is a simplified version of that proposed by Markusic et al (1996) as reported by Suhadolc and Panza (1996). The highest seismicity is found in the southernmost of these zones (DRIN), which marks an area of some complexity on the flank of the Subpelagonic Massif, characterised by significant NE-SW trending faults intersecting with the characteristic Dinaric regional trend.

The Internal Dinarides represent a palaeo-subduction zone between the Adria Plate and the Pannonian Fragment (or microplate). As such, it includes deeper seismicity than is found elsewhere in the Dinarides. For the most part, this system, which preserves the Dinaric NW-SE trend, is modelled as a single zone (IDIN).

The Adriatic Group

The area of the Adriatic Sea and the Italian Peninsula is of marginal concern to the seismic hazard of the Pannonian area, and has therefore been treated in a very simplistic way. The Adriatic Sea has been treated as a single zone (ADRI), with Italy being divided into two zones, one occupying the eastern coast (EITA) and the remainder of this corner of the overall study area as a residual zone (APPN). Obviously, if one were concerned with sites in this area, a finer zonation would be required.

The Carpathian Group

The tectonics and seismicity of this area are discussed by Mândrescu et al (1997). The zones in this group cover the area around the Carpathian Mountains. In the north are three large zones of low seismicity: the Apuseni Mountains (APUS), the Transylvanian Depression (TRAD) and the Eastern Carpathians (ECAR). To the east of these lies the Moldavian Platform, which is of such low seismicity as to be almost aseismic, and is therefore excluded from the model. The same is true of the stable platform area in Ukranian territory to the north. The Transylvanian Depression zone is unusual in that, although it contains a few historical earthquakes, the last recorded event in the catalogue is in 1902.

South of these three low-seismicity zones are three zones of much higher activity. These are the eastern extremity of the South Carpathians (SECM), the well-known Vrancea source, a subduction feature producing intermediate focus earthquakes (VRAN), and a transitional zone around the main Vrancea source to the east and south (TRVR) wherein the seismicity is more diffuse. The Vrancea earthquakes are of key importance to the seismic hazard of the region. As can be seen in Table 2, when seismicity rates are adjusted for area, the activity in the VRAN zone, square kilometre for square kilometre, is fifteen times greater than the next highest zone.

South again lies the Moesian Platform, a stable area of very low seismicity (MOEP). It would be possible to drop this zone from the model in the same way that the Moldavian Platform was omitted; however, it is not so peripheral, and so is kept in. The seismicity is too low to allow sensible analysis, so semi-arbitrary values were assigned (these were actually imported from a neighbouring zone, SSPM, and reduced by not adjusting them for area).

The area from the Vrancean region to the Danube delta is more active than the territory to north and south, and also coincides with the Sfântul Gheorghe Fault. This is modelled as a separate zone (SFGF), as is the course of the deep Intramoesian Fault, which also appears to be an active feature (IMOF). The Intramoesian Fault continues southeastwards into the NE corner of Bulgaria (the Shabla region) where large earthquakes (7 Ms) have been observed. Modelling the entire fault as a single zone implies that such large earthquakes could also occur further north. This is conjectural and open to dispute; the significance of these two faults and the connection between the IMOF and the Shabla earthquake (if any) is poorly understood, and was, for example, a subject of much discussion at the 1st International Workshop on Vrancea Earthquakes, Bucharest, 1-4 November 1997. The area between the Sfântul Gheorghe Fault and the Intramoesian Fault is very inactive and has not been included in the model.

The Balkan Group

The zones of this group largely represent the Balkan Terrane, the western outer zone of the Protomoesian microcontinent, and parts of the Thracian Massif (Haydoutov and Yanev 1997). This area is characterised by a structural trend that continues as N-S from the Banat Depression, and swings around the Moesian Platform, becoming NW-SE and finally almost E-W in central Bulgaria.

This group includes the western part of the South Carpathian Mountains (SWCM) on the NW flank of the Moesian Platform. To the SW of this zone is a complex area of relatively high seismicity, including part of the South European Variscan Suture (IBAR). To the east and southeast runs the Starayna Planina Meganticlinorum, south of the Moesian Platform and the Cis-Balcanic folded zone. The pattern of seismicity has been modelled as two relatively low-activity zones (BALT, SSPM) and several local concentrations of higher activity (SKOP, SOFI, MARI, PLOV, ESPM and SZAG) of which three could be considered as continuing further south beyond the extent of the area covered by the earthquake catalogue. Finally, the Luda Kamcija Synclinorum appears to be a significant feature and occupies a single zone (LUKA).
 
 


Attenuation











The remaining topic to be discussed as respects input is the attenuation relationship to be used. Since a well-known study for Europe (or more specifically, Southern Europe) exists in Ambraseys and Bommer (1991), which includes strong motion data from Italy, Greece and former Yugoslavia, it was considered to be adequate to the purposes of this study to adopt this formula, which is

log10 A = -0.87 + 0.217 Ms - log10 R - 0.00117 R (1)

where A is peak horizontal ground acceleration in g, and R is hypocentral distance in km.

A more recent study also exists (Ambraseys 1995) which draws upon an expanded database. The reasons for not using this study are (a) it does not operate on hypocentral distance, but fault rupture distance with a fixed depth factor; (b) there is little difference between the results of the two studies in practice.

While equation (1) is adequate for most of the study area, it reflects attenuation from crustal earthquakes and is not appropriate for intermediate-focus events from the Vrancea seismic zone. Furthermore, it is well known that the release of seismic energy from Vrancea events is markedly directional (e.g. Mândrescu et al 1988) in nature. Thus, the use of a single isotropic attenuation equation such as that in equation (1) cannot give realistic results.

The attenuation of strong ground motion from Vrancea intermediate-depth earthquakes is studied by Lungu et al (1997) from existing accelerogram data. Such data exist from three earthquakes only (1977, 1986, 1990) and for the first of these, only one data point (Bucharest) exists. Lungu et al (1997) present three directional analyses - for the direction of Moldova (NE), Chernavoda (SE) and Bucharest (S) - and one general case using all data irrespective of azimuth. The data point for Bucharest from the 1977 earthquake is added to the data sets for the Moldova and Chernavoda directions as otherwise there would be only two earthquakes and a regression would not be possible (Lungu pers. com.).

The following equations are given:

ln A = 5.571 + 0.937 Ms -1.256 ln R - 0.0069 h (2)

ln A = 6.470 + 0.923 Ms - 1.403 ln R - 0.007 h (3)

ln A = 8.136 + 0.876 Ms - 1.657 ln R - 0.0076 h (4)

ln A = 4.150 + 0.913 Ms - 0.962 ln R - 0.006 h (5)

where A is peak horizontal ground acceleration in cm/sec2. These four equations relate to the complete data set, the Bucharest direction, the Cernavoda direction and the Moldova direction respectively.

There are two problems with this suite of equations. In the first case, there is no term for anelastic attenuation. Where one would expect to find this, there is instead a term for depth. This means that differences in the azimuthal attenuation that relate to anelastic attenuation cannot be modelled properly. The second problem is that as the equations were all regressed separately, (and perhaps partly because there was no anelastic term in the model used) the equations are not compatible at short distances. It is expected that attenuation is slowest in the NE direction and much more rapid in the SE direction, and this is borne out by the coefficients for ln R (i.e. the geometric spreading term) in equations (4) and (5). But since the constant in equation (5) is much higher than that for equation (4), the net effect is that, for earthquakes at the critical depth of 90 km, for epicentral distances out to about 130 km, calculated ground motions are higher to the SE of the Vrancea zone than they are to the NE. The hazard map (not shown) that results from using these equations has, therefore, two ellipticalities: the contours of highest hazard are oriented with the long axis towards the SE, while the contours of lowest hazard have their long axis aligned to the NE. This is quite contrary to what is seen in isoseismal maps (Mândrescu et al 1988), and would not be acceptable as a realistic map of hazard in the region.

There was not the possibility, within the scope of this project, to attempt to recalculate the regressions from the original data, which were not available. For this project a difficulty therefore arose: none of the published attenuation curves appeared to give reasonable results, and there were no data available against which to construct or test a new set of curves. In these difficult circumstances a stopgap solution was necessary, and after some experimentation the following course of action was adopted:

  1. Equations (2)-(5) were converted to give results in g rather than cm/sec2.
  2. The magnitude terms and geometric spreading terms were left untouched.
  3. The depth term was removed from equations (3)-(5) and an anelastic attenuation term (on R) introduced in its place. This was set to an arbitrary small value for equations (3) and (4), and zero for equation (5).
  4. The constant term in equations (3)-(5) was adjusted until each of them produced calculated accelerations at short epicentral distances that corresponded with those derived from equation (2).
In effect, the attenuation curve from the complete data set was used as an anchor for the other three, ensuring that the absolute values are sufficiently true to the original data set. The retention of the second two terms preserves the relative azimuthal variations in attenuation, and the addition of an arbitrary small anelastic term is a gesture towards what is probably the case, that anelastic attenuation is less in the north-easterly direction. The new coefficients are as follows:

ln A = -1.15 + 0.923 Ms - 1.403 ln R - 0.0004 R (6)

ln A = 0.33 + 0.876 Ms - 1.657 ln R - 0.0004 R (7)

ln A = -3.1 + 0.913 Ms - 0.962 ln R (8)

for the directions of Bucharest, Cernavoda and Moldova respectively, where A is now given in g.

It was assumed that attenuation to the NW is similar to that to the SE, i.e. equation (7) was used in this direction. This seems to borne out by isoseismal maps of the major Vrancea events. For intermediary directions, the software used interpolated linearly between equations (6)-(8). It is to be emphasised that these equations are only intended as stopgaps, incorporating necessarily arbitrary decisions in order to obtain realistic results in the frame of the present project; future work should involve a detailed study of the original data to provide definitive equations.
 
 


Results











Hazard calculations were carried out for the area 14o-28o E, 43o-49o N. This is deliberately slightly smaller than the area covered by the catalogue, in order to leave a border zone free, in which the hazard might be affected by seismicity occurring outside the catalogue area. Results were calculated at a grid interval of 0.25 degrees of longitude (about 20 km) and 0.20 degrees of latitude (about 22 km).

The final hazard maps are shown in Figures 3, 4, 5 & 6.

The hazard maps show that several areas within the region have enhanced hazard. The first of these, obviously enough, is the Vrancea region, which generates the highest hazard values on each map. The highest contour values for return periods of 100, 475, 1000 and 3000 years are 0.25g, 0.4g, 0.5g and 0.65g respectively. The shape of the contours for the Vrancea hazard are strongly influenced by the attenuation, extending strongly to the NE, but the simple shape provided by the attenuation variation is distorted to the S and W by the presence of adjoining source zones which contribute some additional hazard.

The second highest hazard in the region, perhaps surprisingly, is found to the west of Zagreb. Unlike Vrancea, this is not an area commonly associated with large earthquakes. The largest historical earthquake in the ZAGB source zone is the 9 November 1880 earthquake with magnitude 6.2 Ms. However, there are a relatively large number of events occurring within this zone compared to its area. This can be seen in Table 2, where the area-adjusted seismicity rate is actually greater than the zone for Friuli. In addition to this, the seismicity is shallow, with over 60% of the earthquakes which have depth information being at less than 10 km focal depth. The hazard comes, therefore, from earthquakes around magnitude 5 Ms being frequent over a small area, and possibly generating high accelerations through the natural scatter of acceleration values as determined by the parameter ?. For a return period of 475 years, the hazard value reaches 0.35 g.

The hazard is moderately high over the Dinarides, especially along the Dalmatian coast around the town of Split, but less in the area of the Velebit Planina. The hazard drops away from the coast, but rises again in a broad area between Belgrade and Nis, with hazard values of over 0.25 g at the 475 year level.

A strong SW-NE lineation of higher hazard runs through eastern Austria along the Mur-M?rz Valley and into Slovakia. Here, the highest contour values for return periods of 100, 475, 1000 and 30000 years are 0.1g, 0.2g, 0.25g and 0.35g respectively. To the south of this belt, a patch in north-central Hungary, NE of Lake Balaton, shows similar or higher hazard levels. A similar degree of hazard is also shown in part of Northern Bulgaria.

The places with the lowest hazard are the central part of the Hungarian-Romanian border (near Debrecen), the Czech-Austrian border and the extreme western part of the Ukraine.

The shape of the hazard contours is influenced by decisions made in the modelling process, and it is as well to be aware of where these may be subject to judgement, such that alternative interpretations might have a significant effect on the hazard. Such decisions have to be made with respect to the purpose of the hazard study, and a generalised hazard mapping study, as here, requires a different approach from a site hazard mapping study as might be required for providing engineering parameters for a particular structure. This is discussed to good effect by Page and Basham (1985). In order to show hazard at a broad, regional scale, the assumption is generally made that large earthquakes such as the 14 June 1913 earthquake (M=7) in northern Bulgaria are not restricted to the exact locations in which they occurred historically, but that other related features exist within the same geological structures which could produce similar earthquakes. The 1913 earthquake is therefore placed within the Eastern Starayna Planina Meganticlinorum zone, which is sufficiently broadly defined as to encompass some 7500 sq km. If this zone were drawn more tightly around the seismicity within this zone (which also includes events in 1875, 1914 and 1986) the effect would be to increase the hazard by concentrating it within a more constricted area. For a site-specific hazard study in this area, the effect of such an interpretation would have to be considered in the interests of conservatism. In a study such as the present one, which serves to present a regional outlook on hazard, there is not such a requirement to apply conservative options.

There are a number of other areas where this applies; perhaps the most important is the Shabla region in NE Bulgaria, where the large 31 March 1901 earthquake occurred. As discussed previously, in this study it is assumed that this earthquake is not fated to recur in one place only. An alternative interpretation would produce a peak of high hazard in the Shabla area. This possibility would need to be considered in the case of applications where the hazard in this area particularly was important. A less important case, but typical, is presented by Banja Luka, where events of magnitude 5.3 occurred in 1969 and 1981, with smaller events in 1897, 1935, 1950, 1970 and 1977. Modelling this activity as a small active zone would produce a peak of hazard at this location, the significance of which would be debatable. This argument can also be directed against the peak of hazard on the Slovenian-Croatian border, which appears counter-intuitive; however, it does seem that the seismicity here matches the tectonics and that the zonation is justified.

Finally, considering that hazard in the eastern part of the region is dominated by the Vrancea seismic source, there is still some uncertainty about the directional attenuation from intermediate focus earthquakes. There are three points to be made in this respect. Firstly, as discussed above, the attenuation relations used in this study are partly artificially produced, and it is expected these could be improved with further work and more data. Secondly, not much is known about attenuation to the W, SW and NW of Vrancea. Judging from isoseismal maps, it is likely that attenuation to the SW is less than simulated in this study, but more data is necessary before this can be quantified. Thirdly, some Vrancea earthquakes are known to have isoseismal elongation perpendicular to the usual axis, and also fault plane solutions different from the norm. These different events have not been factored into the hazard calculations. Probably they result from tearing within the subducting slab; they are therefore likely to be consistently smaller than "normal" Vrancea earthquakes (limited rupture width) and so not taking account of them in the hazard calculations should not have too great an effect.
 
 


Conclusions











This study presents four maps of seismic hazard at different return periods for a wide area including all of Hungary, Slovenia and Croatia, most of Romania and Serbia, and parts of Bulgaria, Slovakia, Austria and other adjoining countries. The hazard in the eastern part of the region is dominated by a single source of seismicity: the intermediate focus Vrancea activity. The western part is more complex, with the hazard being greatest along the Dalmatian coast and the Croatian-Slovenian border.

These maps are not intended for use in generating values for specific sites for engineering purposes, and are also not intended to supplant national maps of seismic hazard in any of the countries covered. National maps may implement specific policies for local purposes which are not considered here.

These maps can be used as a general guideline to regional distributions of seismic hazard in terms of peak ground acceleration values over the general area. They can also be used as a basis for comparative studies of seismic hazard in the area using different methodologies.
 
 


Acknowledgements











I would like to thank my various colleagues on the Copernicus Project "Quantitative Seismic Zoning of the Circum Pannonian Basin" (EC Project CIPA-CT94-0238) for helpful discussions at various stages of this work; especially to Mircea Radulian, Mladen Zivcic, Peter Suhadolc and Giuliano Panza. I would also like to thank Paul Henni for his assistance with preparing the maps in this report, and also Dan Lungu, Gottfried Grünthal, Max Stucchi for assistance at various stages.This work was supported by the European Community and the Natural Environment Research Council and is published with the permission of the Director of the British Geological Survey (NERC).
 
 


References











Ahorner, L. and Rosenhauer, W., Seismiche risikoanalyse, In Naturkatastrophen und katastrophenvorbeugung, Bericht zur Int. Decade Natural Disaster Reduction, (Deutsche Forschungsgemeinschaft, Verlag, Weinheim 1993), 177-190.

Ambraseys, N.N., and Bommer, J.J., (1991) The attenuation of ground accelerations in Europe, Eq. Eng. and Structural Dyn., 20, 1179-1202.

Ambraseys, N.N., (1995) The prediction of earthquake peak ground acceleration in Europe, Eq. Eng. and Structural Dyn., 24, 467-490

Bus, Z., Szeidovitz, G., and Vaccari, F., (1999), Synthetic seismogram based deterministic seismic zoning of the Hungarian part of the Pannonian Basin, Pageoph., in press.

Coppersmith, K.J. and Youngs, R.R., Capturing uncertainty in probabilistic seismic hazard assessments within intraplate tectonic environments, In Proc. Third US Nat. Conf. Eq. Eng, (Charleston, 1986) 1, 301-312.

Gardner, J.K. and Knopoff, L., (1974) Is the sequence of earthquakes in Southern California, with aftershocks removed, Poissonian?, Bull. Seism. Soc. Am., 64, 1363-1367.

Grünthal, G., Bosse, C., Musson, R.M.W., Gariel, J.-C., de Crook, T., Verbeiren, R., Camelbeeck, T., Mayer-Rosa, D. and Lenhardt, W., Joint seismic hazard assessment for the central and western part of GSHAP Region 3 (Central and Northwest Europe), In Seismology in Europe, (ed. Thorkelsson, B.) (Icelandic Met. Office, Reykjavik, 1996) 339-342.

Gutdeutsch, R. and Aric, K., Tectonic block models based on the seismicity in the East Alpine-Carpathian and Pannonian area, In Geodynamics of the Eastern Alps, ( ed. Flügel, H.W. and Faupl, P.) (Deuticke, Vienna, 1987).

Haydoutov, I. And Yanev, S., (1997), The Protomoesian microcontinent of the Balkan Peninsula - a peri-Gondwanaland piece, Tectonophysics, 272, 303-313.

Herak, M., Herak, D. and Markusic, S., (1996), Revision of the earthquake catalogue and seismicity of Croatia 1908-1992, Terra Nova, 8, 86-94.

Johnson, C.E. and Koyanagi, R.Y., (1988), A Monte-Carlo approach applied to the estimation of seismic hazard for the state of Hawaii, Seism. Res. Letters, 59 no 1, 18.

Lungu, D., Cornea, T., Aldea, A. and Zaicenco, A., Basic representation of seismic action in Design of structures in seismic zones (ed. Lungu, D., Mazzolani, F. and Savidis, S.) (Bridgeman Ltd, Timisoara, 1997) 9-60.

Mândrescu, N., Anghel, M. and Smalbergher, V., The Vrancea intermediate-depth earthquakes and the peculiarities of the seismic intensity distribution over the Romanian territory, In Recent seismological investigations in Europe (ed. Neresov, I.L. et al), Proc. XIX ESC General Assembly, (Moscow, 1988), 59-65.

Mândrescu, N. and Radulian, M., (1996), Characterisation of seismogenic zones of Romania, EEC Technical Report, Project CIPA-CT94-0238.

Mândrescu, N., Popescu, E., Radulian, M., Utale, A. and Panza, G., (1997) Seismicity and stress field characteristics for the seismogenic zones of Romania, EEC Technical Report, Project CIPA-CT94-0238.

Markusic, S., Costa, G., Vaccari, F., Suhadolc, P. and Herak, M., (1996) Deterministic seismic zoning of the Croatian territory derived from complete synthetic seismograms, in preparation (referred to by Suhadolc and Panza 1996, vide infra).

Markusic, S., Suhadolc, P., Herak, M., and Vaccari, F., (1999) Seismic hazard elements in Croatia based on deterministic modelling, Pageoph., in press.

Musson, R.M.W., (1994), A catalogue of British earthquakes, BGS Technical Report No WL/94/04.

Musson, R.M.W., An earthquake catalogue for the Circum-Pannonian Basin, In Seismicity of the Carpatho-Balcan Region, Proc. XV Congress of the Carpatho-Balcan Geol. Ass., (ed. Papanikolaou, D. and Papoulia, J.,) (Athens, 1996), 233-238.

Musson R.M.W. and Winter P.W., (1997), Seismic hazard maps for the UK, Natural Hazards, 14, 141-154.

Page, R.A. and Basham, P.W., (1985), Earthquake hazards in the offshore environment, US Geological Survey Bulletin no 1630.

Radulian, M., Vaccari, F., Mandrescu, N., Panza, G.F. and Moldoveanu, C.L. (1999), Seismic hazard of Romania: deterministic approach, Pageoph., in press.

Reasenberg, P., (1985), Second-order moment of Central Californian seismicity, 1969-1982, J. Geophys. Res., 90, 5479-5495.

Reasenberg, P. and Jones, L.M., (1989), Earthquake hazard after a mainshock in California, Science, 243, 1173-1176.

Reiter, L., Earthquake hazard analysis, (Columbia UP, New York, 1990).

Rosenhauer, W., Methodological aspects encountered in the Lower Rhine area seismic hazard analysis, In, Seismicity and seismic risk in the offshore North Sea area, (ed. Ritsema, A.R. and Gürpinar) (Reidel, Dordrecht, 1983) 385-396.

Royden, L.H. and Horváth, F. (editors) (1988), The Pannonian Basin. Am. Assoc. Pet. Geol., Mem., 45, 27-48.

Shapira, A., (1983), Potential earthquake risk estimations by application of a simulation process, Tectonophysics, 95, 75-89.

Suhadolc, P., (1996), Structural models in the Eastern Alps and Bulgaria, EEC Technical Report, Project CIPA-CT94-0238.

Suhadolc, P. and Panza, G., (1996), Focal mechanisms and seismogenetic zones, EEC Technical Report, Project CIPA-CT94-0238.

Zivcic, M. Cecic, I. and Poljak, M., (1996), Seismicity and geodynamics of Slovenia, EEC Technical Report, Project CIPA-CT94-0238.

Zivcic, M. and Poljak, M., (1997), Seismogenetic areas of Slovenia, EEC Technical Report, Project CIPA-CT94-0238.

Zivcic, M., Suhadolc, P. and Vaccari, F., (1999), Seismic zoning of Slovenia based on deterministic hazard computations, Pageoph, in press.
 
 


Figure captions









Figure 1   Seismicity of the study area. Earthquakes identified as foreshocks or aftershocks are not shown.

Figure 2   The seismic source zone model used in this study. See Table 2 for details of each zone.

Figure 3   Contour map of horizontal peak ground acceleration values with return period of 100 years.

Figure 4   Contour map of horizontal peak ground acceleration values with return period of 475 years.

Figure 5   Contour map of horizontal peak ground acceleration values with return period of 1,000 years.

Figure 6   Contour map of horizontal peak ground acceleration values with return period of 3,000 years.