In this subsection, several relevant seismic parameters will be introduced as they are indicators of the seismicity of an area. These parameters are the maximum magnitude, the b-value, and the normalized annual rate for events exceeding M3.

#### 3.2.2. The b-value

The

b-

value is the parameter most studied in seismology and corresponds to the slope of the FMD in a log–log plot. The majority of the authors stated its relationship with the physics of the studied area. It has been deeply studied and is usually analyzed in PSHA works, in prediction, in locating asperities, periodical tidal loading, and energetic characterization [

4,

8,

11,

16,

25,

92,

93,

94,

95].

Although firstly, the least square solution was usually employed to calculate its value, currently, there is a consensus in considering the Maximum-Likelihood-Estimate (MLE) as the best approach to obtain it. This is due to the fact that it does not present interdependency between variables [

2]. Over time, a considerable amount of different methods have been suggested for its calculations, such as those found in [

13,

96,

97,

98,

99,

100]. One of the most employed formulae for MLE was proposed by [

98,

99] for binned data.

where

M_{c} is the cut-off magnitude,

$\overline{M}$ is the average magnitude of the earthquakes whose magnitude is larger than or equal to

M_{c}, and Δ

M is the binning interval of the magnitude.

In this research, the bin interval is 0.1. The solution proposed by Kijko and Smit [

13], through the exposed MLE expression, has been applied. It permits taking into account a longer temporal extent and different magnitude–year of completeness pairs. Moreover, it is said to be simple, manageable, and it is not based on iterations [

13].

The method suggested by Kijko and Smit [

13] is based on dividing the catalog into more coherent

s sub-catalogs, each of a different level of completeness, and with its corresponding year of completeness. It is particularly indicated for incomplete or inhomogeneous catalogs. For every sub-catalog, the MLE proposed by [

98,

99] is used. Later,

b-

value is estimated as a weighted solution as:

where

b_{i} is the

b-

value of each of the

s sub-catalogs,

n_{i} is the sample size of the sub-catalogs, and

n is the total number of events considered (

n = n_{1} + n_{2} + …+ n_{s}).

Each sub-catalog has a known but different level of completeness, ${M}_{min}^{1},{M}_{min}^{2},\dots ,{M}_{min}^{s}$, and it spans ${t}_{1},{t}_{2},\dots ,{t}_{s}$ years.

Finally, from now on, the obtained

b-

value after the correction proposed by [

101], to minimize the overestimation produced for small samples, will be noted as

b-

value, or only

b.

Besides, despite the method proposed by Kijko and Smit [

13] gives the expression to calculate the standard deviation, that suggested by Shi and Bolt [

102] has been preferred as it considers the real dispersion of the sample:

where

n is the total samples.

As stated above, the method requires pairs of values to estimate the

b-

value. These have been previously calculated and are shown in

Table 2.

Once the formulae have been defined, the minimum number of events to generate a representative

b-

value must be established. It is a very controversial issue, and there is no general academic agreement: extraordinarily, Dominique and Andre [

103] considered only six events: Bender [

98] or Skordas and Kulhánek [

104] chose 25 events; Mousavi [

11] and Amorèse et al. [

28] proposed 50; González [

60] suggested 60, and Roberts et al. [

105] established 200. A thorough study regarding the relationship between error and number of events for different

b-

values can be found in Nava et al. [

24].

When a continuous representation is prepared, another crucial parameter is how the geographical space is divided. This split can lead to completeness issues, as the number of events for every cell may be small, particularly when a grid division is employed [

106], as in [

7,

11,

107]. The minimum number of events is a trade-off between accuracy and coverage, whereas cell size is a trade-off between coverage and resolution [

11]. In this work, using a GIS tool, two different grid sizes have been established as in [

7,

11]. Considering these trade-offs, a 0.5 × 0.5° grid was selected with 100 events as a minimum value for the most active area; and a 1 × 1° grid has been considered for the whole area, with a minimum of 50 events. Besides, to reduce the border effect, four overlapped grids have been defined (the original; one shifted half of cell size to the south; another moved half to the east; and finally, displaced in south and east), as was done previously in [

7,

15].

For every cell, the average geographical coordinates of the epicenters have been estimated. Thus, seismic parameters have been assigned to this location. Later, the method proposed by Kijko and Smit [

13] has been adopted to compute the

b-

value. This has been for every cell of every grid. Finally, to conduct the spatial analysis, an interpolation by the Inverse Distance Weighting (IDW) algorithm has been applied. It must be noted that where the minimum number of events has not been reached, its associated

b-

value has not been considered in the interpolation, so caution should be taken when analyzing these areas.

Different color maps have been produced with the conditions above exposed to represent the b-value distribution. The location of the points used to generate the maps is shown.