# Statistical analysis of time-dependent earthquake occurrence and its impact on hazard in the low seismicity region Lower Rhine Embayment

L. Faenza, S. Hainzl, F. Scherbaum, C. Beauval

<figure><img src="image_1.png" /></figure>

## To cite this version:

L. Faenza, S. Hainzl, F. Scherbaum, C. Beauval. Statistical analysis of time-dependent earthquake occurrence and its impact on hazard in the low seismicity region Lower Rhine Embayment. Geophysical Journal International, 2007, 171 (2), pp.797-806. 10.1111/j.1365-246X.2007.03564.x. hal-00407725

HAL Id: hal-00407725

https://hal.science/hal-00407725v1

Submitted on 2 Oct 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Statistical analysis of time-dependent earthquake occurrence and its impact on hazard in the low seismicity region Lower Rhine Embayment

Licia Faenza,<fnref n="1" /> Sebastian Hainzl,<fnref n="1" />,<fnref n="2" /> Frank Scherbaum<fnref n="1" /> and Céline Beauval<fnref n="3" />

Accepted 2007 July 27. Received 2007 July 26; in original form 2007 May 9

## SUMMARY

The time-dependence of earthquake occurrence is mostly ignored in standard seismic hazard assessment even though earthquake clustering is well known. In this work, we attempt to quantify the impact of more realistic dynamics on the seismic hazard estimations. We include the time and space dependences between earthquakes into the hazard analysis via Monte Carlo simulations. Our target region is the Lower Rhine Embayment, a low seismicity area in Germany. Including aftershock sequences by using the epidemic type aftershock-sequence (ETAS) model, we find that on average the hypothesis of uncorrelated random earthquake activity underestimates the hazard by 5–10 per cent. Furthermore, we show that aftershock activity of past large earthquakes can locally increase the hazard even centuries later. We also analyse the impact of the so-called long-term behaviour, assuming a quasi-periodic occurrence of main events on a major fault in that region. We found that a significant impact on hazard is only expected for the special case of a very regular recurrence of the main shocks.

**Key words:** long- and short-time earthquake behaviour, low seismicity regions, Monte Carlo simulations, probabilistic hazard analysis, statistical analysis.

## 1 INTRODUCTION

Seismic risk assessment of urbanized areas is a major challenge to societies exposed to earthquakes. In a general framework, the seismic risk is usually defined as the product of three factors: the ‘hazard’, the ‘vulnerability’ and the ‘element at risk’. The ‘hazard’ is the probability of any particular area to be shaken by an earthquake during the exposure time; the ‘element at risk’ consists of the social and economical quantification of the objects exposed to risk. It can be expressed in terms of density of inhabitants, as capital value (land, buildings, etc.), or as productive capability (factories, power plants, highways, etc.). The ‘vulnerability’ quantifies the predisposition of an object to undergo any damage.

For the seismic hazard assessment, the quantification of the hazard can be related to a variety of parameters indicating the intensity of shaking. Therefore, the expected outcome of a hazard assessment is the evaluation of the probability of exceeding a selected set of ground motion intensity values during an assumed exposure time (SSHAC 1997).

Often the first step of a hazard analysis is the statistical modelling of the past earthquakes that have occurred in the region under study. This requires the evaluation of the seismicity rates, the identification of the seismo-tectonic sources, the quantification of ground motion relations and the study of earthquake spatio-temporal distribution

(Reiter 1990). In this paper, we will focus on the latter as a crucial aspect.

In the most common approach to hazard assessment, the temporal occurrence of earthquake is assumed to be Poissonian (Cornell 1968; McGuire 1976), hence imposing a random temporal behaviour of the seismicity. Recent studies on the other hand (e.g. Kagan & Jackson 2000; Kagan et al. 2003; Cinti et al. 2004) show clear deviations from a Poissonian behaviour in seismic catalogues. While the presence of cluster activity after a main shock is worldwide accepted (Ogata 1988), the so-called ‘long-term’ behaviour is still a question of debate, and there are different methodologies to approach it (an overview of the problem can be found on Working Group On California Earthquake Probabilities 2003).

Moreover, in low seismicity areas the scarce information on active faults and the poverty of the data available make in general the analysis and the discrimination between competitive models more complicated.

In a recent paper (Beauval et al. 2006), the impact of non-Poissonian occurrence on hazard has been studied for the seismic region of southwest France. While Beauval et al. (2006) studied the inclusion of aftershock sequence into seismic hazard, in this paper we further generalize this approach to include also the long-term behaviour of cyclic occurrence of main events and we provide a study of the role of real past historical events in nowadays hazard. The

<footnote n="1">¹Institute of Earth Sciences, University of Potsdam, Karl-Liebknecht Str. 24, 14476 Potsdam – Golm, Germany. E-mail: licia@geo.uni-potsdam.de</footnote>

<footnote n="2">²GeoForschungsZentrum, Telegrafenberg, 14473 Potsdam, Germany</footnote>

<footnote n="3">³IRD Géosciences Azur 250, rue Albert Einstein, Sophia Antipolis, 06560 Valbonne, France</footnote>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Target area: Lower Rhine Embayment, Germany. The two sites Cologne and Aachen. Seismicity data come from the Leydecker (2005) catalogue, complete for $M_L \ge 2.0$ since 1974, for this region. The box indicates the study area. The line represents the location of the synthetic fault used for the simulation. The diamonds indicate the location of the Düren earthquake of 1756 and the Roermond earthquake of 1992.</figcaption></figure>

new target region is the Lower Rhine Embayment (LRE) (Fig. 1), one of the regions of highest seismic risk in Germany (Hinzen & Oemisch 2001; Tyagunov et al. 2006). The short-term clustering is modelled through the epidemic type aftershocks-sequence model (ETAS, Ogata 1988) which has become a popular technique to model earthquake interaction in a statistical sense. Via a Monte Carlo technique, the ETAS model is applied on timescales often used in seismic hazard studies with 50 yr of exposure time (Swiss Earthquake Hazard, Giardini et al. 2004; The Seismic Hazard Map for Italy, MPS Working Group 2004). The 50 yr exposure time is used in earthquake engineering because it is supposed to represent the mean lifetime of conventional buildings. Using the same model, we also estimate the potential contribution of aftershocks from the Düren earthquake of 1756 February and the Roermond earthquake of 1992 April for the hazard in the following 50 yr. For the long-term behaviour, we perform a synthetic test on the impact of possibly cyclic activity on a single major fault. In this case, we quantify the impact for varying degrees of knowledge about the past fault activity.

## 2 THE METHODOLOGY

Following Beauval *et al.* (2006), the evaluation of the probability of non-exceedance of a particular ground motion value during the exposure time is evaluated using a Monte Carlo technique. The use of synthetic catalogues analysis is one way to calculate the probabilistic seismic hazard (Rosenblueth 1964; Musson 1999; Smith 2003; Giardini *et al.* 2004).

In the conventional approach, the hazard at a site is evaluated by considering the ground motion from all possible damaging earthquakes that can occur in a region. If a Poissonian model is assumed, it is easy to express the exceedance probability in term of the return periods (Cornell 1968; Ang & Tang 1975). In fact, in a time-independent study, the annual rate is sufficient to express earthquake occurrence. In this case, the seismic hazard can be computed in terms of annual rates of exceedance of selected acceleration levels ($A^*$) at the site of interest, using the relation (Ang & Tang 1975): $P(A \ge A^*) = 1 - \exp(-\lambda t)$, where $\lambda$ is the annual rate of the target event and $t$ is the time of interest. Using a Taylor series expansion for $\lambda t$ small, and truncated at the first term, the exceedance probability becomes approximated equal to the annual rate $\lambda$ times $t$. Considering a time-dependent earthquake occurrence, the seismic hazard cannot be expressed in term of annual rate of earthquake occurrence any longer. In Beauval *et al.* (2006), the Monte Carlo techniques are applied in the study of the temporal behaviour of earthquakes, in particular to estimate the impact of the Poisson hypothesis on seismic hazard, versus a time-dependent behaviour of seismicity, either in the short- or long-term. This approach is further utilized in the present study.

For the evaluation of the hazard, a sufficient large number of synthetic catalogues ($N$) of time duration $t$ is generated. The time duration $t$ is set to 50 yr, a commonly used value for the exposure time (e.g. The Seismic Hazard Map for Italy, MPS Working group 2004). Each catalogue can be considered as a representation of seismicity in the time period of interest. For each catalogue, a set of ground motions is evaluated. The probability of non-exceedance of a level $A^*$ of ground motion at a specific site in the time $t$ is computed by
counting the intervals in which A* did not occur

$$P(A^*; t) = \lim_{N \to \infty} \frac{1}{N} \sum_{i=1}^{N} H(A^* - A_{\max,i}), \quad (1)$$

where N is the number of catalogues of time duration t; H is the Heaviside function and Amax,i is the maximum ground motion value occurred at a site of interest during the ith catalogue of time duration t. The complement of P(A*; t) is the probability that A* is exceeded at least once in the time period t. An alternative way to estimate the probability of non-exceedance is to consider directly the empirical density function of the maximum acceleration values Amax,i. Each probability of non-exceedance is identical to the corresponding percentile of this distribution (further details in Beauval et al. 2006).

## 2 SHORT-TERM CLUSTERING

We focus on the LRE in Germany using the so-called Leydecker catalogue (2005), maintained by the BGR, Federal Institute for Geosciences and Natural Resources. We have chosen this catalogue because it is the most updated and with the lowest threshold magnitude for our target areas. We found that the catalogue is complete for this study area after 1974 for $M_L \ge 2.0$ by using the cumulative number of events as a function of time and the Gutenberg–Richter relation; our results are in agreement also with the finding of Schmedes et al. (2005). In Fig. 1, the distribution of the instrumental seismicity, the box defined for the analysis and the locations of the two target sites (Aachen & Cologne) are shown; 191 events have been analysed. We are conscious that this time interval is dominated by the 1992 Roermond earthquake, but this time window is the best data available we have at the moment for this regions. In the first step in our analysis, we study the statistical patterns that represent the recorded seismicity, applying three different techniques.

### 3.1 Coefficient of variation

The coefficient of variation (CV) is a non-parametric test of the measure of the dispersion. The random variables considered are the interevent times, which are the time intervals between two subsequent events. The CV is expressed as the ratio of the standard deviation to the mean. If the CV is less then 1, the distribution corresponds to a quasi periodic behaviour in the temporal domain; if the CV is equal to 1 the distribution is random in the temporal domain, since the standard deviation and the mean have the same value; and, finally, if the CV is larger than 1 the distribution has a cluster behaviour (Cox & Lewis 1966). Our results suggest that the events occur in clusters, since the CV is equal to 1.5.

### 3.2 Study of the hazard-rate function

The second technique we apply is a non-parametric estimation of the so-called hazard-rate function which represents the instantaneous conditional probability of occurrence of an event at time t upon survival (non-occurrence of an event) until time t. The hazard-rate function unambiguously describes the statistical distribution of the earthquake point process (e.g. Kalbfleisch 1985; Faenza 2005). In some aspects it is preferable to other statistical functions, such as the cumulative or the density function, since a simple analysis of its trend versus time can provide useful insights on the physics of the process.

The Tanner & Wong (1984) technique, which we apply in our study, approaches the problem of estimating the hazard-rate function directly by smoothing the empirical rate. The random variables considered are the interevent time and the censoring time. The consideration of the censoring time, which represents the time elapsed between the most recent event and the end of the catalogue, becomes very important in time-dependent analysis. The Kernel estimator of the hazard-rate function $\lambda(t)$ is

$$\lambda(t) = \sum_{i=1}^{n} \frac{\delta_i}{n - i + 1} K_{\Theta}(t - y_i), \quad (2)$$

where $y_i$ is the ith ordered random variable of the system, $\delta_i$ is an indicator associated to $y_i$, that is $\delta_i = 1$ in case of interevent data, $\delta_i = 0$ in case of censoring; $K_\Theta$ is the kernel function depending on a positive smoothing vector $\Theta$; and n is the total number of data points. The smoothing vector $\Theta$ is evaluated using a modified-likelihood criterion. The details of the algorithm can be found in Tanner & Wong (1984). Numerical checks have been done to verify the consistency of the results and to evaluate the performance of the non-parametric estimators in comparison to the one of the parametric estimators at different setting.

In Fig. 2, the shape of the hazard-rate function versus time for our data set (events with $M_L \ge 2.0$ since 1974) is shown. The decreasing trend of the function stands for a clustered behaviour of the seismicity (Kalbfleisch 1985; Faenza 2005) for a time duration of a few years. These results of a non-parametric analysis of the data motivated the application of the ETAS model in the following.

### 3.3 Epidemic type aftershock-sequence

The inclusion of the aftershock activity into the analysis is done through the ETAS model. The ETAS model is a stochastic marked point process for the representation of the occurrence of earthquakes of size larger than or equal to a threshold magnitude, in a region and in a period of time (Ogata 1988). The advantage of this triggering model is that its application on real data does not require the discrimination of events; it works by considering the seismicity as the superposition of seismicity induced by previous events on the background. The background events represent the tectonic loading and are modelled as a Poisson process with a constant rate $\mu$ while the cascade of aftershocks is described by the empirical Omori–Utsu law (Utsu et al. 1995), in which the number of events is proportional to $(t + c)^{-p}$, with t being the time from the occurrence of the main shock. The magnitude is randomly selected from a Gutenberg–Richter relation, for both the tectonic and triggered events; the productivity of the sequence is proportional to $K10^{\alpha M}$, with M being the magnitude and $\alpha$ and K two constants. We do not describe the characteristics of this parametric model in detail: its complete description can be found in Ogata (1988, 1998). For the joint inversion of the five parameters for our catalogue ($M_L \ge 2.0$ since 1974), the maximum-likelihood method (Ogata et al. 1993) is used; we obtain: $\mu = 1.35 \text{ yr}^{-1}$, $K = 0.0083$, $\alpha = 0.70$, $p = 0.98$ and $c = 0.5310^{-5} \text{ yr}$. The parameter b of the Gutenberg–Richter relation is taken from Schmedes et al. (2005) and it is equal to $0.96 \pm 0.03$.

Application of the Akaike's Information Criterion (AIC Akaike 1974) to compare the goodness-of-fit between the Poisson model and the ETAS one for this catalogue yields that the ETAS model gives a significantly better description of the instrumental seismicity ($AIC_{ETAS} - AIC_{Poisson} = -574$).
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Hazard-rate function as a function of time using the algorithm of Tanner & Wong (1984). The decreasing trend stands for a cluster behaviour of the seismicity in the temporal domain, for a time duration of a few years (Kalbfleisch 1985; Cinti *et al.* 2004; Faenza 2005).</figcaption></figure>

## 3.4 Generation of the synthetic catalogues

Following the previous analysis, we conclude that the ETAS model is a good way to represent the seismicity. We will refer to the catalogues generated by using the ETAS model as the ETAS-catalogues.

For the generation of these catalogues, an ‘inverse’ method is used according to Felzer *et al.* (2002). In each realization, events are simulated sequentially: first the time, then the magnitude and the epicentral coordinates. The maximum magnitude is set to $M_w$ 7.0 following Schmedes *et al.* (2005), in which is given the probability density function of the maximum events in the LRE. The conversion from $M_w$ and $M_L$ is taken from Grünthal & Wahlström (2003). The inclusion of earthquakes with long periods is, therefore, controlled by a Gutenberg–Richter relation calibrated on the seismicity of this area. For the spatial distribution of the tectonic events, a Gaussian filter with correlation distance equal to 20 km is applied (Frankel 1995). For the selection of the independent events, we used the declustering procedure proposed by Zhuang *et al.* (2002) and we find that 40 per cent of the events are identified as independent. Remarkably, by using the technique of Hainzl *et al.* (2006) to estimate the aftershock fraction, we get the same result. The filter is applied in a grid with cells of 10 × 10 km. For the location of the clustered events we follow an isotropic power-law distribution $\rho(r) = cr^{-n}$ with $n = 1.36$, in agreement with Helmstetter *et al.* (2003) and Felzer & Brodsky (2006). The depth of the events is set to 10 km for simplicity. The threshold magnitude for the cascade process is set equal to 2.0.

## 4 CONSIDERATION OF MAIN SHOCK CYCLES

One of the goals of this paper is to study the impact of a possible quasi-periodic recurrence of characteristic events on major faults. For this purpose, we perform a synthetic test for a generic fault which

is in its position, length and strike in agreement with observed surface signatures (Hinzen 2004), see Fig. 1. The magnitude of the characteristic earthquake on that fault is assumed to be $M_w$ 7.0 (see Schmedes *et al.* 2005) with a recurrence time distribution according to the Brownian Passage Time (BPT) distribution (Matthews *et al.* 2002). Palaeoseismologic studies indicate that similar large events occurred in the past (Camelbeech & Meghraoui 1998; Camelbeech *et al.* 2000). We choose the BPT distribution since it is the statistical representation of Reid’s elastic rebound theory, but we are aware that our choice is subjective. In this model, the occurrence of major earthquakes is controlled by a steady tectonic loading perturbed by a Brownian motion representing stochastic stress fluctuations, for example, due to earthquake interactions. Events happen in this model when the critical failure threshold is reached. From a statistical point of view, the distribution has two parameters: the mean and the aperiodicity. In this study, the mean recurrence time is deduced from the Gutenberg–Richter relation obtained in the work of Schmedes *et al.* (2005); it is 9500 ± 3900 yr for the assumed magnitude $M_w$ 7.0 event. The aperiodicity is the shape parameter of the distribution because its value changes the width of the distribution. Smaller values of the aperiodicity represent more regular temporal behaviour of the sequence, with a nearly symmetric density function where the central value is close to the mean; while larger values produce more random sequences, with a density function skewed to the right and sharply peaked at a value left of the mean (Matthews *et al.* 2002). One of the problems of a low seismicity area, as formerly mentioned, is the scarce knowledge of the activity on individual faults, therefore, it is not possible to estimate the aperiodicity for such regions, and we decide to work with three distinct values: 0.3, 0.5 and 0.7; in order to range different possibilities in cycling earthquake occurrence and to agree with the Working Group On California Earthquake Probabilities (2003).

We called these catalogues BPT+ETAS-catalogues; for their generation, we follow the same procedure as before for the ETAS-catalogues. However, now the characteristic on-fault earthquakes
<figure><img src="image_5.png" /><figcaption>**Figure 3.** The empirical density function of the $A_{max,i}$ in 50 yr based in 20 000 catalogues, for the ETAS-catalogues in the two sites of Cologne and Aachen. The solid line represents the 50 per cent percentile; the dashed line the 90 per cent percentile at Cologne: maximum acceleration value of 0.036 and 0.09 g, respectively. For Aachen, the corresponding values are 0.049 and 0.12 g.</figcaption></figure>

are added to the Poissonian background activity before aftershock sequences are calculated.

# 5 APPLICATION AND RESULTS

As target locations we choose the cities of Cologne and Aachen, two highly urbanized areas in Germany. For the calculation of the hazard, a sufficient large number of catalogues (20 000) of 50 yr duration was generated. For each $M_s \ge 4.0$ earthquake, the peak ground acceleration value at the site of interest was calculated by using the Berge-Thierry *et al.* (2003) ground motion relation. $M_s = 4.0$ is the threshold magnitude for this ground motion relation. The conversion from $M_L$ to $M_s$ follows Ambraseys (1990). The log (PGA)-value was chosen randomly from a Gaussian density function with a standard deviation of 0.2923 and truncation at three standard deviations.

In Fig. 3, the density function of maximum acceleration in 50 yr for the two sites is shown for the ETAS catalogues. The probability of non-exceedance can be deduced directly from this figure and it corresponds to the percentile of the distribution (Beauval *et al.* 2006). For instance, 90 per cent of probability of non-exceedance corresponds to 0.09 g at Cologne and 0.12 g at Aachen which is in general agreement with previous estimations (Grünthal & Wahlström 2006).

## 5.1 Impact of the short-term clustering

Under the Poissonian hypothesis, only independent events are considered for analysis. In our study, we have used two different approaches for declustering: in the first case, following the ETAS model philosophy, the main shock is defined as the first event in the cluster, independently of its magnitude value (case 1); in the second case, the main shock is chosen as the largest event in the cluster and

we will refer to this latter case as the ‘perfect declustering’ (case 2). The declustered catalogues are called Poissonian-catalogues.

Once the ETAS-catalogues and their corresponding Poissonian-catalogues have been generated, the impact of the time-independent hypothesis can be evaluated for the two sites (Fig. 4). Its quantification is done by computing for every percentile the difference between the PGA-values calculated for the Poissonian-catalogue and the ETAS-catalogue, normalized over the value for the ETAS-catalogue.

As shown in Fig. 4, the results are similar for the two locations. We refer to case 1 and case 2 as the lower and upper bound of the impact of the Poissonian hypothesis. Here we remark that the difference between case 1 and case 2 is in the declustering algorithm. Case 2 yields a systematic but not very high impact equal to 8 at 90 per cent probability of non-exceedance in 50 yr. This implies that perfect declustering yields a systematic underestimation of the hazard of 8 per cent which is in agreement with the result Beauval *et al.* (2006) for the seismicity in southwest France.

## 5.2 Impact of historic events

In a previous study, Ebel *et al.* (2000) showed that aftershock sequences of historic main events can locally dominate the present seismicity. To quantify the impact of ongoing aftershock activity of larger historic events, we selected two prominent earthquakes: the Düren earthquake of February, 1756 with $M_L = 6.4$ and the Roermond earthquake of 1992 April with $M_L = 5.9$ (Leydecker 2005), see Fig. 1. The quantification of their aftershock hazard was done by using Monte Carlo simulations. For this purpose, we generated random ETAS-type aftershock sequences triggered by a main shock with the same characteristics as the observed historic earthquakes in magnitude and location. We assume that the ETAS parameters for these sequences are the same we derived for the seismic
<figure><img src="image_6.png" /><figcaption>**Figure 4.** The impact of the Poisson hypothesis at the two locations. In case 1, the main shock is defined as the first event in the sequence; in case 2, the main shock is defined as the largest event in the cluster.</figcaption></figure>

catalogue. The hazard was then calculated on basis of 20 000 cata-
logues of 50 yr. Each catalogue begins at an elapsed time matching
the occurrence of the event, namely 250 yr in the first and 15 yr
in the second case. The contribution of the specific aftershock
sequence was calculated as the ratio of the hazard estimated from
the aftershock sequence alone to the one estimated from the ETAS-
catalogue.

The results are shown in Fig. 5. To consider the uncertainty linked
to the magnitude for the Düren event, we consider as lower and up-

per estimates the values M_L = 5.9 and M_L = 6.9. In Fig. 5 the
results are reported as dash and dotted lines, respectively. Remark-
ably, the contribution to the hazard driven by the 1756 event is still
large for high probabilities of non-exceedance, where it reaches the
maximum value of 20 per cent for the location of Aachen, while for
Cologne the contribution if significantly less. The difference for the
two sites results from the location of the event, which is close to
the city of Aachen; the location errors in Leydecker (2005) ranges
between ±5 km for this events, see Fig. 1. Furthermore, it is found

<figure><img src="image_7.png" /><figcaption>**Figure 5.** Top panel: the contribution of ongoing aftershock sequences triggered by the 1992 Roermond event, for the location of Aachen and Cologne. Lower panels: the contribution of ongoing aftershock sequences triggered by the 1756 Düren event. To consider the uncertainties linked to its magnitude, the values *M*<sub>L</sub> = 5.9 and 6.9 are considered as lower and upper estimations. They are shown in figure as dash and dotted lines.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Plot of the hazard-rate function for the BPT distribution for three values of the aperiodicity ($a = 0.3, 0.5$ and $0.7$). The vertical lines represent the cases we study: elapsed time equal to 10 per cent of the mean (black); elapsed time equal to the mean recurrence time (light grey) and elapsed time equal to 10 yr (dark grey).</figcaption></figure>

that the aftershocks of the Roermond earthquake, occurred in 1992, will contribute more than 10 per cent to the hazard at 90 per cent probability of non-exceedance in the next 50 yr.

## 5.3 Impact of the quasi-periodic on-fault main shocks

In a time-dependent perspective, as for example, a quasi-periodic recurrence of a characteristic event on a fault, the time elapsed since the most recent event plays a crucial role for the hazard function. For the BPT distribution, three different periods can be identified. The probability of the occurrence of the next characteristic event is almost zero immediately after the last event, representing, from a physical point of view, the unloading of the fault after the shock and its inability to generate another one. With increasing elapsing time, the conditional probability raises till it reaches its maximum approximately at the mean recurrence time. This is illustrated in Fig. 6 for the three values of aperiodicity. In this study, to quantify the impact on hazard of the cycling of main events at different elapsed times, the impact was computed for four distinct cases. In particular, we take into consideration: (i) the elapsed time equals to the mean of the distribution (i.e. 9500 yr); (ii) the elapsed time is 10 per cent of the mean recurrence time (i.e. 950 yr); (iii) without knowledge of the elapsed time and (iv) for sensitivity purpose only, the elapsed time is equal to 10 yr. In Fig. 6, the shape on the hazard-rate function for the three values of aperiodicity is plotted.

Once again, the impact is computed as the normalized difference between PGA-value of the BPT+ETAS-catalogues and the ones of the ETAS-catalogue. The estimation are done for 50 yr of exposure time. Here we note that the ETAS-catalogues are composed by background tectonic events plus triggered seismicity (aftershocks). The results for Cologne and Aachen are similar, hence only the one for Cologne are shown. Fig. 7 displays the impact of each of the four-cases described in the previous paragraph. In case i), the im-

pact on hazard is different for the three values of the aperiodicity; in particular, only the function with a very regular behaviour (aperiodicity equal to 0.3, see above and Fig. 6) gives a high impact, while the impact is negligible for the other two cases. The impact is driven by the increase of the probability of occurrence for the next BPT event in the next 50 yr (Fig. 6). In case that the elapsed time is 950 yr (case ii) or in the case that we have no knowledge of the elapsed time (case iii), the impact is almost null. The results for case iv) indicates that the impact at 90 per cent probability of non-exceedance is 10 per cent regardless of the aperiodicity. This kind of behaviour can be attributed to the aftershock activity of a large event, since from Fig. 6 it can be seen that the probability for a characteristic event is almost negligible for an elapsed time of 10–60 yr. Thus even in the time span 10–60 yr after a strong earthquake, the seismic hazard is significantly enlarged due to the ongoing aftershock sequence of the past main shock.

Numerical checks have been done with aperiodicity equal to 0.1, 0.4 and for the elapsed time equal to the mode of the distribution, showing the consistency of our results. In fact, the inclusion of the temporal cycling of main events is significant for hazard study only for very regular behaviour of the occurrence (e.g. small aperiodicity) at an elapsed time equal to the mean recurrence time of the distribution.

### 5.3.1 Implications for present hazard estimates

As a final application, we calculate the impact of the hypothesis of quasi-periodic recurrence of a $M_w$7 event on the given fault for the present hazard estimation for Cologne and Aachen. On the basis of the historic records available so far for this area, we have to assume that the elapsed time for an event of $M_w \ge 7.0$ is at least 1200 yr. The impact is evaluated considering BPT-ETAS catalogues with the elapsed time larger than 1200 yr, and leaving the aperiodicity as a
<figure><img src="image_9.png" /><figcaption>**Figure 7.** The impact of the cycling main shock recurrence on hazard for different elapsed times and for three different values of the aperiodicity. Solid lines: BPT with aperiodicity 0.3; dashed lines: BPT with aperiodicity 0.5; dotted lines: BPT with aperiodicity 0.7. Colour code: light grey for the elapsed time equal to the mean recurrence time; black for no-knowledge of the elapsed time and dark grey for the elapsed time equal to 10 yr. The results for elapsed time equal to 10 per cent of the mean are not shown since they are similar to the one without knowledge on the elapsed time.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 8.** The impact on present day hazard in LRE where we only know that the elapsed time is ≥1200. The aperiodicity is a free parameter: *a* = 0.3 (solid lines); *a* = 0.5 (dashed lines) and *a* = 0.7 (dotted lines).</figcaption></figure>

free parameter. The result, shown in Fig. 8, indicates that the impact
of such an assumption is almost negligible. We point out that this
is a first and simple model to represents the LRE seismicity. We
think that future works need to be address in this direction in order
to better represents this tectonic structure.

# 6 SUMMARY AND CONCLUSIONS

Our study aims at quantifying the impact of the time-variability of
seismicity on seismic hazard in low seismicity regions. Our specific
target region is an area in Germany, the LRE.
At first we focus on the short-term behaviour by including after-
shock activity into the hazard study via the ETAS model. To test the
applicability of this model, we have performed three statistical tests:
the analysis of the CV, the non-parametric study of the hazard-rate
function and the comparison of the AIC for the ETAS model and
the Poisson one, Sections 3.1, 3.2 and 3.3. Our analysis shows that
neglecting aftershocks leads to an underestimation of the hazard of
8 per cent at 90 per cent probability of non-exceedance in 50 yr. Here
we point out that, although the impact is not very high and taking
into account that each other ingredient within the chain of the hazard
computation has further uncertainties, this result is important since
the time-dependent behaviour leads to a systematic underestimation
of the estimated values. Moreover, the ongoing aftershock sequence
of the Roermond event occurred in 1992 still contributes 10–15 per
cent to the hazard at the level of 90 per cent of probability of non-
exceedance. Even the Düren earthquake which occurred 250 yr ago
still contributes about 20 per cent to the present hazard for the city
of Aachen at the level of 95 per cent probability of non-exceedance.

The second test was done with the purpose to examine the effect
of the repeating occurrence of main shocks. We study the behaviour
of a generic fault that follows the BPT distribution. Especially, the
influence of the aperiodicity and the elapsed time on the hazard
estimation has been evaluated. Knowing the fault and the statistic of
the earthquake time recurrence times (i.e. the elapsed time and the
BPT distribution) the impact on hazard ranges from 5 to 10 per cent
for the 90 per cent probability of non-exceedance in 50 yr. For the
two special cases of an ongoing aftershock sequence on the one hand
and an elapsed time equal to the mean recurrence time on the other
hand, the increase of the hazard value is more than 10 per cent at the
level of 90 per cent probability of non-exceedance. However, after
the complete decay of the aftershock activity and for aperiodicity
values larger then 0.3, the impact on hazard is negligible for any
elapsed time.

It is important to note that, given our limited knowledge of the
seismicity of this area, the assumption of a quasi-periodic main
shock occurrence on a fault does not change the present hazard
significantly. Thus on the basis of this study, we can deduce that
the implementation of a ‘long-term’ dependence of the seismicity
via the BPT distribution does not substantially affect the hazard
estimation for the LRE region.

ACKNOWLEDGMENTS

This work was supported by the Deutsche Forschungsgemeinschaft
(SCHE280/14-2). We thank Dr G. Grünthal and Dr K. G. Hinzen
for reviewing a preliminary version of the paper and for the fruitful
comments. We thank the editor Andrew Curtis and two anonymous
reviewer for valuable suggestion.

REFERENCES

Akaike, H., 1974. A new look at the statistical model identification, *IEEE Trans. Auto. Control*, **19**, 716–723.

Ambraseys, N.N., 1990. Uniform magnitude re-evaluation of European earthquakes associated with strong motions, *Earthq. Eng. Struct. Dyn.*, **19**, 1–20.

Ang, A. & Tang, W., 1975. *Probability Concepts in Engineering Planning and Design*, Vol. I: Basic Principles, John Wiley & Sons, UK.

Beauval, C., Hainzl, S. & Scherbaum, F., 2006. Probabilistic seismic haz-
ard estimation in low-seismicity regions considering non-Poissonian seis-
mic occurrence, *Geophys. J. Int.*, **164**, 543–550, doi:10.1111/j.1365-
246X.2006.02863.x

Berge-Thierry, C., Cotton, F., Scotti, O., Griot-Pommere, D.A. & Fukushima, Y., 2003. New empirical response spectral attenuation laws for moderate European earthquakes, *J. Earthq. Eng.*, **7**, 193–222.

Camelbeech, T. & Meghraoui, M., 1998. Geological and geophysical evidence for large paleo-earthquakes with surface faulting in the Roer Graben (northwest Europe), *Geophys. J. Int.*, **132**, 347–362.

Camelbeech, T., Alexandre, P., Vanneste, K. & Meghraoui, M., 2000. Long-
term seismicity in regions of present day low seismic activity: the example
of western Europe, *Soil Dyn. Earthq. Eng.*, **20**, 405–415.

Cinti, F.R., Faenza, L., Marzocchi, W. & Montone, P., 2004. Probability map of the next $M \ge 5.5$ earthquakes in Italy, *Geochem. Geophys. Geosyst.*, **5**, Q11003, doi:10.1029/2004GC000724.

Cornell, C.A., 1968. Engineering seismic risk analysis, *Bull. seism. Soc. Am.*, **58**, 1583–1606.

Cox, D.R. & Lewis, P.A.W., 1966. *The Statistical Analysis of Series of Events*, Methuen, London.

Ebel, J.E., Bonjer K-P. & Oncescu M.C., 2000. Paloseismicity: seismic-
ity evidence for past large earthquakes, *Seism. Res. Lett.*, **71**, 283–
294.

Faenza, L., 2005. Analysis of the spatio-temporal distribution of large earth-
quakes, *PhD thesis*, Università degli Studi di Bologna, *Alma Mater Stu-
diorum*, Bologna, Italy.

Felzer, K.R. & Brodsky, E.E., 2006. Decay of aftershock density with distance indicates triggering by dynamic stress, *Nature*, **441**, doi:10.1038/nature04799.

Felzer, K.R., Becker, T.W., Abercrombie, R.E., Ekström, G. & Rice, J., 2002. Triggering of the 1999 Mw 7.1 Hector Mine earthquake by aftershocks of the 1992 Mw 7.2 Landers earthquake, *J. geophys. Res.*, **107**, B9, 2190, doi:10.1029/2001JB000911.

Frankel, A., 1995. Mapping Seismic Hazard in the Central and Easter United States, *Seism. Res. Lett.*, **66**, 8–21.

Giardini, D., Wiemer, S., Fäth, D. & Deichmann, N., 2004. Seis-
mic hazard assessment of Switzerland, 2004, Report of ETH, Zurich
(http://www.earthquake.ethz.ch/research/Swiss / 'Hazard).

Grünthal, G. & Wahlström, R., 2003. An Mw based earthquake catalogue for central, northern and northwestern Europe using a hierarchy of magnitude conversions, *J. Seimol.*, **7**, 507–531.

Grünthal, G. & Wahlström, R., 2006. New generation of probabilistic seismic hazard assessment of the Cologne/Aachen considering the uncertainties of the input data, *Nat. Hazards*, **38**, 159–176.

Hainzl, S., Beauval, C. & Scherbaum, F., 2006. Estimating background ac-
tivity based on inter-event-time distribution, *Bull. seism. Soc. Am.*, **96**,
313–320.

Helmstetter, A.M, Ouillon, G. & Sornette, D., 2003. Are aftershocks of
large California earthquakes diffusing?, *J. geophys. Res.* **108**, 2483,
doi:10.1029/2003JB002593.

Hinzen, K.G., 2004. The use of engineering seismological models to interpret archeoseismolgical findings in Tolbiacum, Germany, a case study, *Bull. seism. Soc. Am.*, **95**, 521–539.

Hinzen, K.G. & Oemisch, M., 2001. Location and magnitude from seismic intensity data of recent and historic earthquakes in the Northern Rhine Area, Central Europe, *Bull. seism. Soc. Am.*, **91**, 40–56.

Kagan, Y.Y. & Jackson, D.D., 2000. Probabilistic forecasting of earthquakes,
*Geophys. J. Int.*, **143**, 438–453.

Kagan, Y.Y., Rong, Y.F. & Jackson, D.D., 2003. Probabilistic forecasting of seismicity, in *Earthquake Science and Seismic Risk Reduction*, Chapter 5.2, pp. 185–200, eds Mulargia, F. & Geller, R.J., Kluwer, Dordrecht.

Kalbfleisch, J.D., 1985. *Probability and Statistical Inference*, 2nd edn,
Springer-Verlag, New York.

Leydecker, G., 2005. Earthquake Catalogue for the Federal Republic of Germany and Adjacent Areas for the Years 800-2003, Datafile, http://www.bgr.de/quakecat. Federal Institute for Geoscience and Natural Resources

Matthews, M.V., Ellsworth, W.L. & Reasenberg, P.A., 2002. A Brownian Model for the recurrent earthquakes, *Bull. seism. Soc. Am.*, **92**, 2233–2250.

McGuire, R.K., 1976. Fortran computer program for seismic risk analysis,
*US Geological Survey open-File Report 76-67*.
MPS Working Group, 2004. Redazione della mappa di pericolosità sismica prevista dall'Ordinanza PCM 3274 del 20 marzo 2003. Rapporto Conclusivo per il Dipartimento della Protezione Civile, INGV, Milano-Roma, 2004 aprile, 65 pp. + 5 appendixes. http://zonesismiche.mi.ingv.it

Musson, R.M.W., 1999. Probabilistic seismic hazard maps for the North
Balkan Region, *Annali di Geofisica*, **42**, 1109-1124.

Ogata, Y., 1988. Statistical models for earthquake occurrence and residual analysis for point processes, *J. Am. Stat. Assoc.*, **83**, 9–27.

Ogata, Y., 1998. Space-time point-process models for earthquake occurrence, *Ann. Inst. Statis. Math.*, **50**, 379–402.

Ogata, Y., Matsuúra, R.S., Katsura, K., 1993. Fast likelihood computation of the epidemic type aftershock-sequence model, *Geophys. Res. Lett.*, **20**, 2143–2146.

Reiter, L., 1990. *Earthquake Hazard Analysis: Issues and Insight*, Columbia University Press, New York.

Rosenblueth, E., 1964. Probabilistic design to resist earthquakes, *J. Eng. Mech. ASCE*, **90**, no. EM5, 189–220.

Schmedes, J., Hainzl, S., Reamer, S.K., Scherbaum, F. & Hinzen, K.G., 2005. Moment release in the Lower Rhine Embayment, Germany: seismological perspective of the deformation process, *Geophys. J. Int.*, doi:10.111.j.1365-246X.2005.02525.x

Senior Seismic Hazard Analysis Committee (SSHAC), 1997. Recommendations for probabilistic seismic hazard analysis-guidance on uncertainty and use of experts, U.S. Nuclear Regulatory Commission, NUREG/CR-6372.

Smith, W.D., 2003. Earthquake hazard and risk assessment in New Zeland by Monte Carlo methods, *Seism. Res. Lett.*, **74**, 298–304.

Tanner, M.A. & Wong, W.H., 1984. Data-based Non-parametric Estimation of the Hazard Function with Application to a Model Diagnostics and Exploratory Analysis, *J. Am. Stat. Ass.*, **79**, 174–182.

Tyagunov, S., Grünthal, G., Wahlström, R., Stempniewski, L. & Zschau, J., 2006. Seismic risk mapping of Germany, *Nat. Hazard Earth Syst. Sci.*, **6**, 573–596.

Utsu, T., Ogata, Y. & Matsuúra, R.S., 1995. The centenary of the Omori formula for the decay law of aftershock activity, *J. Phys. Earth*, **43**, 1–33.

Working Group on California Earthquake Probabilities, 2003. Earthquake Probabilities in the San Francisco Bay Region (U.S. Geologic Survey, Denver), *USGS Open File Report* 03–214.

Zhuang, J., Ogata, Y. & Vere-Jones, D., 2002. Stochastic declustering of space-time earthquake occurrences, *J. Am. Stat. Assoc.*, **97**, 369–380.