## 3 Results and Discussion
Typically, these spatiotemporal variations could not be modeled using only linear functions with a few parameters, as they could not catch the complex nonlinear variations present in the seismicity phenomenon. However, it is still possible to use several types of parameterized linear functions, such as the model in Equation (1), where the interaction effects of the original covariates are included through auxiliary covariates. \begin{table}
\begin{tabular}{c c c c} \hline \hline
**Coef.** & **Mode** & **Mean** & **95\% CI** \\ \hline \(\beta_{(1)0}\) & 0.8136 & 0.7004 & (0.4677, 0.8079) \\ \(\sigma\) & 0.5710 & 1.1110 & (1.0855, 1.1697) \\ \(\kappa\) & \(-\)0.5900 & \(-\)0.5930 & (\(-\)0.6014, \(-\)0.5900) \\ \(\sigma_{1}^{2}\) & 26.1363 & 83.3487 & (26.2751, 810.0938) \\ \hline \hline \end{tabular}
\end{table} 
Table 2: Estimates and 95% credible intervals of the non-stationary GEV model for the maxima of earthquakes. Figure 4: Estimated spatial smoothing of the location parameter for the years: (**a**) 1995; (**b**) 2000; (**c**) 2005; (**d**) 2010; (**e**) 2015; (**f**) 2018. Several results were obtained by analyzing the spatial smoothing set for the location parameter over time, as shown in Figure 4a-f. Firstly, we could observe that each of the functions estimated through time showed several peaks or maximums, which correspond to the places where the most intense earthquakes occur. Secondly, it can be observed that the estimated functions for the years 1995, 2000 and 2005 have similar characteristics. However, in the year 2010, the patterns began to change, until later, for the years 2015 and 2018, we observed different patterns compared to the first few years. Furthermore, the intensities of earthquakes increase over time, which shows that the risk of extreme events increases with time. The flexibility of the proposed model allows studying the temporal variations in each geographical location on the study area. Figure 5 shows the temporal dynamics at six geographical locations in which we observed a period of increased seismic activity between the years 1998 and 2002 and a period of greater seismic activity between 2010 and 2016 that remained at high levels until the last years of the period of time studied. Consistently with these results, in 2017 and 2018, three of the strongest earthquakes of the studied period were recorded. A strong earthquake of \(M_{w}8.2\) hit the southern coast of Mexico, on 7 September 2017, approximately 87 kilometers southwest of Pijijapan in the state of Chiapas, being the second strongest recorded in the history of the country. Two weeks later, on 9 September 2017, a second earthquake of magnitude \(M_{w}7.2\), struck the states of Puebla and Morelos, in the center of the country, causing severe material and human losses. A few months later, on 16 February 2018, a strong earthquake of \(M_{w}7.2\) hit the south coast of Mexico again, approximately 37 km northeast of Pinotepa de Don Luis in the State of Oaxaca. The probability of large earthquakes on the coasts of Oaxaca and Guerrero is greater compared to other regions, due to the proximity they have with boundary plates. However, the Puebla and Morelos earthquake hit the central zone of the country, far from the limits between tectonic plates. Recent research suggests that these earthquakes should not be unexpected. An entropy change of seismicity under time reversal, considered as a non-seismic precursor of earthquakes, was identified several months before the occurrence of the major earthquake [34; 35]. In this study, we developed a model of extreme values with censoring data for the earthquakes on the Pacific coast of southern Mexico using a robust and solidly-based theory. An important application of the results of our study is through a risk map or return level map.