# The GPS velocity field of the Aegean. New observations, contribution of the earthquakes, crustal blocks model

Pierre Briole, Athanassios Ganas, Panagiotis Elias, Dimitar Dimitrov

► To cite this version:

Pierre Briole, Athanassios Ganas, Panagiotis Elias, Dimitar Dimitrov. The GPS velocity field of the Aegean. New observations, contribution of the earthquakes, crustal blocks model. Geophysical Journal International, 2021, 226 (1), pp.468-492. 10.1093/gji/gab089. hal-03409600

HAL Id: hal-03409600

https://hal.science/hal-03409600

Submitted on 9 Aug 2022

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_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# The GPS velocity field of the Aegean. New observations, contribution of the earthquakes, crustal blocks model

Pierre Briole<sup>1</sup>,<sup>2</sup> Athanassios Ganas,<sup>2</sup> Panagiotis Elias<sup>1,3</sup> and Dimitar Dimitrov<sup>4</sup>

Accepted 2021 March 2. Received 2021 February 23; in original form 2020 July 27

## SUMMARY

We calculate and analyse the coordinate time-series of 282 permanent GPS stations located in Greece and 47 in surrounding countries. The studied period is 2000–2020. The average GPS time-series length is 6.5 yr. The formal velocity uncertainties are rescaled to be consistent with the velocity scatters measured at 110 pairs of stations separated by less 15 km. We remove the effect of the crustal earthquakes of $M_w \ge 5.3$. We quantify and model the post-seismic deformations. Two relaxation times are usually needed: one short of some weeks and one long of 1 yr or more. For the large $M_w = 6.9$ events of Samothraki 2014 and Methoni 2008, the post-seismic deformation equals or exceeds the coseismic one. We detect at three stations a deformation transient in May 2018 that may correspond to a slow earthquake beneath Zakynthos and northwest Peloponnese, with equivalent magnitude 5.8. The density and accuracy of the velocities make it possible to better quantify several characteristics of the deformation in the Aegean, in particular: (i) the transition from the Anatolian domain, located in the southeast, to the European domain through the western end of the North Anatolian fault; (ii) the north–south extension in the western Aegean; (iii) the east–west extension of the western Peloponnese; (iv) the clockwise rotation of the Pindos; (v) the north–south extension in central Macedonia. Large parts of the central Aegean, eastern Peloponnese and western Crete form a wide stable domain with internal deformation below 2 nstrain yr$^{-1}$. We build a kinematic model comprising 10 crustal blocks corresponding to areas where the velocities present homogeneous gradients. The blocks boundaries are set to fit with known localized deformation zones, for example, the rift of Corinth, the North Anatolian fault and the Katouna fault. When the velocity steps are clear but not localized, for example, through the Peloponnese, the boundary line is arbitrary and represents the transition zone. The model fits the velocities with a root-mean-square deviation of $\pm0.9$ mm yr$^{-1}$. At the boundaries between blocks we compare the predicted and observed deformations. We find shear rates of 7.4 and 9.0 mm yr$^{-1}$ along the Movri and Katouna faults, 14.9 and 8.7 mm yr$^{-1}$ along the North Anatolian fault near Lemnos and near Skopelos respectively, extension of 7.6, 1.5 and 12.6 mm yr$^{-1}$ across the Gulf of Patras, the Trichonis Lake and the Ambracian Gulf. The compression across western Epirus is 3.7 mm yr$^{-1}$. There is a dextral transtensional movement of 4.5 mm yr$^{-1}$ between the Amorgos and Astypalea islands. Only the Ionian Islands region shows evidence of coupling along the subduction interface.

**Key words:** Space geodetic surveys; Transient deformation; Europe; Earthquake ground motions; Fractures, faults, and high strain deformation zones; Kinematics of crustal and mantle deformation.

<footnote n="1"><sup>1</sup>UMR 8538 CNRS - École Normale Supérieure - PSL Research University, F-75005 Paris, France. E-mail: briole@ens.fr</footnote>

<footnote n="2"><sup>2</sup>National Observatory of Athens, Institute of Geodynamics, GR-11810 Athens, Greece</footnote>

<footnote n="3"><sup>3</sup>National Observatory of Athens, IAASARS, GR-15236 Penteli, Greece</footnote>

<footnote n="4"><sup>4</sup>Bulgarian Academy of Sciences, National Institute of Geophysics, Geodesy and Geography, BG-1113 Sofia, Bulgaria</footnote>
# 1 INTRODUCTION

Since its origin, space geodesy has been used in Greece for a wide range of geophysical applications (Marinou et al. 2016). These include sea level monitoring (Becker et al. 2002), quantifying of plate tectonics (Le Pichon et al. 1995) and study of the seismicity (e.g. Jackson et al. 1994; Meyer et al. 1996; Bernard et al. 1997). The first technique used was the Satellite Laser Ranging (SLR; Aardoom & Van Gelder 1984; Christodoulidis et al. 1985; Smith et al. 1994) followed soon by the Global Positioning System (GPS; Noomen et al. 1996; Clarke et al. 1998, Cocard et al. 1999). The Doppler Orbitography and Radiopositioning Integrated by Satellite (DORIS) system was used for the calibration of satellite altimetry and GPS measurements (Willis et al. 2013). Among the first GPS works, it is worth mentioning Billiris et al. (1991) who compared the coordinates of first-order geodetic points measured by triangulation at the end of the 19th century with new coordinates obtained by GPS. This work has been the first to quantify the extension in central Greece across the Gulf of Corinth. Synthetic Aperture Radar Interferometry (InSAR) has been widely used for the study of earthquakes most of the time coupled with GPS (e.g. Briole et al. 2000, 2015).

After the appearance of a small number of permanent GPS stations in the 1990s, their number increased during the 2000s with the deployment of stations throughout the territory (Ganas et al. 2008) by the National Observatory of Athens (NOA) and in greater density in the Corinth rift by the Centre National de la Recherche Scientifique (CNRS). At the end of the 2000s, three nationwide, permanent networks were created for the needs of land use planning and supporting the work of the surveyors doing geodesy in Greece: the Hepos (https://www.hepos.gr) network funded by the European Union and managed by the Greek State through its Cadasstre Department, and two networks operated by the Metrica SA (HxGN SmartNet, https://wwwmetrica.gr) and TreeComp (Uranus network, http://www.uranus.gr) companies. Despite its status, as of today, the first network is not freely open to the science community (only a portion of its data is available on request) whereas we have access to the data from the two others, which represent 180 stations.

We analysed the GPS data acquired at 329 stations (Fig. 1, Supporting Information Table S1), 282 located in Greece, 28 in Turkey, 12 in Bulgaria, 6 in Albania and 1 in Northern Macedonia. In Greece, 96 stations belong to TreeComp, 84 to Metrica and 102 to various Greek, European and foreign research centres: 30 to the French CNRS, 18 to the National Observatory of Athens (NOA), 13 (not any more in use) to the British COMET (Billiris et al. 2014), 8 to the National Technical University of Athens (NTUA), 6 to the Charles University of Prague (CUP), 6 to the University of Corado (UNAVCO), 5 to the Italian Istituto Nazionale di Geofisica e Vulcanologia (INGV), 5 to the Technological Educational Institute of Crete (TEIC) and 11 to various others. Five of those academic stations, AUT1, DYNG, KATC, PAT0 and TUC2 (see Fig. 1 for locations), are included in both the ITRF2014 solution (Altamimi et al. 2016) and the EUREF network (https://www.epncb.oma.be/). The Turkish stations are part of the General Command of Mapping (GCM) network. The data cover the period 2000–2020 and especially 2010–2020 because of the larger number of GPS stations.

With respect to the previous velocities fields (e.g. Clarke et al. 1998; Kahle et al. 2000; Kreemer & Chamot-Rooke 2004; Nyst & Thatcher 2004; Serpelloni et al. 2005; Reilinger et al. 2006; Hollenstein et al. 2006; Floyd et al. 2010; Reilinger et al. 2010;

Müller 2011; Nocquet 2012; Pérouse et al. 2012; Ganas et al. 2013; Müller et al. 2013; Chousianitis et al., 2013, 2015; Vernant et al. 2014; Métois et al. 2015; England et al. 2016; Devoti et al. 2017; Bitharis et al. 2019; Ganas et al. 2019, D’Agostino et al. 2020), this one is much denser. It is more accurate because it is based entirely on data from permanent stations, and with longer time periods of acquisition. Moreover, it is improved, thanks to the correction of the deformations induced by the earthquakes.

In the first part of the paper, we describe our method for processing the data. Details on the processing are presented in Supporting Information Section S2.

In the second part, we quantify the coseismic and post-seismic contributions to the GPS velocity field of the Aegean area for the period 2000–2020. Indeed, the inaccurate correction for the coseismic and post-seismic displacements leads to biased estimates of secular velocities. This is visible, for example, at the EUREF station PAT0 (location in Fig. 1) where the seismic part of the velocity field accounts (for the 12 yr of observation) for 0.67 and 0.85 mm yr⁻¹ in east and north respectively (Supporting Information Fig. S24). In Supporting Information Section S3, we review the models of coseismic and post-seismic deformations of the earthquakes listed in Table 1. In Supporting Information Section S4, we discuss some other signals found in the time-series.

In the third part, we discuss the secular velocity field obtained after the removal of the contribution of the earthquakes.

In the fourth part, we present a crustal blocks kinematic model tailored to fit this secular velocity field. This model fits the velocities with a root-mean-square (r.m.s.) deviation of ± 0.9 mm yr⁻¹.

Fig. 2 shows on a map the location of the various sites mentioned in the text, and Supporting Information Table S2 contains a glossary of the site names, with some additional details.

# 2 GPS DATA PROCESSING

In the following we use the acronym GPS and not GNSS (Global Navigation Satellite System) because we have processed and analysed only the GPS data acquired by the receivers.

## 2.1 Calculation

We processed all the data except those from Turkey since they are not publicly available but only the time-series produced by the Nevada Geodetic Laboratory (NGL, http://geodesy.unr.edu; Blewitt et al. 2018). Details on the processing method are available in Supporting Information Section S2.

We analyse the velocities only for stations with an acquisition period greater than 2 yr and velocity accuracy better that 1 mm yr⁻¹, which sums up to 309 stations out of the 329. We provide in the additional material (Supporting Information Table S7) the coordinates of the GPS stations in ITRF2014 at the epoch 2020.0. For the 20 stations with short acquisition period, these coordinates are established by using, for their velocity, an interpolation of the velocities available at the neighbouring stations. Those 20 extra stations, having a short acquisition period, were kept because, in several cases, they were valuable for the measurement and modelling of coseismic deformations. Therefore, they contribute in this manner to the quality of the final velocity field.

We model the annual term of the coordinate time-series using a simple sine function of annual period. In this way only two parameters are needed, an amplitude and a phase. The time-series of coordinates of the 329 stations are available in ASCII format in the
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Location of the 329 GPS stations (see Table 1). The 47 magenta stations correspond to GPS stations outside Greece: in Albania (6), Northern Macedonia (1), Bulgaria (12) and Turkey (28). The 96 green stations belong to the Uranus network (http://uranus.gr) operated by the Tree Company corporation. The 84 orange stations belong to the Metrica company (https://metrica.gr). The 102 red stations are owned and operated by academic, or assimilated, institutions. They are denser in southwestern Greece, in the Corinth rift (inset) and in Santorini.</figcaption></figure>

Supporting Information. We applied a Gaussian filter $F(t) = \sigma^{-1}(2\pi)^{-1/2} \exp(-t^2/2\sigma^2)$ to the time-series presented in some of the figures, using the characteristic time $\sigma = 2.24$ d. This filter, which is involving weighted data acquired in the window of one week before and after each epoch, reduces the noise of the time-series (for both the horizontal and the vertical components) by a factor close to 2.

## 2.2 Assessing the uncertainties

In the case where the time-series does indeed show a linear trend, without marked fluctuations, which is almost always the case for horizontal components but not for verticals, the uncertainty depends on (i) the total duration of the observations, (ii) the mean dispersion of the daily solutions, (iii) the uncertainties related to the correction of possible jumps in the series and (iv) the uncertainties related to the correction of the coseismic and post-seismic motions.

In order to establish a realistic law for the uncertainties, we studied how the velocity differences between neighbouring stations evolve as a function of their distance. This was possible because the data set contains several clusters of neighbouring stations. We selected 110 pairs of stations located within 15 km of each other and analysed the scatter of their velocities in east, north and vertical (Supporting Information Table S3). We used those scatters to scale the formal uncertainties provided by the calculation so they fit with the observed scatters. Fig. 3 shows the alignment of the scatters and the uncertainties for the east, north and vertical components. Supporting Information Fig. S7 shows the histogram of the uncertainties in east and vertical. In the way that it is made, our estimation of the uncertainties does not take into account the velocity deviations that may exist in some cases between fairly close (typically a few km) stations due to the activity of local active faults during the period 2000–2020, or to any other local effect, for example, hydrological effects or site instability. This is the case, for example, at Aigion where the three nearby stations of AIGI, VALI and EGIO
<table><caption>**Table 1.** Significant earthquakes, from 2001 to 2020, considered for the calculation of coseismic displacements in the GPS network. Magnitude, strike, dip and rake are from G-CMT (Ekström *et al.* 2012). In the case of Zakynthos 2006 and Gulpinar 2017, the various focal mechanisms of the sequence are very similar, and we used the average angles. Longitude, latitude and depth are those inferred from the modelling of the GPS data. Time constant of post-seismic deformation is provided when a post-seismic signal was observed and quantified. The Samos earthquake of 2020 October 30 is mentioned and plotted because we mention it in the discussion but we do not evaluate the fault parameters of the co- and post-seismic displacements.</caption><thead><tr><th rowspan="2">Name</th><th rowspan="2">M<sub>w</sub></th><th rowspan="2">Date</th><th rowspan="2">Long. (°)</th><th rowspan="2">Lat. (°)</th><th rowspan="2">Depth (km)</th><th rowspan="2">Strike (°)</th><th rowspan="2">Dip (°)</th><th rowspan="2">Rake (°)</th><th colspan="2">Fast post-seismic</th></tr><tr><th>(d)</th><th>(d)</th></tr></thead><tbody><tr><td>Skyros</td><td>6.4</td><td>26/07/2001</td><td>24.34</td><td>38.95</td><td>9</td><td>148</td><td>76</td><td>-1</td><td></td><td></td></tr><tr><td>Lefkada</td><td>6.2</td><td>14/08/2003</td><td>20.63</td><td>38.83</td><td>9</td><td>18</td><td>59</td><td>-174</td><td></td><td></td></tr><tr><td>Zakynthos (six events)</td><td>5.0-5.7</td><td>3-19/4/2006</td><td>20.95</td><td>37.68</td><td>14.5</td><td>3</td><td>29</td><td>122</td><td></td><td></td></tr><tr><td>Methoni</td><td>6.9</td><td>14/02/2008</td><td>21.80</td><td>36.17</td><td>26.5</td><td>332</td><td>6</td><td>120</td><td>10</td><td>449</td></tr><tr><td>Movri</td><td>6.4</td><td>08/06/2008</td><td>21.52</td><td>37.94</td><td>22</td><td>209</td><td>83</td><td>164</td><td>103</td><td></td></tr><tr><td>Efpalion</td><td>5.3</td><td>18/01/2010</td><td>21.94</td><td>38.44</td><td>5.5</td><td>296</td><td>49</td><td>-66</td><td></td><td></td></tr><tr><td>Efpalion</td><td>5.3</td><td>22/01/2010</td><td>21.98</td><td>38.43</td><td>5.5</td><td>296</td><td>51</td><td>-59</td><td></td><td></td></tr><tr><td>Western Crete</td><td>6.4</td><td>12/10/2013</td><td>23.37</td><td>35.30</td><td>22</td><td>339</td><td>3</td><td>130</td><td></td><td></td></tr><tr><td>Cephalonia</td><td>6.1</td><td>26/01/2014</td><td>20.42</td><td>38.18</td><td>5.5</td><td>19</td><td>76</td><td>176</td><td>85</td><td></td></tr><tr><td>Cephalonia</td><td>6.0</td><td>03/02/2014</td><td>20.41</td><td>38.24</td><td>3</td><td>12</td><td>57</td><td>157</td><td></td><td></td></tr><tr><td>Samothraki</td><td>6.9</td><td>24/05/2014</td><td>25.68</td><td>40.28</td><td>5</td><td>73</td><td>85</td><td>-177</td><td>65</td><td>1295</td></tr><tr><td>Lefkada</td><td>6.5</td><td>17/11/2015</td><td>20.45</td><td>38.76</td><td>5</td><td>22</td><td>64</td><td>179</td><td>138</td><td>471</td></tr><tr><td>Gulpinar (four events)</td><td>5.2-5.4</td><td>6-12/02/2017</td><td>26.13</td><td>39.53</td><td>5.5</td><td>125</td><td>47</td><td>-89</td><td></td><td></td></tr><tr><td>Lesvos</td><td>6.3</td><td>12/06/2017</td><td>26.31</td><td>38.89</td><td>6.5</td><td>110</td><td>47</td><td>-87</td><td>4</td><td>140</td></tr><tr><td>Kos</td><td>6.6</td><td>20/07/2017</td><td>27.49</td><td>36.92</td><td>6.5</td><td>278</td><td>36</td><td>-82</td><td></td><td>220</td></tr><tr><td>Zakynthos</td><td>6.8</td><td>25/10/2018</td><td>20.70</td><td>37.40</td><td>6.5</td><td>11</td><td>28</td><td>165</td><td>21</td><td>183</td></tr><tr><td>Preveza</td><td>5.4</td><td>05/02/2019</td><td>20.61</td><td>38.93</td><td>10</td><td>2</td><td>81</td><td>180</td><td></td><td></td></tr><tr><td>Eratini</td><td>5.3</td><td>30/03/2019</td><td>22.27</td><td>38.29</td><td>10</td><td>292</td><td>39</td><td>-62</td><td></td><td></td></tr><tr><td>Magoula</td><td>5.3</td><td>19/07/2019</td><td>23.51</td><td>37.99</td><td>11</td><td>103</td><td>55</td><td>-89</td><td></td><td></td></tr><tr><td>Kanallaki</td><td>5.7</td><td>21/03/2020</td><td>20.61</td><td>39.26</td><td>7</td><td>337</td><td>39</td><td>119</td><td></td><td></td></tr><tr><td>Samos</td><td>7.0</td><td>30/10/2020</td><td>26.68</td><td>37.76</td><td>17</td><td>270</td><td>37</td><td>-95</td><td></td><td></td></tr></tbody></table>

(location in the inset of Fig. 1) have differential velocities of ~1 mm yr⁻¹ due to the local tectonic activity as shown by Elias & Briole (2018). However, this particular context is minor in our set of 110 pairs. Therefore, we believe that our uncertainties are well scaled and robust.

## 3 COSEISMIC AND POST-SEISMIC DEFORMATIONS

In the first subsection below, we review and quantify the coseismic and post-seismic contributions to the GPS velocity field of Greece and the Aegean area for the period 2000–2020. In the second subsection, we review the other phenomena found to disrupt some of the GPS time-series. The correction of those signals is needed to estimate a secular velocity field of the whole region.

### 3.1 Earthquakes disrupting the GPS time-series

Table 1 and Fig. 2 contain and display the list of earthquakes (or swarms for Zakynthos 2006 and Gulpinar 2017) occurred in Greece or close to Greece having produced a displacement larger than 1 mm at one at least of our 329 stations. To make this list, we extracted from the European Mediterranean Seismological Centre (EMSC, https://www.emsc-csem.org) catalogue, the earthquakes of magnitude ≥ 5.3 having occurred at a depth of less than 40 km between January 2000 and July 2020 in the box 20°–29° in longitude and 34°–42° in latitude. For each of those earthquakes, we used the Global Centroid Moment Tensor Catalog (G-CMT,

https://www.globalcmt.org/; Ekström *et al.* 2012) focal mechanism or, if unavailable, the NOA focal mechanism. We assume the following values for the length $L = 24.7 (M \mu^{-1})^{(1/3)}$, with width $W = 0.625 L$ and slip $D = M/(\mu L W)$, where M is the GCMT (or NOA) seismic moment (in N m) and μ the rigidity of the crust assumed to be equal to $3 \times 10^{10}$ Pa. This value of rigidity, currently used for the crust, is giving on average a good fit between the seismic and geodetic moments for the population of earthquakes analysed. Those scaling relations are used by Briole (2020) for the a priori estimate of fault parameters of crustal earthquakes in the Mediterranean aimed to produce a priori synthetic GPS displacements and interferograms after events, and to initialize inverse models. They are simple and, for the required need, sufficiently consistent with those published by Wells and Coppersmith (1994) and Thingbaijam *et al.* (2017), the former being based on models, while the second being constrained from field observations. Indeed, what is needed to estimate displacements in the intermediate and far field is not the detail of the fault geometry and slip, but the use of the correct focal mechanism and seismic moment provided by seismology. Moreover, the a posteriori analysis of the faults we modelled in this paper shows that their geometric parameters sometimes changed by more than 50 per cent with respect to their a priori value predicted by the above relations. This emphasizes the fact that these scale relations for seismic faults have a statistical value but cannot be used to precisely prescribe the geometric parameters of a particular event.

For each event, the predicted values were then used to calculate with the model of Okada (1992) the theoretical displacements at the 329 stations. This protocol eliminated all earthquakes occurred in the central and western Hellenic arc and Karpathos basin, several
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Toponymy (coordinates and description are in Supporting Information Table S2), topography, bathymetry, major faults from Ganas (2020), and G-CMT focal mechanisms of the events analysed in the text.</figcaption></figure>

with magnitude larger than 6, except the $M_w = 6.4$ earthquake of western Crete of 2013 October 12. Some moderate earthquakes of continental Greece such as the $M_w = 5.4$ and 5.3 Kallidromon earthquakes of 2013 August 7 and September 16 (Ganas et al. 2014), the $M_w = 5.3$ Amfilochia earthquake of 2014 October 24 and the $M_w = 5.3$ northern Evia earthquake of 2015 June 9 (Ganas et al. 2016a) were also eliminated after checking carefully the time-series of the nearby GPS stations.

The events of Table 1 are analysed one by one in Supporting Information Section S3.

## 3.2 Other transients disrupting the time-series

We identified transients in the time-series of several stations. The most noteworthy are described in Supporting Information Section S4. Most of them correspond to earthquakes. The most interesting of those events is a possible slow-slip event recorded at several stations close to the island of Zakynthos five months before the earthquake of 2018 October 25. We also re-analysed the GPS data

of Santorini. Thanks to the longer time-series, it is possible to revisit the unrest period, to correct its effect and to assess a secular velocity of the stations located on this volcanic island. This knowledge is important for the analysis of the volcanic and local tectonic deformations.

## 3.3 Separating transient and secular velocity fields

Among the five major earthquakes of the period 2000–2020, the ones that most disrupt the GPS time-series are those of Methoni 2008 February 14, Samothraki 2014 May 24, Lefkada 2015 November 17, Zakynthos 2018 October 25 and Samos 2020 October 30. The last one is not analysed in this paper.

Figure 4 shows the part of the velocity field induced by coseismic and post-seismic displacements during the sampled period at the GPS stations (only the 129 stations where the velocity is biased by more than 0.6 mm yr$^{-1}$ are plotted) and Supporting Information Table S7 contains the numerical values of this contribution of the earthquakes to the velocity field.
<figure><img src="image_5.png" /><figcaption>Figure 3. Estimating the uncertainties (in mm yr⁻¹) in the east, north and vertical axes by using the observed scatters between neighbouring stations as a metrics.</figcaption></figure>

The field of vectors shown in Fig. 4 is not a stationary field. It depends on the earthquakes that occured during the investigated period and on the duration of the available time-series, the shorter being more biased. Fig. 4 as a whole may however exhibit some characteristics likely to be relatively stationary in the long-term.

This is in particular the case of the overall southwestern velocity pattern of the stations of western Greece around the western Gulf of Corinth and the Ionian Islands. In that area, the contribution of the earthquakes to the velocity field exceeds 10 mm yr⁻¹ towards the southwest at stations not necessarily located in the near field of any recent earthquake. The most characteristic is the case of the EUREF station PAT0 (Supporting Information Fig. S24) where the velocity is affected with no abrupt anomaly in the time-series except for the Zakynthos earthquake of 2018.

The contribution of the earthquakes is also significant around Lemnos and Lesvos. The long duration of the post-seismic relaxation of the May 2014 Samothraki earthquake induces velocity anomalies at all stations located in northeastern Greece between Thessaloniki and Alexandroupoli. The correction of this effect is needed to properly measure and interpret the north-south extension of northern Greece and southern Bulgaria.

In Fig. 4, we distinguish two wide areas not affected by large earthquakes from 2000 to 2020, that is, below the threshold of 0.6 mm yr⁻¹: northwestern Greece and the Cyclades region.

While the time-series of the stations located close to the epicentres show clear offsets (which are actually the offsets used to constrain the models), it is not the case of the time-series of GPS stations located at larger distances. For many GPS stations located at intermediate distance (say 50 km), the coseismic and post-seismic effects can be present and significant in the time-series but hidden in the noise or in the yearly fluctuations. This is the case in particular for the post-seismic effects that are smooth, especially when the relaxation time is long. The corrections of the velocities can reach

~1 mm yr⁻¹ at stations located at several tens of km of the epicentre of the closest large event.

Supporting Information Table S6 shows how the velocity of PAT0 and the IGS stations AUT1, DYNG, KATC and TUC are biased by the earthquakes. Those biases are 0.2 and −0.2 mm yr⁻¹ at AUT1, −0.3 and −0.5 mm yr⁻¹ at DYNG, −0.7 and −0.9 mm yr⁻¹ at PAT0 (Supporting Information Fig. S24), 0.0 and −0.8 mm yr⁻¹ at KATC, and −0.4 and −0.6 mm yr⁻¹ at TUC2, in east and north, respectively.

Supporting Information Table S6 shows also a comparison of the velocities of four fundamental stations of Greece computed by five groups: LTK-EUREF, ROB-EUREF, UGA-CNRS (all three provided by the online EPOS portal https://gnssproducts.epos.ubi.pt/), NGL and ourselves. The velocities differ from one solution to another with the r.m.s. deviation of the solutions notably larger than the claimed velocity uncertainties in the case of LTK-EUREF and ROB-EUREF, the others being more conservative. Our solution before applying the correction of co- and post-seismic is relatively consistent with their average. The solution of NGL is closer to our final solution because NGL includes, like us, a routine correction for the coseismic motion of the largest earthquakes, yet using standard parameters provided by the seismological centres, and not an accurate modelling of the observed offsets as we do here.

The velocities in Supporting Information Table S1 are corrected for the effects of the earthquakes and other effects discussed in the next section.

## 3.4 Weight of the post-seismic relaxations

It appears that post-seismic relaxation is not the exception but the rule. In several cases, we can distinguish two distinct relaxation times: a first fast component, with a relaxation time of a few weeks
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Velocity biases corresponding to the coseismic and post-seismic signals for the period 2000–2020. The amplitudes of those biases are dependent on the duration of the available time-series. The projections of the faults of the earthquakes considered in Table 1 are plotted in red for the coseismic and green for the post-seismic (parameters of the faults in Table 2; use also the supporting KML file) and the six main ones labelled. The ellipses correspond to uncertainties at 95 per cent confidence level. The circle in the velocity scale bar corresponds to an uncertainty of 0.5 mm yr⁻¹. Only the 129 stations where the velocity is biased by more than 0.6 mm yr⁻¹ are drawn. The location of the five IGS and EUREF reference stations of AUT1, DYNG, KATC, PAT0 and TUC2 is shown.</figcaption></figure>

and sometimes even a few days (e.g. in Lesvos 2017); a second slower component, with a relaxation time ranging from several months to several years.

Although an exponential law allows to fit well the post-seismic relaxation, we found several cases in which it was needed to combine two exponentials, for example, in the case of the Samothraki 2014 earthquake (Supporting Information Fig. S2), one short relaxation time (65 d) and one long relaxation time (1295 d). The need to involve two time constants arises also for the earthquakes of Methoni 2008, Lefkada 2015 and Zakynthos 2018. For the other earthquakes, with magnitude less than $M_w = 6.5$, the amplitude of their post-seismic signals is less and the dual process not clearly apparent, even if it may also exist, and not needed to fit the data well below the level of noise of the time-series. Table 1 indicates the post-seismic relaxation times found for the studied earthquakes.

Figure 5 shows the percentage of the seismic moment released during the post-seismic phase for nine large events, that is, the quantity $100 M_p / (M_w + M_p)$ where $M_w$ is the coseismic moment and $M_p$ the post-seismic moment. For two of the largest earthquakes, $M_p$ is comparable or larger than $M_w$, up to three times larger for Methoni 2008 (see values in Table 3), while around the magnitude

$M_w = 6$ the average weight of $M_p$ is around 15 per cent. Our sample is however insufficient to attempt to establish a relationship between the magnitude and the weight of the post-seismic component.

As the relaxation time can be of several years, especially for the large events, the post-seismic effects are often smoothed in time and can be hard to distinguish from the secular velocity. Fig. 6 shows how, according to our fault parameters (Table 2) and relaxation times (Table 1), the seismic moment is released as a function of time. This ‘breathing’ impacts the velocities of the GPS stations that are affected by it, and it could well be invisible when using sparse campaign data instead of continuous data.

When modelling the post-seismic signal, we find that the short-term relaxation can always be modelled with a slight extension of the initial fault plane: laterally, upward and downward. For example, for Zakynthos 2018 the rapid post-seismic can be modelled by an extension of the slip on a vertical fault on top of the coseismic one connected at the depth of 3 km below the seafloor to the surface (Ganas et al. 2020). For long-period relaxation it is the case for Samothraki 2014 only, with the GPS data fitted by a 4.5 km downward extension of the slip on the fault plane, from 8.5 to 13 km depth. The cases of Lefkada 2015 and Zakynthos 2018 are different, with the post-seismic slip occurring on planes with
<table><caption>**Table 2.** Parameters of the modelled coseismic and post-seismic dislocations.</caption><thead><tr><th>Name</th><th>Date</th><th>Long (°)</th><th>Lat (°)</th><th>Depth (km)</th><th>Fault length (km)</th><th>Fault width (km)</th><th>Strike-slip (mm)</th><th>Dip slip (mm)</th><th>Fault azimuth (°)</th><th>Fault dip angle (°)</th></tr></thead><tbody><tr><td>Skyros<sup>0-1</sup></td><td>26/07/2001</td><td>24.349</td><td>38.950</td><td>5.2</td><td>16.6</td><td>7.6</td><td>1500</td><td>0</td><td>148</td><td>76</td></tr><tr><td>Lefkada<sup>2</sup></td><td>14/08/2003</td><td>20.601</td><td>38.840</td><td>4.5</td><td>16</td><td>10</td><td>-750</td><td>0</td><td>18</td><td>59</td></tr><tr><td>Zakynthos<sup>0-1</sup></td><td>3-19/4/2006</td><td>20.921</td><td>37.683</td><td>13</td><td>10</td><td>7</td><td>-290</td><td>-385</td><td>3</td><td>25</td></tr><tr><td>Methoni<sup>0-3</sup></td><td>14/02/2008</td><td>21.700</td><td>36.089</td><td>24.5</td><td>30</td><td>23</td><td>-631</td><td>-901</td><td>310</td><td>8</td></tr><tr><td>Methoni<sup>0-3</sup></td><td>Coseismic + 10 d</td><td>21.705</td><td>36.089</td><td>24.5</td><td>34</td><td>26</td><td>-631</td><td>-901</td><td>310</td><td>8</td></tr><tr><td>Methoni<sup>0-3</sup></td><td>Post-seismic 449 d</td><td>21.189</td><td>36.271</td><td>22.5</td><td>130</td><td>15</td><td>-631</td><td>-901</td><td>310</td><td>8</td></tr><tr><td>Movri<sup>4</sup></td><td>08/06/2008</td><td>21.523</td><td>37.937</td><td>17.5</td><td>15</td><td>8.7</td><td>-1200</td><td>0</td><td>25</td><td>90</td></tr><tr><td>Efpalion<sup>5</sup></td><td>18-22/01/2010</td><td>21.954</td><td>38.416</td><td>4.1</td><td>8.8</td><td>5</td><td>0</td><td>182</td><td>-75</td><td>38</td></tr><tr><td>Western Crete<sup>0</sup></td><td>12/10/2013</td><td>23.320</td><td>35.275</td><td>20</td><td>18</td><td>12</td><td>-840</td><td>-1000</td><td>325</td><td>12</td></tr><tr><td>Cephalonia<sup>6</sup></td><td>26/01/2014</td><td>20.406</td><td>38.180</td><td>2.5</td><td>16</td><td>6.5</td><td>-1160</td><td>0</td><td>14</td><td>70</td></tr><tr><td>Cephalonia<sup>6</sup></td><td>03/02/2014</td><td>20.394</td><td>38.243</td><td>2</td><td>15</td><td>4</td><td>-300</td><td>-300</td><td>14</td><td>30</td></tr><tr><td>Samothraki<sup>0</sup></td><td>24/05/2014</td><td>25.683</td><td>40.283</td><td>0.8</td><td>101</td><td>8.5</td><td>-900</td><td>0</td><td>66</td><td>90</td></tr><tr><td>Samothraki<sup>0</sup></td><td>Coseismic + 1295 d</td><td>25.631</td><td>40.268</td><td>0.8</td><td>110</td><td>13</td><td>-900</td><td>0</td><td>66</td><td>90</td></tr><tr><td>Lefkada<sup>7</sup></td><td>17/11/2015</td><td>20.518</td><td>38.595</td><td>0.15</td><td>23</td><td>10</td><td>-1200</td><td>0</td><td>18</td><td>71</td></tr><tr><td>Lefkada<sup>0</sup></td><td>Post-seismic 138 d</td><td>20.529</td><td>38.598</td><td>2.0</td><td>25</td><td>12</td><td>-100</td><td>0</td><td>18</td><td>71</td></tr><tr><td>Lefkada<sup>0</sup></td><td>Post-seismic 471 d</td><td>20.643</td><td>38.758</td><td>13.5</td><td>51</td><td>22</td><td>-29</td><td>-21</td><td>18</td><td>15</td></tr><tr><td>Gulpinar<sup>8</sup></td><td>6-12/02/2017</td><td>26.136</td><td>39.549</td><td>3.3</td><td>6</td><td>6</td><td>0</td><td>280</td><td>110</td><td>45</td></tr><tr><td>Lesvos<sup>0</sup></td><td>12/06/2017 + 4 d</td><td>26.325</td><td>38.910</td><td>4</td><td>16.1</td><td>6.7</td><td>98</td><td>1250</td><td>111</td><td>45</td></tr><tr><td>Lesvos<sup>0</sup></td><td>Coseismic + 140 d</td><td>26.332</td><td>38.915</td><td>2.8</td><td>16.8</td><td>7.8</td><td>98</td><td>1250</td><td>111</td><td>45</td></tr><tr><td>Kos<sup>9</sup></td><td>20/07/2017</td><td>27.471</td><td>36.874</td><td>2.5</td><td>14</td><td>12.5</td><td>525</td><td>1961</td><td>283</td><td>37</td></tr><tr><td>Zakynthos<sup>10</sup></td><td>25/10/2018</td><td>20.627</td><td>37.415</td><td>3</td><td>26</td><td>16.5</td><td>-1900</td><td>-350</td><td>9</td><td>25</td></tr><tr><td>Zakynthos<sup>10</sup></td><td>Post-seismic 21 d</td><td>20.627</td><td>37.415</td><td>0</td><td>26</td><td>3</td><td>-1900</td><td>0</td><td>9</td><td>90</td></tr><tr><td>Zakynthos<sup>0</sup></td><td>Post-seismic 183 d</td><td>20.876</td><td>37.353</td><td>12</td><td>26</td><td>78</td><td>0</td><td>-90</td><td>310</td><td>12</td></tr><tr><td>Preveza<sup>0</sup></td><td>05/02/2019</td><td>20.608</td><td>38.929</td><td>9.2</td><td>4.2</td><td>2.7</td><td>-467</td><td>0</td><td>2</td><td>81</td></tr><tr><td>Eratini<sup>0</sup></td><td>30/03/2019</td><td>22.269</td><td>38.283</td><td>9</td><td>4</td><td>3</td><td>0</td><td>363</td><td>291</td><td>39</td></tr><tr><td>Magoula<sup>0</sup></td><td>19/07/2019</td><td>23.544</td><td>38.127</td><td>5</td><td>3.8</td><td>2.4</td><td>0</td><td>500</td><td>103</td><td>50</td></tr><tr><td>Kanallaki<sup>11</sup></td><td>21/03/2020</td><td>20.592</td><td>39.247</td><td>7</td><td>4.7</td><td>4.7</td><td>-87</td><td>-492</td><td>327</td><td>39</td></tr><tr><td>Transient of May 2018<sup>0</sup></td><td>10/5/2018</td><td>21.408</td><td>37.595</td><td>22.</td><td>24</td><td>13</td><td>-90</td><td>0</td><td>310</td><td>12</td></tr></tbody></table>

The index in the first column is the source of the information: ⁰ this study, ¹ Müller (2011), ² Ilieva et al. (2016), ³ Howell et al. (2017a), ⁴ Serpetsidaki et al. (2014), ⁵ Elias (2013), ⁶ Briole et al. (2015), ⁷ Ganas et al. (2016b), ⁸ Ganas et al. (2018), ⁹ Ganas et al. (2019), ¹⁰ Ganas et al. (2020), ¹¹ Valkaniotis et al. (2020). A value preceded by a zero means that we calculated a model starting initially from the indicated reference. The coordinates (longitude, latitude and depth) are those of the centre of the top edge of the rectangular fault. Positive values of slip mean left lateral and normal, negative values right lateral and reverse, respectively. The zero reference for the fault azimuth is the north. The models were made with assuming the value of 3 × 10¹⁰ Pa for the rigidity of the medium.

<figure><img src="image_7.png" /><figcaption>**Figure 5.** Percentage of the seismic moment released during the post-seismic phase for nine events of magnitude larger than 6. The parameters used for the earthquakes are presented in Tables 2 and 3.</figcaption></figure>
<table><caption>**Table 3.** Contribution of the post-seismic signal to the total seismic moment for the major events of magnitude larger than 6. The largest events are those that also trigger the post-seismic signals with major relative weight. For several events we could not separate accurately the coseismic and the fast post-seismic, so we give only one value of the total moment of the two phases.</caption><thead><tr><th>Earthquake</th><th>M<sub>w</sub></th><th>Coseismic moment (× 10<sup>16</sup> N m)</th><th>Fast post-seismic moment (× 10<sup>16</sup> N m)</th><th>Time constant (d)</th><th>Slow post-seismic moment (× 10<sup>16</sup> N m)</th><th>Time constant (d)</th><th>Total moment (× 10<sup>16</sup> N m)</th><th>Weight of the post-seismic</th></tr></thead><tbody><tr><td>Skyros 2001</td><td>6.4</td><td>568</td><td></td><td></td><td></td><td></td><td></td><td>Unknown</td></tr><tr><td>Lefkada 2003</td><td>6.2</td><td>360</td><td></td><td></td><td>No evidences</td><td></td><td></td><td>&lt;20 per cent</td></tr><tr><td>Methoni 2008</td><td>6.9</td><td>2277</td><td>640</td><td>10</td><td>6435</td><td>449</td><td>9352</td><td>~75 per cent</td></tr><tr><td>Movri 2008</td><td>6.4</td><td>470</td><td></td><td>103</td><td>Small, seen at RLSO only</td><td></td><td></td><td>~10 per cent</td></tr><tr><td>Cephalonia 2014 (two events)</td><td>equiv 6.2</td><td>438</td><td>minor</td><td></td><td>Not seen in GPS time-series</td><td></td><td></td><td>&lt;10 per cent</td></tr><tr><td>Samothraki 2014</td><td>6.9</td><td>2318</td><td></td><td>65</td><td>1543</td><td>1295</td><td>3861</td><td>~45 per cent</td></tr><tr><td>Lefkada 2015</td><td>6.5</td><td>828</td><td>90</td><td>138</td><td>121</td><td>471</td><td>1039</td><td>20 per cent</td></tr><tr><td>Lesvos 2017</td><td>6.3</td><td>406</td><td></td><td>4</td><td>87</td><td>140</td><td>493</td><td>~25 per cent</td></tr><tr><td>Kos 2017</td><td>6.6</td><td>1066</td><td></td><td>Not seen</td><td>Data not available</td><td>220</td><td></td><td>&lt;25 per cent?</td></tr><tr><td>Zakynthos 2018</td><td>6.8</td><td>2486</td><td>445</td><td>21</td><td>548</td><td>183</td><td>3479</td><td>~30 per cent</td></tr><tr><td>Samos 2020</td><td>7.0</td><td>3200</td><td>Not estimated</td><td></td><td></td><td></td><td></td><td>Not estimated</td></tr><tr><td>Total</td><td></td><td>14417</td><td>1175</td><td></td><td>8734</td><td></td><td></td><td></td></tr></tbody></table>

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Evolution of the seismic moment released by the events listed in Table 3</figcaption></figure>

low dip angle and not coplanar with the coseismic plane yet adjacent to it. In all studied cases, the azimuth of the post-seismic slip vector remains little different from the azimuth of the coseismic slip.

# 4 ANALYSIS OF THE VELOCITY FIELD

The secular velocity field represented with respect to stable Europe is plotted in Fig. 7. The three profiles in red pass through the Euler pole of rotation of Anatolia with respect to Europe (AT-EU) as defined and located by Le Pichon & Kreemer (2010) at 31.96°E, 32.02°N. The radial and tangential vector amplitudes along those profiles are discussed in the Supporting Information. Those profiles are useful in particular because they display well the gradients of velocities.

## 4.1 North-south extension in central Macedonia

Fig. 7 shows that there is a gradient of north-south velocity between the stations located in central Macedonia and those located in Thrace. The latter located in Bulgaria and east of the Strymon river (LOVE, SOFI, ..., KAVA, THAS) present an average rate of

extension of $9 \times 10^{-9}$ strain yr$^{-1}$. The former (TEIS, PROV, ..., MARM, AFY0) present an average rate of extension of $42 \times 10^{-9}$ strain yr$^{-1}$; thus, more than four times larger. This extension is accompanied by the occurrence of large earthquakes, such as, during the 20th century, the Ierissos 1932 (Pavlides & Tranos 1991) and the Thessaloniki 1978 (Stiros & Drakos 2000) events.

## 4.2 North-south extension in the eastern Aegean

Fig. 8 shows the vectors obtained after removing a counterclockwise rotation of $1.56^\circ$ Ma$^{-1}$ (e.g. $27.3 \times 10^{-9}$ rad yr$^{-1}$) around a pole located at 31.96°E, 32.02°N (Fig. 7). In this representation, that cancels the rotation of Anatolia, the vectors along the coast of Turkey are oriented almost north-south with increasing length from north to south. The major step of velocity is located between the islands of Lesvos and Chios (see also Supporting Information Fig. S27) between which we measure a north-south extension of 8 mm yr$^{-1}$. This is the area of the $M_w = 6.3$ earthquake of 2017 June 12 (Papadimitriou et al. 2018) and the large $M_w = 6.8$ earthquake of 1949 July 23 (Melis et al. 2020). Further south there is another step between the stations KALY and KALU (Kalymnos) and the stations
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Velocities of the 309 selected stations, corrected of the effect of the earthquakes (Fig. 4, Supporting Information Table S1), plotted with respect to stable Europe. The three profiles along which the velocities are plotted in Supporting Information Figs S25, S27 and S29 are shown in red. The main segment of the North Anatolian fault is plotted from the west of Lemnos to the east of the Marmara Sea. The labelled stations are those of the profile of Supporting Information Fig. S25. The pole of rotation Anatolia–Eurasia, located at the 31.96°E, 32.02°N, is the one proposed by Le Pichon & Kreemer (2010).</figcaption></figure>

DATC and TILO (Tilos island). We interpret the velocity step between Kalymnos and Tilos as corresponding to extension across the east-west Gökova Gulf and its westward continuation south of Kos (Ganas *et al.* 2019).

## 4.3 The wide stable area of central Aegean, eastern Peloponnese and western Crete

Fig. 9 is a plot of the velocities shifted by 15.6 and -25.2 mm yr⁻¹ with respect to the velocity of stable Europe in east and north respectively. It shows a vast area with little deformation in the central Aegean, stretching from the central and southeastern Peloponnese in the west to the islands of Ikaria, Leros and Kalymnos in the east.

To the north, the transition is abrupt at the south shore of the Gulf of Corinth and extends to the east from the region of Thiva to the south of Evia, as suggested previously by Briole *et al.* (2000), Avallone *et al.* (2004) and Chousianitis *et al.* (2013).

The western part of Crete can be included to this stable domain. Indeed, within our bars of uncertainties, there is minor deformation between the Lefka mountains of Crete, west of VAM0 (8.2 ± 0.3 and -11.3 ± 0.3 mm yr⁻¹ in east and north) and the stable central

Aegean and eastern Peloponnese region (7.9 ± 0.1 and -11.7 ± 0.1 mm yr⁻¹).

## 4.4 The southeastern Aegean

The domain including the eastern part of Crete, Karpathos, Rodos, Astypalea and Tilos constitutes a region where distributed deformation occurs. As discussed in the previous section, to the west this domain starts in central Crete where there is a velocity step of ~1.5 mm yr⁻¹ between the Lefka and the Psiloritis mountains (Fig. 9 and Supporting Information Fig. S31). There is another step of ~1.5 mm yr⁻¹ visible mostly in the east velocity that concerns the most eastern stations, SITI (Sitia) and ZKRO (Zakros). This step is consistent with the location and mechanism of deformation of the northeast/southwest Ierapetra fault zone and the Ptolemy trench (Caputo *et al.* 2010). Our velocity field for Crete is, on the whole, in line with the one of Saltogianni *et al.* (2020).

If we now consider Crete as a single body, the ITRF2014 velocities of the centre of the island are 8.7 ± 0.15 and -12.4 ± 0.15 mm yr⁻¹ in east and north, respectively. With respect to the stable
<figure><img src="image_10.png" /><figcaption>**Figure 8.** Residual velocities after removing a rigid rotation of the Anatolian plate assuming the pole of rotation at 31.96°E, 32.02°N and a rate of counterclockwise rotation of 1.56° Ma⁻¹ (e.g. 27.3 × 10⁻⁹ rad yr⁻¹).</figcaption></figure>

central Aegean region this represented a velocity towards southeast of 0.8 ± 0.2 and -0.7 ± 0.2 mm yr⁻¹ in east and north, respectively. Therefore, the whole island can be considered as rotating clockwise at a rate of 11 × 10⁻⁹ rad yr⁻¹ (i.e. 0.6° Ma⁻¹). This overall rotation is consistent with the azimuth and amplitude of the vectors observed at Karpathos further east. The overall rate of east-west extension of Crete, when considered as a single body, is 5.9 × 10⁻⁹ strain yr⁻¹ (see the slope of the east velocity in Supporting Information Fig. S31) which means that the net extension from the east to the west coast is around 1.5 mm yr⁻¹.

To the north, this southeastern Aegean domain is bounded by a gradient of velocity of around 4 mm yr⁻¹ between Astypalea and Naxos (Supporting Information Fig. S29) with the deformation more likely to be localized in the NE-SW Amorgos basin, south of the eponymous island (Nomikou et al. 2018).

To the northeast, this domain terminates south of the island of Kos and in the Gökova gulf that is known to be active with a localized and relatively fast extension of ~4 mm yr⁻¹ (Ganas et al. 2019) and where the M<sub>w</sub> = 6.6 Kos earthquake occurred in 2017.

## 4.5 East-west extension in the western Peloponnese

In the western Peloponnese, Fig. 9 shows a significant extension that starts west of the line joining Kalavrita, Tripoli and Monemvasia. This extension, almost east-west, has been described by Lyon-Caen et al. (1988) after the M<sub>w</sub> = 5.9 1986 September 13 Kalamata earthquake [see also Stiros & Kontogianni (2008) for geodetic observations of that earthquake] and by Armijo et al. (1992) who studied the historical Sparta earthquake of 464 B.C., which they estimated around M = 7.2. It was also analysed using
<figure><img src="image_11.png" /><figcaption>**Figure 9.** Velocities shifted by 15.6 and -25.2 mm yr⁻¹ with respect to the velocity of the stable Europe in east and north respectively, highlighting the stability of a region constituted of the central Aegean, eastern Peloponnese and western Crete.</figcaption></figure>

GPS data, in particular by Floyd *et al.* (2010) and Chousianitis *et al.* (2015). The extension captured by the GPS stations can be analysed in Supporting Information Fig. S30. It concerns mostly the domain located between Sparta and the western coast (Pylos) with ~4 mm yr⁻¹ east-west extension localized there in a band of ~50 km. To the east of Sparta there is also extension, but it is much smaller, with amplitude ~2 mm yr⁻¹ occurring in a wider band of ~100 km.

## 4.6 Rotation of the Pindos region and collision with the Apulian platform

Fig. 10 shows the vectors of northern Greece after removing a clockwise rotation of 2.6° Ma⁻¹ (see Fig. 7) with respect to stable Europe. This rotation minimizes almost completely the vectors in an extended area starting in central Greece around the latitude 38.5° and bounded to the east by the line formed by the Olympus mountain range. East of this range Fig. 10 shows that the rotation rapidly fails to minimize the vectors. To the west the rotation does not minimize the vectors along the western coast of Epirus (e.g. station GARD) and in the islands of Paxoi and Corfu. This is expected as we know that there is a significant and localized compression in the mountain ranges along the Epirus coast (e.g. Baker *et al.* 1997; Ganas *et al.* 2013; D'Agostino *et al.* 2020) with the frequent occurrence of significant earthquakes, for example, the Kanallaki earthquake of 2020 March 21 (Table 1).

Using our GPS velocities, Valkaniotis *et al.* (2020) have established that, at the latitude of Kanallaki and Paxoi, the convergence between the Apulian platform (using the velocity of D'Agostino *et al.* 2020) and the central Pindos region (using our velocity) is

8.9 mm yr⁻¹ in the N228°E azimuth. This azimuth is significantly different than the N247°E obtained by Louvari *et al.* (2001) by analysing the slip vectors of the earthquakes in northwestern Greece. Concerning the rate of convergence, we find an extra 3.7 mm yr⁻¹ with respect to the value of the convergence between Apulian and stable Europe at that location (5.2 mm yr⁻¹). This extra convergence is induced by the compression imposed by the motion of the Anatolian micro plate, itself resulting in the rotation of the Pindos region (Valkaniotis *et al.* 2020). In consequence, we expect a larger rate of convergence to the south of this zone, around Preveza and a lower rate of convergence to the north around Corfu and the continuation of this decrease until Vlora. Valkaniotis *et al.* (2020) evaluated at 60 km the width of the deforming zone (defined as the zone where 80 per cent of the total convergence is accommodated), a value consistent with that found by Kassaras *et al.* (2016) from their analysis of a large database of focal mechanisms (see top section in their fig. 11).

The thickness of the brittle (locked/seismic) area found by Valkaniotis *et al.* (2020) is 11 km, a value consistent with the depth of the fault of the Kanallaki earthquake. They showed also that the central axis of the plate boundary runs along the southwestern coast of Corfu and the northeastern coast of Paxoi, which means that half of the strain is accommodated offshore in a band of ~30 km.

## 4.7 Distributed deformation in central Greece

This is an east-west shaped area that roughly coincides with the azimuth projected inland of the major branch of the North Anatolian fault (Hatzfeld *et al.* 1999). It is a region with distributed extension and the largest strain rates observed in Europe, especially
<figure><img src="image_12.png" /><figcaption>**Figure 10.** Velocities in and around the Pindos region plotted after counter-clockwise rotation of 2.6° Ma⁻¹ of the vectors of Fig. 7 around a pole of rotation Pindos-Eurasia located 10 km west of Permet (Albania) represented by the large red circle. The small red circle is the location of the pole obtained when using also the stations of central Greece located south of the latitude 38.5° and north of the Gulf of Corinth.</figcaption></figure>

in the western Gulf of Corinth (Clarke *et al.* 1998; Avallone *et al.* 2004) where strong historical seismicity is recorded (Ambraseys & Jackson, 1990, 1997) and large faults well exposed in the field, for example, Roberts & Ganas (2000) and Goldsworthy & Jackson (2001). There are several grabens that localize most of the deformation:

(1) the graben of the Ambracian Gulf, south of Arta where we measure an extension larger than 10 mm yr⁻¹ between the stations of ARTA and ARTU and the stations AMFI and ABEL;

(2) the graben of the Trichonis lake with an oblique extension of ~2 mm yr⁻¹ measured between the station of AGRI, AGRU, PVOG and the stations MESO and RETS, and where moderate seismicity is recorded with the occurrence of events of swarms like the one of 2007 (Kiratzi *et al.* 2008);

(3) the Gulf of Patras with an extension rate of ~7 mm yr⁻¹ measured between RETS and MESO to the north and MAOR and ANOA to the south;

(4) the rift of Corinth, extensively studied for decades (e.g. Ambraseys & Jackson 1997) where the deformation reaches its acme in the western part of the rift west of the island of Trizonia, with an intense level of microseismicity recorded almost at all times;

(5) the graben of Evia with major normal faults north-dipping (Ganas & Papoulia 2000) and showing less seismicity during the 20th century with respect to the Gulf of Corinth.

The whole area of central Greece is rotating with respect to the central Aegean domain. Here, taking the stations located south of the latitude of ~38.5° (defined in the previous section as the approximate limit of the Pindos domain) and the northern shore of the Gulf of Corinth, we find a clockwise rotation of ~4° Ma⁻¹ for this domain, thus more than the 2.6° Ma⁻¹ of the Pindos domain. This value is consistent with the determination of Chousianitis *et al.* 2015. The larger value of 7 ± 0.5° Ma⁻¹ found by Avallone *et al.* (2004, see their fig. 3) is due to the fact that they use mostly stations located very close to the Gulf of Corinth and therefore affected by the extension (increasing from east to west) along the rift.

## 4.8 The South Ionian Islands and surroundings

This south Ionian domain, at the northwestern termination of the Hellenic megathrust, is an area of high seismicity and strong internal deformation (Ganas *et al.* 2020; Haddad *et al.* 2020). It is bounded
by the Cephalonia fault zone to the west and by the major strike-slip faults of Movri and Katouna to the east. The deformation of this relatively small area (Kahle et al. 1993; Louvari et al. 1999) is complex as a significant part of the deformation corresponds to the accumulation of strain in the vicinity of the major right lateral Cephalonia fault zone where the recurrence time of earthquakes is short (Briole et al. 2015).

The left lateral strike-slip motion along the Katouna fault is larger than 6 mm yr⁻¹ according to the GPS as inferred by the relative velocities of KTCH (and to some extent RGNI, yet it should be localized in the deforming zone of the fault) with respect to MESO and RETS.

In our velocity field, we can see also a difference in the behaviour of the vectors when moving from north to south. While compression is visible in Lefkada, Cephalonia and Zakynthos, it disappears further south at Strofades. This suggests the existence of two different domains. To the north, the plate interface is dominantly locked and produce earthquakes along the subduction megathrust, and along the Cephalonia fault zone. The GPS velocities at the southern stations of Cephalonia KIPO, KEFA, ARGO, VLSM and SISS, which are roughly aligned perpendicular to the Cephalonia fault, show well a gradient of velocities characteristic of strain accumulation during the interseismic period. The observed gradient is consistent with the one evaluated by Briole et al. (2015) using the velocities of campaign stations located in the south of Cephalonia. To the south of Zakynthos, in particular at Strofades (STRF), the convergence appears almost unlocked as it is also observed in the entire offshore southern Peloponnese and south of Crete (see Section 4.4), with a population of subduction earthquakes small compared to what it would be in the case of a fully coupled zone. This strong unlocking has been quantified and discuss by various authors since the early works of the 1990s, in particular Ekström & England (1989).

## 4.9 The vertical velocities

Fig. 11 shows the vertical velocities at a selection of 177 points having an uncertainty on that velocity lower than 1 mm yr⁻¹. We do not see clustering in the distribution of the velocities with areas showing preferential up or down movement within the bars of uncertainties. Supporting Information Fig. S33 shows the cumulative distribution and best Gaussian fit calculated with the whole set of 322 stations and with the subset of 177 stations having velocity uncertainty lower than 1 mm yr⁻¹. The best Gaussian fit is calculated for the velocities comprised between −2.7 and 1.4 mm yr⁻¹ only, which corresponds to 90 per cent of the population of stations, the data at the edge of the distribution (6 per cent on the negative velocity side and 4 per cent on the positive velocity side) showing clear departure from a Gaussian distribution. Most of the outliers of the negative velocity side are likely to correspond to local subsidence that can have various causes, natural or anthropogenic, for example, GEYB that is mentioned in Supporting Information Section S4.7. The median of the 177 best velocities is −0.38 mm yr⁻¹. This overall subsidence of the whole area is slightly less but relatively consistent with the observations and models of Serpelloni et al. (2013) who find (e.g. in their fig. 8) a subsidence of Greece and western Turkey in the range −0.4 to −0.8 mm yr⁻¹, increasing towards southwest, with average −0.6 mm yr⁻¹. Our geodetic results are not sufficient to investigate local aspects of the vertical deformations in geological timescales, and therefore we do not compare our velocities with the velocities predicted by Howell et al. (2017b) or with the observations of Robertson et al. (2019).

# 5 MODEL WITH TEN CRUSTAL BLOCKS

Using the analysis made in the previous section, we will now build and assess a kinematic model composed of ten crustal blocks. This model is intended to reproduce at the first order the velocity field. The applications of a deformation model can concern geophysical considerations but also various practical needs of surveyors working in Greece. To set up the model requires (i) deciding the number of blocks, (ii) deciding the boundaries between the blocks, (iii) deciding the parameters that will characterize each block, (iv) evaluating those parameters and (v) assessing the performance of the model. Let us also point out that there was already such a block model in the article of Nyst & Thatcher (2004), yet coarser because of the limited available GPS velocities at that time, but it shares already several conceptual characteristics with ours.

## 5.1 Definition of the crustal blocks

In the previous sections, we have seen that there is a relatively well defined and large stable region located in the central Aegean Sea and the eastern Peloponnese. To the southeast the islands of Karpathos, Rodos and Astypalea belong to a different structure that is moving towards the southeast with respect to the previous one and rotating as we have seen for southeast Crete and even more in Karpathos. Those two structures will constitute our blocks CA (Central Aegean) and SA (Southern Aegean). We have seen that there is a large east-west extension in the western Peloponnese mostly located west of the Taygetus mountain but with some minor expression even east of it. Our third block, called WP (West Peloponnese), will contain this east-west extension as well as most of the variable deformation related to the subduction in the West Hellenic arc.

To the northeast of the CA block, we define a block that we call NA (Northern Aegean) that contains all the islands located between the north of Samos and Lemnos, including Skyros and Skopelos to the west. We have seen in Section 4.2 that this domain is deforming internally. However, in our model we consider it as a single block despite the strong internal deformation between Lesvos and Chios. The reason is that we do not have enough GPS stations and long enough time-series to split the zone in two and retrieve robust parameters for both.

North and northwest of the North Anatolian fault we have seen that there are three domains with relatively distinct velocity fields: the largest one is the one that we will call PI (Pindos) and that we used to estimate the rotation of this region. East of it we have seen that central Macedonia, east of the line of the Olympus mountain present a different kinematics with minor rotation and large north-south extension. Further east we have seen that this large extension ends at the longitude of the elbow of the North Anatolian fault which corresponds to the region of Thrace, east of the Strymon river. Therefore, this defines two other blocks that we call CM (Central Macedonia) and TH (Thrace). To the west of the Pindos block we have seen that, despite the limited number of stations, the GPS vectors show well the compression in the mountain ranges located along the coast in the west of Epirus. We define a block there that we call WE (Western Epirus).

Between the Pindos area and the Peloponnese there is the area of central Greece, north of the Gulf of Corinth and including the centre and north of Evia. We define a block there that we call CG (Central Greece). The limit between the PI and CG blocks is not well defined because of the limited number of GPS points and the distributed nature of the deformation. However Fig. 10 indicates that KARP (Karpenissi), LAMI and LAMA (Lamia) belong to the
<figure><img src="image_13.png" /><figcaption>**Figure 11.** Vertical velocities of the 177 points with uncertainty on this velocity lower than 1 mm yr⁻¹. Black and white dots mean that the velocity is below -2 mm yr⁻¹ or above 2 mm yr⁻¹, respectively. The size of the dots is inversely proportional to the uncertainty. The 68 points where the uncertainty is below 0.7 mm yr⁻¹ have their code plotted.</figcaption></figure>

former while ANOC and LIDO belong to the latter, so we draw the boundary from the Trichonis Lake where we know that there is active faulting and north-south extension to the northern entrance of the Gulf of Evia.

Our last block, which fits with the IAB (Ionian Islands-Akarnania Block) of Pérouse et al. (2017), contains the South Ionian islands (Zakynthos, Cephalonia, Lefkada) and the regions of north-west Peloponnese (around Killini) and the western part of Aetolia-Acarnania, west of the Katouna fault, we call it the SI block.

Supporting Information Fig. S34 presents the location of the ten blocks together with the focal mechanisms (those for depths shallower than 35 km) published by the National and Kapodistrian University of Athens for the period 2003–2018.

The major fault zones and zones of localized deformation of Greece constitute boundaries between the blocks, this is the case in particular of the North Anatolian fault (segments i and j), the Corinth rift and its extension through the south of Evia (segments c and m), the Amorgos basin fault zone (segment k), the Movri-Patras fault (segment a), the Katouna fault (segment e), the narrow and deforming mountain ranges of WE (segment h). The Cephalonia fault (segment l) is also an obvious boundary with the Apulian platform. To the north, the boundary with stable Europe is arbitrary but can be used to verify that our model is not implying a gap of velocity there which can be verified in Fig. 13. To the east the boundary inside Turkey is also arbitrary. To the south and west the boundaries of the model are at sea.

## 5.2 Boundaries of the blocks

The boundaries of the blocks are plotted in Fig. 12. Details of the criteria and location of 21 of the boundaries segments (coded using the letters a–u) are given in Table 4. For practical reasons, to simplify further analysis, we minimized the number of edges and maximize the boundaries located at sea (e.g. at the north entrance of the Evia gulf, south of Pilion) where the boundary matters for seismology and geology but cannot be well constrained by GPS.

## 5.3 Parameters of the blocks

For each block, we define seven parameters. Two are the coordinates of the barycentre of the GPS stations of the block. Two are the horizontal components of the estimated velocity of this barycentre calculated using the average of the velocities of the various GPS stations of the block without weighting by their uncertainties. We found insignificant differences (below 0.1 mm yr⁻¹) when weighting each velocity by its uncertainty. The fifth parameter is the rotation of the
<figure><img src="image_14.png" /><figcaption>**Figure 12.** Boundaries of the crustal blocks model HelVel2020. The letters in yellow (a–u) correspond to boundaries that are discussed in Table 4 and in the text. The residual velocities after applying the blocks model are plotted (compare with the original vectors plotted at the same scale in Figs 7 and 8). Among the 265 stations used for determining the parameters of the model, 119 with post-fit residual lower than 1 mm yr⁻¹ are not plotted, 80 with post-fit residual between 1 and 2 mm yr⁻¹ are plotted with green arrow, 32 with post-fit residual between 2 and 3 mm yr⁻¹ are plotted with purple arrow and 34 with post-fit residual larger than 3 mm yr⁻¹ are plotted with black arrows.</figcaption></figure>

block determined for all block in the same way we did for Pindos in the previous section. The last two parameters are a homogeneous (linear) component of strain in the east and north components estimated by minimization of the residuals obtained after the rotation. The parameters of the ten blocks are listed in Table 5. In two of the blocks not all stations were used for the minimization. This is the case of CA where nine stations located in a radius of 10 km around the city of Patras were removed because of the large discrepancy of those stations that are located right in the boundary between two or more blocks (Elias & Briole 2018; Papadopoulos et al. 2019). This is also the case of TH where six stations of Bulgaria were included for further analysis but not used for the minimization. The total number of stations used for the assessment of the model are therefore 249 (for a total of 265 included in the blocks analysis). The number of stations used in each block is reported in Table 5.

The background concept of our crustal blocks model is close to the concept of blocks model developed by Simpson *et al.* (2012) and refined in particular by Savage & Simpson (2013a,b),

Savage & Well (2015) and Savage (2018). Both are using the concept of Euler-vector clustering, with the same theoretical background, but instead of inverting for the two coordinates of a Euler pole defined with respect to the ITRF2014 and a rotation with respect to this pole, we invert for the velocity of our blocks in the ITRF2014 and an internal rotation of the block (with respect to its barycentre). Both methods are mathematically equivalent, and we believe that ours is more intuitive in terms of understanding of the local relative displacements and differential rotations of the blocks. Concerning the internal deformation of the block, in their formalism Savage & Simpson (2013a) make the assumption that the upper crust is elastic. They use the corresponding equations for uniform strain and rotation rates (their eq. 4), so they invert for four parameters ε<sub>E</sub>, ε<sub>N</sub>, ε<sub>EN</sub> and Ω. However, later in their article (section 39)) they discuss the fact that the upper crust might be inelastic and discuss the impact for the rotation rate Ω. Savage & Simpson (2013b) make the same observation, and in their section [24] indicate that the strain accumulation they calculate need not be inelastic. In the case of the Aegean, there is a distributed faulting
<table><caption>**Table 4.** Definition of the edges of the ten blocks model. The codes are plotted in Fig. 12.</caption><thead><tr><th>Code</th><th>Long.</th><th>Lat.</th><th>Nature of the discontinuity / data supporting placing a block boundary there</th></tr></thead><tbody><tr><td>a</td><td>21.53</td><td>37.95</td><td>Right-lateral strike-slip illustrated by the 2008 Movri earthquake (Serpetsidaki et al. 2014).</td></tr><tr><td>b</td><td>21.56</td><td>38.19</td><td>Extension across the Gulf of Patras involving east-west faults located south of the Gulf according to Flotté et al. (2005).</td></tr><tr><td>c</td><td>22.14</td><td>38.31</td><td>Zone of the fastest extension of the rift of Corinth around Aigion.</td></tr><tr><td>d</td><td>22.80</td><td>38.18</td><td>Eastern rift of Corinth.</td></tr><tr><td>e</td><td>21.19</td><td>38.70</td><td>Left-lateral strike-slip across the Katouna-Stamna fault zone (Pérouse et al. 2017).</td></tr><tr><td>f</td><td>21.49</td><td>38.56</td><td>Extension across the Trichonis lake (see also u).</td></tr><tr><td>g</td><td>20.88</td><td>38.96</td><td>Extension across the Ambracian Gulf.</td></tr><tr><td>h</td><td>20.44</td><td>39.53</td><td>Compression in western Epirus where moderate to strong seismicity occurs episodically, for example, on March 21, 2020 with the M<sub>w</sub> = 5.7 Kanallaki earthquake (Valkaniotis et al. 2020).</td></tr><tr><td>i</td><td>25.28</td><td>40.18</td><td>Main branch of the North Anatolian fault on which the M<sub>w</sub> = 6.9 Samothraki 2004 earthquake occurred.</td></tr><tr><td>j</td><td>24.02</td><td>39.58</td><td>Main branch of the North Anatolian fault west of the elbow located southeast of the Athos mount range.</td></tr><tr><td>k</td><td>25.88</td><td>26.67</td><td>Extension across the Amorgos basin.</td></tr><tr><td>l</td><td>20.40</td><td>38.71</td><td>Major right lateral fault of Lefkada-Cephalonia.</td></tr><tr><td>m</td><td>23.76</td><td>38.32</td><td>Crosses through southern Evia to reach the eastern termination of the Gulf of Corinth in the region of the Kaparelli fault (1981 earthquake) and other faults.</td></tr><tr><td>n</td><td>27.55</td><td>36.83</td><td>Extension across the Gökova gulf (where the Kos, 2017 earthquake occurred) and further west and south of Kos.</td></tr><tr><td>o</td><td>24.25</td><td>35.26</td><td>Location of the offset of velocity between western Crete (which shows almost no motion with respect to central Greece) and eastern Crete (see Supporting Information Fig. S31).</td></tr><tr><td>p</td><td>22.89</td><td>36.94</td><td>Diffuse zone where the east-west extension in the Peloponnese becomes significant. The western part of the Peloponnese is characterized by distributed extension.</td></tr><tr><td>q</td><td>22.31</td><td>40.34</td><td>Discontinuity observed in the GPS data along the Olympus mountain range and further north towards Naoussa and Florina.</td></tr><tr><td>r</td><td>23.82</td><td>41.11</td><td>Discontinuity observed in the GPS velocities following approximately the course of the Strymon valley.</td></tr><tr><td>s</td><td>24.66</td><td>38.29</td><td>Eastern coast of Evia, Andros and Tinos.</td></tr><tr><td>t</td><td>26.47</td><td>37.75</td><td>Extension across a major fault zone located along the northern coasts of Ikaria and Samos. The residual velocity at the station SAMU (east of Samos) might correspond to strain accumulation in this deforming fault zone. Indeed, SAMU is located in the footwall of the 2020 fault, close to it, thus within the region sensitive to strain accumulation and release during the seismic cycle.</td></tr><tr><td>u</td><td>22.50</td><td>38.81</td><td>Through the northwestern entrance of the Gulf of Evia and south of the Pindos mountains towards the triple junction of Stamna west of Agrinion.</td></tr></tbody></table>

<table><caption>**Table 5.** Characteristics of the ten blocks of our velocity model HelVel2020. The codes of the domains are the following: SA, Southern Aegean; CA, Central Aegean; WP, Western Peloponnese; SI, Southern Ionian Islands; WE, Western Epirus; PI, Pindos; CM, Central Macedonia; TH, Thrace; NA, Northern Aegean; CG, Central Greece. The file with the boundaries of the 10 domains is in the Supporting Information. The longitude and latitude are those of the barycentre of the GPS stations used to characterize the domain. The velocities are the average of the ITRF2014 velocities of the stations considered for each domain (without weighting with the single uncertainty of each station). The uncertainties σ<sub>E</sub> and σ<sub>N</sub> are the r.m.s. deviations of the residual vectors in the block after rotation and corrected for the strain, divided by the square root of the number of stations. <sup>a</sup>There are nine extra stations in this block, located around Patras, not used for the modelling because they are in the deforming zone between the block CA, SI and CG. <sup>b</sup>There are seven extra stations further north in Bulgaria not used for the modelling but present in the final list of residual velocities. The r.m.s. deviations of the residual vectors are listed in the last two columns. The global average scatters for the whole selection of 249 stations (used in the model) is 0.87 and 0.97 mm yr<sup>-1</sup> in east and north, respectively. Fig. 12 shows the map of the residual vectors, corrected from the blocks model, at the 265 selected GPS stations.</caption><thead><tr><th>Block</th><th>Stations used</th><th>Long. (°)</th><th>Lat. (°)</th><th>v<sub>E</sub> (mm yr<sup>-1</sup>)</th><th>v<sub>N</sub> (mm yr<sup>-1</sup>)</th><th>σ<sub>E</sub> (mm yr<sup>-1</sup>)</th><th>σ<sub>N</sub> (mm yr<sup>-1</sup>)</th><th>Ω (° Ma<sup>-1</sup>)</th><th>σ<sub>Ω</sub> (° Ma<sup>-1</sup>)</th><th>ε<sub>E</sub> (nstrain)</th><th>ε<sub>N</sub> (nstrain)</th><th>rms<sub>E</sub> (mm yr<sup>-1</sup>)</th><th>rms<sub>N</sub> (mm yr<sup>-1</sup>)</th></tr></thead><tbody><tr><td>SA</td><td>28</td><td>25.731</td><td>35.538</td><td>10.1</td><td>-13.5</td><td>0.28</td><td>0.21</td><td>0.6</td><td>0.2</td><td>15.1</td><td>4.3</td><td>1.01</td><td>0.99</td></tr><tr><td>CA</td><td>47<sup>a</sup></td><td>23.994</td><td>37.414</td><td>8.1</td><td>-11.5</td><td>0.08</td><td>0.08</td><td>0.2</td><td>0.1</td><td>-0.6</td><td>2.1</td><td>0.69</td><td>0.71</td></tr><tr><td>WP</td><td>21</td><td>22.213</td><td>36.923</td><td>5.0</td><td>-12.0</td><td>0.20</td><td>0.09</td><td>0.1</td><td>0.3</td><td>26.4</td><td>5.6</td><td>0.76</td><td>0.50</td></tr><tr><td>SI</td><td>24</td><td>20.863</td><td>38.263</td><td>15.2</td><td>1.2</td><td>0.21</td><td>0.27</td><td>6.2</td><td>0.8</td><td>-50.5</td><td>28</td><td>1.48</td><td>1.43</td></tr><tr><td>WE</td><td>7</td><td>20.140</td><td>39.485</td><td>21.6</td><td>14.4</td><td>0.23</td><td>0.15</td><td>0.7</td><td>1.1</td><td>7.6</td><td>-32.5</td><td>0.97</td><td>0.49</td></tr><tr><td>PI</td><td>39</td><td>21.874</td><td>39.602</td><td>20.1</td><td>7.1</td><td>0.10</td><td>0.19</td><td>2.6</td><td>0.2</td><td>3.4</td><td>18.2</td><td>0.88</td><td>0.82</td></tr><tr><td>CM</td><td>21</td><td>23.141</td><td>40.649</td><td>23.9</td><td>7.9</td><td>0.10</td><td>0.19</td><td>0.4</td><td>0.3</td><td>-9.0</td><td>36.1</td><td>0.39</td><td>0.72</td></tr><tr><td>TH</td><td>18<sup>b</sup></td><td>25.340</td><td>41.191</td><td>23.2</td><td>12.3</td><td>0.07</td><td>0.10</td><td>-0.4</td><td>0.2</td><td>-3.6</td><td>10.2</td><td>0.41</td><td>0.45</td></tr><tr><td>NA</td><td>18</td><td>25.718</td><td>38.943</td><td>8.2</td><td>-3.6</td><td>0.56</td><td>0.28</td><td>0.5</td><td>0.7</td><td>-33.9</td><td>48.2</td><td>1.09</td><td>2.93</td></tr><tr><td>CG</td><td>26</td><td>22.628</td><td>38.477</td><td>11.7</td><td>-2.4</td><td>0.22</td><td>0.20</td><td>4.0</td><td>0.4</td><td>-16.1</td><td>56.9</td><td>0.83</td><td>0.93</td></tr></tbody></table>

and we do not expect the upper crust to be perfectly elastic. Therefore, in our formalism we calculated three terms only, ε<sub>E</sub>, ε<sub>N</sub> and Ω, thus ignoring ε<sub>EN</sub>. This approximation is also making the discussion simpler for the geological analysis with only two terms of east-west and north-south homogeneous deformations within the blocks.

A major difference is the fact that we have established the boundaries of our blocks based on the analysis of the GPS velocity field,

the geology and fault maps, and the seismicity, without trying to perform a “neutral” clustering analysis like the k-medoids method applied by Savage & Simpson (2013a) or Savage & Simpson (2013b).

## 5.4 Fit to the data

For the 249 stations used for defining the parameters of the blocks we calculated the residuals after applying the deformation of the
<figure><img src="image_15.png" /><figcaption>**Figure 13.** Velocities (vectors) of the ten blocks with respect to Eurasia, and their rates of rotation (see numerical values in Table 5). The areas in grey represent the theoretical location of the blocks after two million years. The gaps and overlaps that exist between some limits of blocks correspond to areas where there is a localized deformation (extension, compression, strike-slip) between the two blocks not taken into account in the blocks model. In particular, the rift of Corinth and its eastern continuation across Evia, are well identified. The areas listed in Table 6 are shown with boxes containing the letters a–m (the same letters as those used in Fig. 12 and Table 5).</figcaption></figure>

block. Fig. 12 shows the residual vectors. The average r.m.s. fit is 0.87 and 0.97 mm yr⁻¹ respectively in east and north.

## 5.5 A programme to calculate the modelled velocities

We wrote a Fortran programme (called helvel2020.f) to calculate the predicted velocity of any point of Greece, along with the block to which it belongs. The programme and its documentation are provided in the supplementary material. The programme can be used to estimate the velocity gradient in several boundaries of blocks, this can be particularly useful where there are no (or only a few) points located close to the analysed boundary.

## 5.6 Discussion of the model

We calculated and mapped in Fig. 13 the motion and internal deformation of the ten blocks considering, arbitrarily, a period of 2 Ma.

We used this duration because it allows at the same time a good visual understanding of Fig. 13, and a correct kinematic analysis, as during Quaternary the main features of active tectonics have not changed much. Several features are visible, in particular, as expected when setting the boundaries, the localized deformation along them. This is the case in particular of the Gulf of Corinth, the North Anatolian fault, the Amorgos basin, the Movri fault zone. We observe also ad hoc deformation in the central Peloponnese, this is because we placed arbitrarily the boundary there while we know from the GPS velocities (Supporting Information Fig. S30) and from seismological analysis (e.g. Kapetanidis & Kassaras 2019) that it is a diffuse boundary.

Table 6 gathers the predictions of the model at various block boundaries. We find shear velocities of 7.4 mm yr⁻¹ and shortening of 2.3 mm yr⁻¹ across the Movri fault zone at the location of the $M_w = 6.4$ 2008 June 8 earthquake in northwest Peloponnese. Given the limited accuracy of the parameters of the SI block (see Table 5) and the amplitude of the residuals, the uncertainty on the component of compression is around ± 2 mm yr⁻¹, therefore we conclude that
<table><caption>**Table 6.** Relative east (dVE) and north (dVN) velocities between blocks at several specific areas of interest. Az, Azimuth of the boundary with respect to north (positive = clockwise); Str, Strike-slip component (positive = left lateral); Ext, Extension component (positive = extension); a components estimated assuming that the fault is oriented east-west in the studied area.</caption><thead><tr><th>Object</th><th>Code in Fig. 13</th><th>Long (°)</th><th>Lat (°)</th><th>dVE (mm yr<sup>-1</sup>)</th><th>dVN (mm yr<sup>-1</sup>)</th><th>Az (°)</th><th>Str (mm yr<sup>-1</sup>)</th><th>Ext (mm yr<sup>-1</sup>)</th></tr></thead><tbody><tr><td>Movri fault near the epicentre of the 2008 earthquake</td><td>a</td><td>21.53</td><td>37.95</td><td>5.2</td><td>5.7</td><td>25</td><td>-7.4</td><td>-2.3</td></tr><tr><td>Extension across the Gulf of Patras</td><td>b</td><td>21.56</td><td>38.19</td><td>-0.1</td><td>7.9</td><td>285</td><td>2.1</td><td>7.6</td></tr><tr><td>Extension of the Gulf of Corinth at the longitude of Aigion</td><td>c</td><td>22.14</td><td>38.31</td><td>2.8</td><td>10.4</td><td>284</td><td>0.2</td><td>10.8</td></tr><tr><td>Extension of the Gulf of Corinth at the longitude of Kiato</td><td>d</td><td>22.83</td><td>38.17</td><td>0.8</td><td>5.4</td><td>284</td><td>-0.7</td><td>5.4</td></tr><tr><td>Strike-slip along the Katouna fault</td><td>e</td><td>21.19</td><td>38.70</td><td>-3.6</td><td>8.3</td><td>336</td><td>9.0</td><td>-0.1</td></tr><tr><td>Extension across the Trichonis lake</td><td>f</td><td>21.49</td><td>38.56</td><td>1.0</td><td>1.5</td><td>90</td><td>1.0<sup>a</sup></td><td>1.5<sup>a</sup></td></tr><tr><td>Extension across the Ambracian Gulf</td><td>g</td><td>20.88</td><td>38.96</td><td>-1.7</td><td>12.6</td><td>90</td><td>-1.7<sup>a</sup></td><td>12.6<sup>a</sup></td></tr><tr><td>Compression across the Paramythia mountains (assuming that they localize the entire compression)</td><td>h</td><td>20.53</td><td>39.45</td><td>2.8</td><td>2.3</td><td>321</td><td>0.0</td><td>-3.7</td></tr><tr><td>North Anatolian fault near Lemnos</td><td>i</td><td>25.02</td><td>40.08</td><td>12.6</td><td>8.2</td><td>65</td><td>-14.9</td><td>2.2</td></tr><tr><td>North Anatolian north of Alonissos (west of the elbow)</td><td>j</td><td>24.02</td><td>39.58</td><td>8.9</td><td>3.1</td><td>48</td><td>-8.7</td><td>-3.7</td></tr><tr><td>Amorgos right lateral basin</td><td>k</td><td>26.12</td><td>36.71</td><td>4.3</td><td>-1.1</td><td>70</td><td>3.7</td><td>2.5</td></tr><tr><td>North end of the Cephalonia fault (north of Lefkada)</td><td>l</td><td>20.60</td><td>38.91</td><td>-2.7</td><td>10.6</td><td>17</td><td>-9.4</td><td>5.7</td></tr><tr><td>Shear component in the transition zone Corinth-Evia</td><td>m</td><td>23.76</td><td>38.32</td><td>0.7</td><td>1.0</td><td>70</td><td>-1.0</td><td>0.7</td></tr></tbody></table>

compression is probable across this fault zone, with an amplitude not well defined in the range ~0.5–4 mm yr⁻¹.

Further north, along the Rio-Patras right lateral fault, we find 12.1 mm yr⁻¹ of right lateral shear and 3.8 mm yr⁻¹ of extension, two values relatively close to those (15 and 6.2 mm yr⁻¹) inferred by Elias & Briole (2018) by modelling the gradients of GPS vectors across the termination of the Corinth rift (see their table 4). To the west of those strike-slip faults, we find an extension of 7.6 mm yr⁻¹ across the Gulf of Patras with a small component of left lateral strike-slip of 2.1 mm yr⁻¹. The relatively low seismicity there, as recorded during the last decades (Haddad et al. 2020), suggests that a significant part of the extension may be accommodated aseismically. Further north, we find a shear rate of 9.0 mm yr⁻¹ on the left-lateral Katouna fault, which is in the upper range of what was predicted by Pérouse et al. (2017) using geological evidence (minimum 4 mm yr⁻¹) and the GPS data available at that time (around 10 mm yr⁻¹). North of this fault, we find a fast extension of 12.6 mm yr⁻¹ with minor right lateral component of 1.7 mm yr⁻¹, across the Ambracian Gulf, an area poorly documented by GPS especially to the north. This fast extension is comparable to that of the western end of the Corinth rift (which is coherent and expected from the kinematics points of view), but it does not produce intense microseismicity like the western Corinth rift does. The compression across WE is estimated at 3.7 mm yr⁻¹. There are not enough GPS stations there, especially in the southern part of this block, between Preveza and Parga along the coast, and along the line Arta-Ioannina in the interior to characterize accurately this deformation. It may occur in a relatively narrow zone of 50 km or less but much more GPS data is needed there to get the required accuracy.

We find 14.9 mm yr⁻¹ of shear with small extension of 2.2 mm yr⁻¹ on the North Anatolian fault north of Lemnos, and 8.7 mm yr⁻¹ with 3.7 mm yr⁻¹ of shortening along the Sporades segment of the fault located between Skopelos and the elbow of the fault west of Lemnos, where the strike of the fault changes from N68°E to N50°E. The extension of the Gulf of Corinth is well fitted by the model, with almost pure extension of 10.8 mm yr⁻¹ at the longitude of Aigion and 5.5 mm yr⁻¹ at the longitude of Kiato, consistent with the velocities given by Elias & Briole (2018), Chousianitis et al. (2013) and Avallone et al. (2004). There is a transtensional movement of 4.5 mm yr⁻¹ between the islands of Amorgos and Astypalea.

Only the Ionian Islands region shows evidence of dominant coupling along the subduction interface. For the Cephalonia fault, if we assume that the central axis of the fault zone is in the middle of the Palliki peninsula, following Briole et al. (2015), and use their velocity of Apulia in front of Cephalonia (23.6 and 18.9 mm yr⁻¹ in east and north), our model gives a velocity of 16.9 and 5.7 mm yr⁻¹ at that point, thus a vector relative to Apulia of 6.7 and 13.2 mm yr⁻¹. Projected along the fault axis azimuth, this corresponds to 14.3 mm yr⁻¹ for the shear velocity along the fault at that latitude of 38.25°N, thus significantly less than the 19.5 mm yr⁻¹ proposed by Briole et al. (2015). There is also a significant component of shortening of 3.8 mm yr⁻¹. At the northern tip of the Cephalonia fault zone, in the north of Lefkada, we measure 9.4 mm yr⁻¹ of shear between the block SI and the southeast of the block WE, thus 4.9 mm yr⁻¹ less than the fault velocity at the south of Cephalonia. This decrease of velocity on the Cephalonia fault zone at its northern termination is consistent with the description made by Valkaniotis et al. (2020) of the location and size of the deforming band between the Apulian platform and the Pindos block. According to them, the off shore extension of the deforming band is ~30 km. Therefore, the progressive attenuation of the slip rate on the Cephalonia fault is expected to occur on this length scale. The amplitude of the attenuation should be ~4.4 mm yr⁻¹ (half convergence rate at the latitude of GARD) plus ~0.4 mm yr⁻¹ due to the north-south increase of convergence rate induced by the rotation of the Pindos block (see Section 5.5), thus a total of ~4.8 mm yr⁻¹ of difference which fits very well with the 4.9 mm yr⁻¹ of the model.

The right lateral shear in the transition zone Corinth-Evia is around 1 mm yr⁻¹.

## 6 CONCLUSIONS

The velocity field presented here for Greece is taking into account precisely the coseismic and post-seismic displacements produced by the earthquakes. This is something that Müller (2011) did already in his thesis, quantifying and correcting the effects of the Strofades 1997 (Hollenstein et al. 2008), Skyros 2001 and Zakynthos 2006 earthquakes, for which there were GPS data available from permanent and campaign stations. Here we extended the analysis to the period 2000–2020 and we added the correction of the post-seismic deformations. Our coordinates of the stations at the epoch 2020.0
(Supporting Information Table S7) and velocity field (Supporting Information Table S1) provide a precise framework for the nesting of denser local geodetic networks, such as the one of the western Corinth rift (http://crlab.eu), and for the nesting of any geodetic work in Greece. After applying this correction for the effect of the earthquakes we constructed a crustal blocks model, comprising ten blocks. This model reproduces the Greek GPS velocity field within 1 mm yr⁻¹ on average on the entire country.

Compared to the duration of the seismic cycles, the time period sampled by GPS is short, from a few years to 20 yr, with an average of 10 yr. This time interval allowed us to draw a number of conclusions about the role of earthquakes based on the population of earthquakes occurring in the period. But what is the representativeness of this population of earthquakes and what would be the situation with a different period? Would a longer period lead to very different conclusions or not so different? With the data already available can we already draw scenarios over periods longer than the 20 yr of GPS? If we consider a period twice as long, we should have integrated important earthquakes such as those of Thessaloniki 1978 (Stiros & Drakos 2000), Volos 1980 (Drakos et al. 2001), Corinth 1981 (Stiros et al. 2007), Kalamata 1986 (Stiros & Kontogianni 2008), Grevena 1995 (Meyer et al. 1996) and Aigion 1995 (Bernard et al. 1997), which affected regions different from those affected by the earthquakes of the period 2000–2020. We note that the 1981 and June 1995 events occurred along a block boundary (CA-CG), while the 1978, 1980, 1986 and May 1995 events, all with extensional kinematics, occurred inside several of our blocks. This highlights the fact that significant amounts of seismic strains are accumulated and released across fault zones located inside the deforming blocks.

We observed a transient signal affecting stations located between Zakynthos and the Peloponnese in May 2018 (Supporting Information Fig. S10). This is the first detection of this type of transient made in Greece. It is also reported by Mouslopoulou et al. (2020) yet their estimate of the magnitude and the fault size is too large and does not fit the GPS data. We modelled this event as slow aseismic slip on the subduction interface. The fact that no other events of this type have been found in all the time-series studied does not mean that there were no others in the period, but simply shows that, if this is the case, their amplitude was such that they passed through the meshes of the net allowed by the GPS.

The velocity field presented in this paper, relatively dense and accurate, gives a more accurate vision of the velocity field of Greece. For example, the fast relative motion across the Katouna fault (Pérouse et al. 2017), is better quantified, as well as the extension across the Amorgos Basin in the area of the 1956 earthquake (Okal et al. 2009). Like Mouslopoulou et al. (2014) we find possible evidence of tectonic activity along the Strymon river valley, which is a site of crustal weakness because of the previous phase of Miocene extension in central-north Macedonia (Dinter 1998), but there is not enough resolution to further discuss this observation in terms of activity of a specific fault and in view of geological data suggesting north-south extension during Quaternary (Tranos 2011).

We find a significant gradient of extensional deformation in the marine basin between Lesvos and Chios, while the gradient is minor between Lemnos and Lesvos. In the latter region, we note the synchronous activation of northeast/southwest right lateral faults (Pavlides & Tranos 1991; Ganas et al. 2014) and east-west normal faults, such as the 2017 Gulpinar sequence (Ganas et al. 2018).

The GPS networks, including most of the 102 academic stations, were not designed with enough care to capture well the vertical

tectonic velocities. A majority of the GPS stations are located on buildings and sedimentary areas possibly affected by deformations in the subsurface of anthropogenic or natural origin. This is a weakness to be addressed when upgrading the scientific GPS networks in the future. Improving the understanding of the geodynamic processes, and the forces involved inside and below the brittle crust (on top of which the GPS observations are made), requires more high-quality data on the vertical.

## ACKNOWLEDGEMENTS

We are grateful to the EUREF, IGS and several academic institutions from Greece (NOA, NTUA, TEIC, AUTH, TUC, NKUA, U. Patras), France (CNRS), Bulgaria (BAS), Italy (INGV), Czech Republic (Charles University), UK (COMET) and USA (UNAVCO) for distributing their GPS data through open channels. We are grateful to the companies TreeComp and Metrica SA (HxGN SmartNet) who provided us with the GPS data of their stations. The data from Albania are distributed by ALBPOS. The time-series of the Turkish stations (GCM network) are those distributed by the Nevada Geodetic Laboratory (Blewitt et al. 2018). We had access to data from three stations of the Greek KTIMATOLOGO SA network. We thank JPL/NASA for the GIPSY 6.4 software that was used to process the GPS data. We thank Nicolas Chamot-Rooke, Simon Bufféral, Hélène Lyon-Caen, Ioannis Kassaras, Sotiris Valkaniotis and Antonio Avallone for fruitful discussions. We thank Vassilis Zorapas (IGME) for providing us borehole water level data. For several figures we used the Generic Mapping Tool (Wessel & Smith 1998) and SRTM and EMODNET topography data. Author contribution statement: The four authors have been deeply involved for 20 yr in the deployment, maintenance and data management of the academic permanent GPS stations from Greece and Bulgaria contributing to this paper. PB and AG gathered the RINEX data. PB processed the GPS data, analysed the time-series, made the blocks model and wrote the paper with the help from co-authors. All authors contributed to the edition and review of the manuscript. We are very grateful to the Editor, one anonymous reviewer and Wayne Thatcher for their very detailed and constructive reviews.

## DATA AVAILABILITY

The GPS data acquired by the academic partners is available through their web portals.

The time-series of the GPS station coordinates are provided in the Supporting Information.

The computer code HelVel2020 that calculates velocities in Greece is provided in the Supporting Information.

The ITRF2014 coordinates of all GPS stations at the epoch 2020.0 are provided in the Supporting Information (Table S1).

## REFERENCES

* Aardoom, L. & Van Gelder, B.G.W., 1984. Satellite laser ranging to measure crustal motion in the eastern Mediterranean area—instrumentation and network design, *Ann. Geophys.*, **2**(3), 249–258.

* Altamimi, Z., Rebischung, P., Metivier, L. & Collilieux, X., 2016. ITRF2014: a new release of the International Terrestrial Reference Frame modeling nonlinear station motions, *J. geophys. Res.*, **121**, 6109–6131.

* Ambraseys, N.N. & Jackson, J.A., 1990. Seismicity and associated strain of central Greece between 1890 and 1988, *Geophys. J. Int.*, **101**(3), 663–708

* Ambraseys, N.N. & Jackson, J.A., 1997. Seismicity and strain in the gulf of Corinth (Greece) since 1694, *J. Earthq. Eng.*, **1**(3), 433–474.
Armijo, R., Lyon-Caen, H. & Papanastassiou, D., 1992. East-west extension and Holocene normal-fault scarps in the Hellenic arc, *Geology*, **20**(6), 491–494.

Avallone, A. et al. 2004. Analysis of eleven years of deformation measured by GPS in the Corinth Rift Laboratory area, *C. R. Geosci.*, **336**(4–5), 301–311.

Baker, C., Hatzfeld, D., Lyon-Caen, H., Papadimitriou, E. & Rigo, A., 1997. Earthquake mechanisms of the Adriatic Sea and Western Greece: implications for the oceanic subduction-continental collision transition. *Geophys. J. Int.*, **131**, 559–594.

Becker, M. et al. 2002. Assessment of height variations by GPS at Mediterranean and Black Sea coast tide gauges from the SELF projects, *Glob. Planet. Change*, **34**(1–2), 5–35.

Berardino, P., Fornaro, G., Lanari, R. & Sansosti, E., 2002. A new algorithm for surface deformation monitoring based on small baseline differential SAR interferograms, *IEEE Transactions on Geosciences and Remote Sensing*, **40**(11),2375–2383, https://doi.org/10.1109/TGRS.2002.803792

Bernard, P. et al. 1997. The $M_s = 6.2$, June 15, 1995 Aigion earthquake (Greece): evidence for low angle normal faulting in the Corinth rift, *J. Seismol.* 1(1997), 131–150.

Billiris, H., England, P., Floyd, M., Nocquet, J.M., Palamartchouk, K., Paradissis, D., Parsons, B. & Raptakis, C., 2014. COMET: continuous GPS measurements from the Aegean from 2002 to present, NERC Earth Observation Data Centre, http://catalogue.ceda.ac.uk/uuid/8aa71e3a6899480baea9d0f7e8e21239

Billiris, H. et al. 1991. Geodetic determination of tectonic deformation in central Greece from 1900 to 1988, *Nature*, **350**(6314), 124–129.

Bitharis, S., Papadopoulos, N., Pikridas, C., Fotiou, A., Rossikopoulos, D. & Kagiadakis, V., 2019. Assessing a new velocity field in Greece towards a new semi-kinematic datum, *Surv. Rev.*, **51**(368), 450–459.

Blewitt, G., Hammond, W.C. & Kreemer, C., 2018. Harnessing the GPS data explosion for interdisciplinary science, *EOS, Trans. Am. geophys. Un.*, **99**, https://doi.org/10.1029/2018EO104623

Briole, P., 2020. Monitor (estimates seismic fault parameters from available centroid moment tensors, version 3.2). *Zenodo*, http://doi.org/10.5281/zenodo.4295872

Briole, P. et al. 2000. Active deformation of the Corinth rift, Greece: results from repeated Global Positioning System surveys between 1990 and 1995, *J. geophys. Res.*, **105**(B11), 2528–2540.

Briole, P. et al. 2015. The seismic sequence of January–February 2014 at Cephalonia Island (Greece): constraints from SAR interferometry and GPS, *Geophys. J. Int.*, **203**(3), 1528–1540.

Caputo, R., Catalano, S., Monaco, C., Romagnoli, G., Tortorici, G. & Tortorici, L., 2010. Active faulting on the island of Crete (Greece), *Geophys. J. Int.*, **183**(1), 111–126.

Casu, F., 2014. SBAS-DInSAR parallel processing for deformation time-series computation, *IEEE JSTARS*, **7**(8),3285–3296, https://doi.org/10.1109/JSTARS.2014.2322671.

Chousianitis, K., Ganas, A. & Evangelidis, C.P., 2015. Strain and rotation rate patterns of mainland Greece from continuous GPS data and comparison between seismic and geodetic moment release, *J. geophys. Res.*, **120**(5), 3909–3931.

Chousianitis, K., Ganas, A. & Gianniou, M., 2013. Kinematic interpretation of present-day crustal deformation in central Greece from continuous GPS measurements, *J. Geodyn.*, **7**, 1–13.

Christodoulidis, D.C., Smith, D.E., Kolenkiewicz, R., Klosko, S.M., Torrence, M.H. & Dunn, P.J., 1985. Observing tectonic plate motions and deformations from satellite laser ranging, *J. geophys. Res.*, **90**(NB11), 9249–9263.

Clarke, P.J. et al. 1998. Crustal strain in central Greece from repeated GPS measurements in the interval 1989–1997, *Geophys. J. Int.*, **135**(1), 195–214.

Cocard, M., Kahle, H.G., Peter, Y., Geiger, A., Veis, G., Felekis, S., Paradissis, D. & Billiris, H., 1999. New constraints on the rapid crustal motion of the Aegean region: recent results inferred from GPS measurements (1993-1998) across the West Hellenic Arc, Greece, *Earth planet. Sci. Lett.*, **172**(1-2), 39–47.

D'Agostino, N. et al. 2020. Active crustal deformation and rotations in the southwestern Balkans from continuous GPS measurements, *Earth planet. Sci. Lett.*, **539**, 116246, https://doi.org/10.1016/j.epsl.2020.116246

Devoti, R. et al. 2017. A combined velocity field of the Mediterranean region, *Ann. Geophys.*, **60**(2), https://doi.org/10.4401/ag-7059

Dinter, D.A., 1998. Late Cenozoic extension of the Alpine collisional orogen, northeastern Greece: origin of the north Aegean basin, *Bull. geol. Soc. Am.*, **110**(9), 1208–1230.

Drakos, A.G., Stiros, S.C. & Kiratzi, A.A., 2001. Fault parameters of the 1980 ($M_w$ 6.5) Volos, central Greece, earthquake from inversion of repeated leveling data, *Bull. seism. Soc. Am.*, **91**(6), 1673–1684.

Ekström, G. & England, P., 1989. Seismic strain rates in regions of distributed continental deformation, *J. geophys. Res.*, **94**(B8), 10 231–10 257

Ekström, G., Nettles, M. & Dziewonski, A.M., 2012. The global CMT project 2004–2010: centroid-moment tensors for 13,017 earthquakes, *Phys. Earth planet. Inter.*, **200–201**, 1–9.

Elias, P. & Briole, P., 2018. Ground deformations in the Corinth rift, Greece, investigated through the means of SAR multitemporal interferometry, *Geochem. Geophys. Geosyst.*, **19**(12), 4836–4857.

Elias, P. 2013. Ground deformation observed in the western Corinth rift (Greece) by means of SAR interferometry, *PhD thesis*, Ecole Normale Supérieure of Paris and University of Patras.

England, P., Houseman, G. & Nocquet, J.M., 2016. Constraints from GPS measurements on the dynamics of deformation in Anatolia and the Aegean, *J. geophys. Res.*, **121**(12), 8888–8916.

Flotté, N., Sorel, D., Müller, C. & Tensi, J., 2005. Along strike changes in the structural evolution over a brittle detachment fault: example of the Pleistocene Corinth–Patras rift (Greece), *Tectonophysics*, **403**, 77–94.

Floyd, M.A. et al. 2010. A new velocity field for Greece: implications for the kinematics and dynamics of the Aegean, *J. geophys. Res.*, **115**, B10403, http://dx.doi.org/10.1029/2009JB007040

Ganas, A., 2020. NOAFAULTS KMZ layer Version 3.0 (2020 update) [Data set], *Zenodo*, http://doi.org/10.5281/zenodo.4304613

Ganas, A., Karastathis, V., Moshou, A., Valkaniotis, S., Mouzakiotis, E. & Papathanassiou, G., 2014. Aftershock relocation and frequency-size distribution, stress inversion and seismotectonic setting of the 7 August 2013 M = 5.4 earthquake in Kallidromon Mountain, central Greece, *Tectonophysics*, **617**, 101–113.

Ganas, A., Kourkouli, P., Briole, P., Moshou, A., Elias, P. & Parcharidis, I., 2018. Coseismic displacements from moderate-size earthquakes mapped by Sentinel-1 differential interferometry: the case of February 2017 Gulpinar earthquake sequence (Biga Peninsula, Turkey), *Remote Sens.*, **10**(7), 1089, https://doi.org/10.3390/rs10071089

Ganas, A., Marinou, A., Anastasiou, D., Paradissis, D., Papazissi, K., Tzavaras, P. & Drakatos, G., 2013. GPS-derived estimates of crustal deformation in the central and north Ionian Sea, Greece: 3-yr results from NOANET continuous network data, *J. Geodyn.*, **6**, 62–71.

Ganas, A., Mouzakiotis, E., Moshou, A. & Karastathis, V., 2016a. Left-lateral shear inside the North Gulf of Evia Rift, Central Greece, evidenced by relocated earthquake sequences and moment tensor inversion, *Tectonophysics*, **682**, 237–248.

Ganas, A. & Papoulia, I., 2000. High-resolution digital mapping of the seismic hazard within the Gulf of Evia rift, Eastern Central Greece using normal fault segments as line sources, *Nat. Hazards*, **22**, 203–223.

Ganas, A. et al. 2008. NOANET: the new permanent GPS network for geodynamics in Greece, in *Geophysical Research Abstracts*, Vol. 10, EGU2008-A-04380. EGU General Assembly 2008, https://www.cosis.net/abstracts/EGU2008/04380/EGU2008-A-04380.pdf.

Ganas, A. et al. 2016b. Coseismic deformation, field observations and seismic fault of the 17 November 2015 M = 6.5, Lefkada Island, Greece earthquake, *Tectonophysics*, **687**, 210–222.

Ganas, A. et al. 2019. The July 20, 2017 M6.6 Kos Earthquake: seismic and geodetic evidence for an active north-dipping normal fault at the western end of the Gulf of Gökova (SE Aegean Sea), *Pure appl. Geophys.*, **176**(10), 4177–4211.

Ganas, A. et al. 2020. The 25 October 2018 $M_w$ = 6.7 Zakynthos earthquake (Ionian Sea, Greece): a low-angle fault model based on GNSS data, relocated seismicity, small tsunami and implications
for the seismic hazard in the west Hellenic Arc, *J. Geodyn.*, **137**,
https://doi.org/10.1016/j.jog.2020.101731

Goldsworthy, M. & Jackson, J., 2001. Migration of activity within normal fault systems: examples from the Quaternary of mainland Greece, *J. Struct. Geol.*, **23**(2–3), 489–506.

Haddad, A., Ganas, A., Kassaras, I. & Lupi, M., 2020. Seismicity and geodynamics of western Peloponnese and central Ionian Islands: insights from a local seismic deployment, *Tectonophysics*, **778**, 228353, https://doi.org/10.1016/j.tecto.2020.228353

Hatzfeld, D., Ziazia, M., Kementzetzidou, D., Hatzidimitriou, P., Pana-
giotopoulos, D. & Makropoulos, K., 1999. Microseismicity and focal
mechanisms at the western termination of the North Anatolian Fault and
their implications for continental tectonics, *Geophys. J. Int.*, **137**(3), 891-
908.

Hollenstein, C., Geiger, A., Kahle, H.G. & Veis, G., 2006. CGPS time-series and trajectories of crustal motion along the West Hellenic Arc, *Geophys. J. Int.*, **164**(1), 182–191

Hollenstein, C., Muller, M., Geiger, A. & Kahle, H.G., 2008. GPS-derived coseismic displacements associated with the 2001 Skyros and 2003 Lefkada earthquakes in Greece, *Bull. seism. Soc. Am.*, **98**(1), 149–161.

Howell, A., Jackson, J., Copley, A., McKenzie, D. & Nissen, E., 2017b. Subduction and vertical coastal motions in the eastern Mediterranean, *Geophys. J. Int.*, **211**(1), 593–620.

Howell, A., Palamartchouk, K., Papanikolaou, X., Paradissis, D., Raptakis, C., Copley, A., England, P. & Jackson, J., 2017a. The 2008 Methoni earthquake sequence: the relationship between the earthquake cycle on the subduction interface and coastal uplift in SW Greece, *Geophys. J. Int.*, **208**(3), 1592–1610.

Ilieva, M., Briole, P., Ganas, A., Dimitrov, D., Elias, P., Mouratidis, A. & Charara, R., 2016. Fault plane modelling of the 2003 August 14 Lefkada Island (Greece) earthquake based on the analysis of ENVISAT SAR interferograms, *Tectonophysics*, **693**, 47–65.

Jackson, J., Haines, J. & Holt, W., 1994. A comparison of satellite laser ranging and seismicity data in the Aegean region, *Geophys. Res. Lett.*, **21**(25), 2849–2852.

Kahle, H.G., Cocard, M., Peter, Y., Geiger, A., Reilinger, R., Barka, A. & Veis, G., 2000. GPS-derived strain rate field within the boundary zones of the Eurasian, African and Arabian Plates, *J. geophys. Res.*, **105**(10), 23 353–23 370

Kahle, H.G., Muller, M.V., Mueller, S. & Veis, G., 1993. The Kephalonia transform-fault and the rotation of the Apulian platform – evidence from satellite geodesy, *Geophys. Res. Lett.*, **20**(8), 651–654.

Kapetanidis, V. & Kassaras, I., 2019. Contemporary crustal stress of the Greek region deduced from earthquake focal mechanisms, *J. Geodyn.*, **123**, 55–82.

Kassaras, I., Kapetanidis, V. & Karakonstantis, A., 2016. On the spatial dis-
tribution of seismicity and the 3D tectonic stress field in western Greece,
*Phys. Chem. Earth* **95**, 50–72.

Kiratzi, A. et al. 2008. The April 2007 earthquake swarm near Lake Trichonis and implications for active tectonics in western Greece, *Tectonophysics*, **452**, 51–65.

Kreemer, C. & Chamot-Rooke, N., 2004. Contemporary kinematics of the southern Aegean and the Mediterranean ridge, *Geophys. J. Int.*, **157**(3), 1377–1392.

Le Pichon, X., Chamot-Rooke, N., Lallemant, S., Noomen, R. & Veis, G., 1995. Geodetic determination of the kinematics of central Greece with respect to Europe—implications for eastern Mediterranean tectonics, *J. geophys. Res.*, **100**(B7), 12 675–12 690.

Le Pichon, X. & Kreemer, C., 2010. The Miocene-to-present kinematic evo-
lution of the eastern Mediterranean and Middle East and Its implications
for dynamics, *Annu. Rev. Earth Planet. Sci.*, **38**, 323–351.

Louvari, E., Kiratzi, A., Papazachos, B. & Hatzidimitriou, P., 2001. Fault-plane solutions determined by waveform modeling confirm tec- tonic collision in the eastern Adriatic, *Pure appl. Geophys.*, **158**, 1613–1637.

Louvari, E., Kiratzi, A.A. & Papazachos, B.C., 1999. The Cephalonia Trans-
form Fault and its extension to western Lefkada Island, *Tectonophysics*,
**308**(1), 223–236.

Lyon-Caen, H. et al. 1988. The 1986 Kalamata (south Peloponnesus) earthquake: detailed study of a normal fault, evidences for east–west extension in the Hellenic arc, *J. geophys. Res.*, **93**(B12), 14 967–15 000.

Marinou, A., Anastasiou, D., Papanikolaou, X., Paradissis, D. & Zacharis, V. 2016. Dionysos satellite observatory and higher geodesy laboratory: history and perspectives, *Bull. Geol. Soc. Greece*, **50**(2), 1091–1099.

Melis, N.S., Okal, E.A., Synolakis, C.E., Kalogeras, I.S. & Kanoglu, U., 2020. The Chios, Greece earthquake of 23 July 1949: seismological reassessment and tsunami investigations, *Pure appl. Geophys.*, **177**(3), 1295–1313.

Meyer, B. et al. 1996. The 1995 Grevena (Northern Greece) earthquake: fault model constrained with tectonic observations and SAR interferometry, *Geophys. Res. Lett.*, **23**(19), 2677–2680.

Mouslopoulou, V., Bocchini, G.M., Cesca, S., Saltogianni, V., Bedford, J., Petersen, G., Giannou & Oncken, O., 2020. Earthquake-squaws, slop-slip and fault-interactions at the western-end of the Hellenic Subduction System precede the M<sub>w</sub> 6.9 Zakynthos earthquake, Greece, *Geochem. Geophys. Geosyst.*, **21**(12), e2020GC009243, https://doi.org/10.1029/2020GC009243

Mouslopoulou, V., Saltogianni, V., Giannou, M. & Stiros, S., 2014. Geodetic evidence for tectonic activity on the Strymon fault system, northeast Greece, *Tectonophysics*, **633**, 246–255.

Métois, M. et al. 2015. Insights on continental collisional processes from GPS data: dynamics of the peri-Adriatic belts, *J. geophys. Res.*, **120**(12), 8701–8719.

Müller, M.D, 2011. Analysis of long-term GPS observations in Greece (1993-2009) and geodynamic implications for the Eastern Mediterranean, *PhD thesis*, ETH Zürich.

Müller, M.D., Geiger, A., Kahle, H.G., Veis, G., Billiris, H., Paradissis, D. & Felekis, S., 2013. Velocity and deformation fields in the North Aegean domain, Greece, and implications for fault kinematics, derived from GPS data 1993–2009, *Tectonophysics*, **597**, 34–49.

Nocquet, J.M., 2012. Present-day kinematics of the Mediterranean: a com-
prehensive overview of GPS results, *Tectonophysics*, **579**, 220–242.

Nomikou, P., Hübscher, C., Papanikolaou, D., Farangitakis, G.P., Ruhnau, M. & Lampridou, D., 2018. Expanding extension, subsidence and lateral segmentation within the Santorini–Amorgos basins during Quaternary: implications for the 1956 Amorgos events, central-south Aegean Sea, Greece, *Tectonophysics*, **722**, 138–153.

Noomen, R., Springer, T.A., Ambrosius, B.A.C., Herzberger, K., Kuijper, D.C., Mets, G.J., Overgaauw, B. & Wakker, K.F., 1996. Crustal deformations in the Mediterranean area computed from SLR and GPS observations, *J. Geodyn.*, **21**(1), 73–96.

Nyst, M. & Thatcher, W., 2004. New constraints on the active tec-
tonic deformation of the Aegean, *J. geophys. Res.*, **209**(B11), B11406,
https://doi.org/10.1029/2003JB002830.

Okada, Y., 1992. Internal deformation due to shear and tensile faults in a
half space, *Bull. seism. Soc. Am.*, **82**, 1018–1040.

Okal, E.A., Synolakis, C.E., Uslu, B., Kalligeris, N. & Voukouvalas, E., 2009. The 1956 earthquake and tsunami in Amorgos, Greece, *Geophys. J. Int.*, **178**(3), 1533–1554.

Papadimitriou, P. et al. 2018. The 12th June 2017 M<sub>w</sub> = 6.3 Lesvos earthquake from detailed seismological observations, *J. Geodyn.*, **115**, 23–42.

Papadopoulos, A., Parcharidis, I., Elias, P. & Briole, P., 2019. Spatio-
temporal evolution of the deformation around the Rio-Patras fault
(Greece) observed by synthetic aperture radar interferometry from 1993
to 2017, *Int. J. Remote Sens.*, **40**(16), 6365–6382.

Pavlides, S.B. & Tranos, M.D., 1991. Structural characteristics of two strong earthquakes in the North Aegean: ierissos (1932) and Agios-Efstratios (1968), *J. Struct. Geol.*, **13**(2), 205–214.

Pérouse, E., Chamot-Rooke, N., Rabaute, A., Briole, P., Jouanne, F., Georgiev, I. & Dimitrov, D., 2012. Bridging onshore and offshore present-day kinematics of central and eastern Mediterranean: implications for crustal dynamics and mantle flow, *Geochem. Geophys. Geosyst.*, **13**, Q09013, https://doi.org/10.1029/2012GC004289.

Pérouse, E. et al. 2017. Transition from collision to subduction in Western Greece: the Katouna-Stamna active fault system and regional kinematics, *Int. J. Earth Sci.*, **106**(3), 967–989.
Reilinger, R. et al. 2006. GPS constraints on continental deformation in the Africa–Arabia–Eurasia continental collision zone and implications for the dynamics of plate interactions. *J. geophys. Res.*, BO5411. https://doi.org/10.1029/2005JB004051

Reilinger, R. et al. 2010. Geodetic constraints on the tectonic evolution of the Aegean region and strain accumulation along the Hellenic subduction zone, *Tectonophysics*, **488**(1–4), 22–30.

Roberts, G.P. & Ganas, A., 2000. Fault-slip directions in central and southern Greece measured from striated and corrugated fault planes: comparison with focal mechanism and geodetic data, *J. geophys. Res.*, **105**(B10), 23443–23462.

Robertson, J., Meschis, M., Roberts, G.P., Ganas, A. & Gheorghiu, D.M., 2019. Temporally constant Quaternary uplift rates and their relationship with extensional upper-plate faults in south Crete (Greece), constrained with CI-36 cosmogenic exposure dating, *Tectonics*, **38**(4), 1189–1222.

Saltogianni, V., Mouslopoulou, V., Oncken, O., Nicol, A., Gianniou, M. & Mertikas, S., 2020. Elastic Fault Interactions and Earthquake Rupture Along the Southern Hellenic Subduction Plate Interface Zone in Greece, *Geophys. Res. Lett.*, **47**(13), https://doi.org/10.1029/2019GL08660

Savage, J.C., 2018. Euler-vector clustering of GPS velocities defines microplate geometry in southwest Japan, *J. geophys. Res.*, **123**, 1954–1968.

Savage, J.C. & Simpson, R.W., 2013a. Clustering of GPS velocities in the Mojave Block, southeastern California, *J. geophys. Res.*, **118**, 1747–1759.

Savage, J.C. & Simpson, R.W., 2013b. Clustering of velocities in a GPS network spanning the Sierra Nevada Block, the Northern Walker Lane Belt, and the Central Nevada Seismic Belt, California-Nevada, *J. geophys. Res.*, **118**, 4937–4947.

Savage, J.C. & Wells, R.E., 2015. Identifying block structure in the Pacific Northwest, USA, *J. geophys. Res.*, **120**, 7905–7916.

Serpelloni, E., Anzidei, M., Baldi, P., Casula, G. & Galvani, A., 2005. Crustal velocity and strain-rate fields in Italy and surrounding regions: new results from the analysis of permanent and non-permanent GPS networks, *Geophys. J. Int.*, **161**(3), 861–880.

Serpelloni, E., Faccenna, C., Spada, G., Dong, D. & Williams, S.D.P., 2013. Vertical GPS ground motion rates in the Euro-Mediterranean region: new evidence of velocity gradients at different spatial scales along the Nubia-Eurasia plate boundary, *J. geophys. Res.*, **118**(11), 6003–6024.

Serpetsidaki, A. et al. 2014. New constraints from seismology and geodesy on the M<sub>w</sub> = 6.4 2008 Movri (Greece) earthquake: evidence for a growing strike-slip fault system, *Geophys. J. Int.*, **198**(3), 1373–1386.

Simpson, R. W., Thatcher, W. & Savage, J.C., 2012. Using cluster analysis to organize and explore regional GPS velocities, *Geophys. Res. Lett.*, **39**(18), L18307, https://doi.org/10.1029/2012GL052755.

Smith, D.E., Kelenkiewicz, R., Robbins, J.W., Dunn, P.J. & Torrence, M.H., 1994. Horizontal crustal motion in the central and eastern Mediterranean inferred from satellite laser ranging measurements, *Geophys. Res. Lett.*, **21**(18), 1979–1982.

Stiros, S. & Kontogianni, V., 2008. Modelling of the 1986 Kalamata (SW Greece) earthquake faulting using geodetic data, *J. Appl. Geod.*, **2**, 179–185.

Stiros, S.C. & Drakos, A., 2000. Geodetic constraints on the fault pattern of the 1978 Thessaloniki (Northern Greece) earthquake (M<sub>s</sub> = 6.4), *Geophys. J. Int.*, **143**(3), 679–688.

Stiros, S.C., Psimoulis, P. & Pitharouli, S., 2007. Geodetic constraints to the kinematics of the Kapareli fault, reactivated during the 1981, Gulf of Corinth earthquakes, *Tectonophysics*, **440**(1–4), 105–119.

Thingbaijam, K.K.S., Martin Mai, P. & Goda, K., 2017. New empirical earthquake source-scaling laws, *Bull. seism. Soc. Am.*, **107**(5), 2225–2246.

Tranos, M.D., 2011. Strymon and Strymonikos Gulf basins (Northern Greece): implications on their formation and evolution from faulting, *J. Geodyn.*, **51**(4), 285–305.

Valkaniotis, S., Briole, P., Ganas, A., Elias, P., Tsironi, V., Fokaefs, A., Partheniou, H. & Paschos, P., 2020. The M<sub>w</sub> = 5.7 Kanallaki earthquake of March 21, 2020 in Epirus, Greece: reverse fault model and seismotectonic implications for Apulia-Eurasia collision, *Geosciences*, **10**, 454, https://doi.org/10.3390/geosciences10110454.

Vernant, P., Reilinger, R. & McClusky, S., 2014. Geodetic evidence for low coupling on the Hellenic subduction plate interface, *Earth planet. Sci. Lett.*, **385**(C), 122–129.

Wells, D.L. & Coppersmith, K.J., 1994. New empirical relationships among magnitude, rupture length, rupture width, rupture area, and surface displacement, *Bull. seism. Soc. Am.*, **84**(4), 974–1002.

Wessel, P. & Smith, W., 1998. New improved version of the Generic Mapping Tool released, *EOS, Trans. Am. geophys. Un.*, **79**, 579.

Willis, P., Mertikas, S., Argus, D.F. & Bock, O., 2013. DORIS and GPS monitoring of the Gavdos calibration site in Crete, *Adv. Space Res.*, **51**, 1438–1447.

## SUPPORTING INFORMATION

Supplementary data are available at GJI online.

**Figure S1.** Greece and surroundings with the location of the 21 earthquakes listed in Table 1, and the seismicity from the NOA catalogue for the period 2000–2020, magnitude ≥ 4 and depth ≤ 35 km. The beachballs show the focal mechanism of the studied events, they are centred at the location of the epicentres.

**Figure S2.** Time-series of the north coordinate of the GPS station LEMN. The large offset corresponds to the nearby Samothraki earthquake of 2014 May 24 (Tables 1 and 2). After the earthquake the velocity of the station has changed. The transient can be modelled with a sum of two exponential terms of time constant 65 and 1295 d, respectively.

**Figure S3.** Additional toponymy for features not mentioned in the main text but in the Supporting Information (coordinates and description of the sites in Table S2).

**Figure S4.** Duration of the time-series of GPS data for the 329 stations. The average duration is 6.5 yr.

**Figure S5.** Completeness of the time-series acquired at the 329 GPS stations. The average for the whole set of stations is 85 per cent.

**Figure S6.** Comparison of our solution with the solution of the Nevada Geodetic Laboratory (NGL), here for the time-series of the station LEMN. The curves are de-trended from the modelled effect of the earthquakes and the average secular velocity. The difference between the slopes is below 0.05 mm yr⁻¹. See the native time-series of the north component of LEMN in Fig. S2.

**Figure S7.** Histogram of east and up velocity uncertainties at the 329 GPS stations.

**Figure S8.** Time-series of the RLSO station showing the coseismic offset of the M<sub>w</sub> = 6.4 2008 June 8 Movri earthquake of amplitude −3.4 mm in east and 9.8 mm in north, and the post-seismic deformation with relaxation time of 103 d and amplitudes −1.1 and 1.4 mm in east and north.

**Figure S9.** Time-series of the stations MET0, AGTH, AIGU and ATHI. The vertical line corresponds to the M<sub>w</sub> = 5.3 Magoula earthquake of 2019 July 19. At AIGU and ATHI we measure an offset at the date of the earthquake with amplitudes that are reported in Table S4.

**Figure S10.** Residual time-series of the east component of various stations corrected for coseismic, post-seismic and secular terms (top) and common mode noise correction (bottom) estimated using the eight top stations (MET0 to PAT0). The solid vertical line is on 2018 October 25, day of the M<sub>w</sub> = 6.8 Zakynthos earthquake. The dashed line corresponds to a transient visible at the stations AMAL and PYRG on 2018 May 10 (see Section S4.1 and see Fig. S11 for the north and up components).

**Figure S11.** North and up components of the same stations whose east component is plotted in Fig. S10. The solid vertical line is on 2018 October 25, day of the M<sub>w</sub> = 6.8 Zakynthos earthquake. The
dashed line corresponds to the transient observed in Fig. S10 at the stations AMAL and PYRG on 2018 May 10. Here, a subsidence is observed at TROP, PYRG and AMAL around 2018 May 10.

**Figure S12.** Time-series of the east and up components of the GPS station NOMI used to assess the linear character of the deformation during unrest and the starting and ending dates of the volcanic unrest. We see also a large annual term of the vertical component that was modelled and removed.

**Figure S13.** Time-series of the east component of the GPS stations of Santorini. (a) Original time-series, (b) time-series corrected from a linear signal starting on 2011 February 12 and ending on 2012 May 3 (446 d). See Table S5 for the amplitude of the slope at the stations.

**Figure S14.** Time-series of the vertical coordinate of the GPS stations VALI and XILI. The secular trends are $2.4 \text{ mm yr}^{-1}$ for VALI and $-6.2 \text{ mm yr}^{-1}$ for XILI. The Zakynthos earthquake of October 2018 (vertical line) produced a subsidence of the two stations that we model with a decreasing exponential starting the day of the earthquake with time constant 36 d for VALI and 16 d for XILI and amplitude 10.2 and 9.6 mm, respectively.

**Figure S15.** Sentinel-1 interferogram of the western Gulf of Corinth made with the scenes acquired on 2013 October 20 and November 13 on the descending track 80, using DIAPASON InSAR Sentinel-1 TOPSAR version 1.1.8 on the GEP-TEP (http://geohazards-tep.eu). One cycle of colour corresponds to 28 mm of displacement along the line of sight. The geoinformation layer is from OpenTopoMap (http://opentopomap.org).

**Figure S16.** InSAR ground velocities in the region of Larissa (Thessaly) in the line of sight of the Sentinel-1 descending track 7, measured for the period 2019 October 16 to 2010 May 31, using the CNR-IREA P-SBAS processor on the GEP-TEP (Berardino et al. 2002; Casu et al. 2014). Red dots are permanent GPS stations. Blue colour on the map indicates range increase from $-5 \text{ cm yr}^{-1}$ (uplift) to $5 \text{ cm yr}^{-1}$ (subsidence).

**Figure S17.** Time-series of the station TRIP. On 2019 January 12 (vertical line) a transient deformation starts with a sudden motion of $-1.9 \text{ mm}$ in east, $0.2 \text{ mm}$ in north and $0.0 \text{ mm}$ in up, followed by a progressive deformation that we model with an exponential of time constant 17 d and amplitude $-3.3 \text{ mm}$ in east, $-9.7 \text{ mm}$ in north and 29.7 mm in up.

**Figure S18.** Sentinel-1 interferogram of the Tripoli area made with the scenes acquired on 2018 December 31 and 2019 January 24 on the ascending track 175. The two GPS stations TRIP and TRIU (time-series in Fig. S17) are located in the centre of Tripoli.

**Figure S19.** Time-series of the horizontal coordinates of the GPS stations KALU and KALY, Kalymnos Island. The vertical line indicates the date of the Lesvos 2017 earthquake. The amplitudes of the offset are $-1.7$ and $3 \text{ mm}$ in east and north for KALU, and $0.2$ and $3.5 \text{ mm}$ for KALY. The relaxation time of the signal is 16 d. See also Fig. S20.

**Figure S20.** Sentinel-1 interferogram of the Kalymnos Island made with the scenes acquired on 2017 May 25 and June 30 on the descending track 36. The vectors correspond to the offset plotted in Fig. S19.

**Figure S21.** Time-series of the three components of the vector KRIN (located on a landslide) – KOUN (the closest permanent station, considered here as a local reference). See also Figs S22 and S23.

**Figure S22.** Ground velocities in the line of sight of the Sentinel-1 ascending track 175 over the area of KRIN calculated for the period 2017 January 1 to 2020 May 15, using the TRE-ALTAMIRA

FASTVEL processor on the GEP-TEP. The amplitude of the velocity of the station KRIN is $101.5, -21.2$ and $-18.7 \text{ mm yr}^{-1}$ in east, north and up, respectively. The approximate boundary of the landslide (according also to the interferogram made with the CNR-IREA P-SBAS processor, see Fig. S23) is plotted as a light red contour. Red indicates motion away from the satellite.

**Figure S23.** InSAR ground velocities around the station KRIN, represented in the line of sight of the Sentinel-1 ascending track 175 for the period January 1, 2017 to May 15, 2020, using the CNR-IREA P-SBAS processor on the GEP-TEP. Red indicates motion away from the satellite.

**Figure S24.** Time-series of the horizontal coordinates of the EUREF (https://www.epncb.oma.be/) station PAT0. Top is corrected for the coseismic and post-seismic effects of the various earthquakes, bottom is not. For example, the Zakynthos 2018 earthquake is clear in the uncorrected east component. The differences in slope between the two solutions are $-0.67 \text{ mm yr}^{-1}$ in east and $-0.85 \text{ mm yr}^{-1}$ in north. East and north ITRF2014 velocities of PAT0: EUREF: $8.35, -6.59 \text{ mm yr}^{-1}$; our uncorrected velocities: $7.43, -6.65 \text{ mm yr}^{-1}$; our corrected velocities: $8.1, -5.8 \text{ mm yr}^{-1}$.

**Figure S25.** Istanbul profile. Velocities plotted with X axis along the line Istanbul - pole AT/EU and Y profile perpendicular. The solid and dashed line represents the best-fitting linear increase of the velocity (using the five points in red) as a function of the distance to the pole. The corresponding rate is $27.3 \times 10^{-9} \text{ rad yr}^{-1}$ counter-clockwise, that is, $1.56^\circ \text{ Ma}^{-1}$.

**Figure S26.** Velocities along a profile Sandanski (Bulgaria)—Rodos. We plotted of a subset of stations located at less than 55 km from the profile. The axis is aligned with the azimuth of the Euler pole of rotation of Anatolia with respect to Eurasia, located at $31.96^\circ \text{E}$ and $32.02^\circ \text{N}$ according to Le Pichon & Kreemer (2010).

**Figure S27.** Sandanski profile (location of the points and vectors in Fig. S26). Velocities of the stations located in the vicinity ($\pm 55 \text{ km}$) of the line Sandanski-Rodos aligned with the Euler pole of rotation of Anatolia with respect to Eurasia according to Le Pichon & Kreemer (2010). The oblique solid line shows the tangential velocities that are expected for a rigid rotation of Anatolia around the pole plotted in Fig. 7. The grey dots show the values at LIMN, LEMN and MOUD after correcting for the expected component of strain accumulation due to the vicinity of the North Anatolian fault.

**Figure S28.** Velocities along a profile Tirana-Karpathos. We plotted of a subset of stations located at less than 55 km from the axis connecting Tirana and Karpathos aligned with the Euler pole of rotation of Anatolia with respect to Eurasia located at $31.96^\circ \text{E}$ and $32.02^\circ \text{N}$ according to Le Pichon & Kreemer (2010).

**Figure S29.** Tirana-Karpathos profile. Velocities of the stations located in the vicinity ($\pm 55 \text{ km}$) of the line. See location of the points and line and the vectors in Fig. S28.

**Figure S30.** Section of the south Peloponnese Pylos-Spetsopoula

**Figure S31.** East and north velocities in Crete with respect to average.

**Figure S32.** North-south profile across central Macedonia (see points in Fig. 10).

**Figure S33.** Cumulative distribution of the whole set of 322 vertical velocities (grey) and the subset of 177 (black) with velocity uncertainties lower than $1 \text{ mm yr}^{-1}$ (see Table S1). Best Gaussian fits are calculated for 90 per cent of the population, which corresponds to velocities comprises between $-2.7$ and $1.4 \text{ mm yr}^{-1}$ in the case of the 177-station distribution. The remaining 6 per cent on the negative velocities side and 4 per cent on the positive velocities side (drawn with + symbol) do not fit with a Gaussian distribution and were not used to calculate the best fit. The median $\mu$ of
the distribution is −0.38 mm yr⁻¹ with the best 177 stations and −0.27 mm yr⁻¹ for the whole set of 322 stations. The corresponding standard deviation σ is 0.79 and .95 mm yr⁻¹, respectively. Two stations (GEYB and PSAT) present anomalous subsidence discussed in the text, GEYB (mentioned in Section S4.7) is built on a recent thick (>10 m) embankment built during the construction of the bridge Rio-Antirio, PSAT is located near (~1 km) the fault plane of the active Psathopirgos fault, on its footwall. Three stations present anomalous uplift, ZAKY and AMAL presumably affected local uplift (insufficiently modelled by our post-seismic model, see Section S3.14) during the post-seismic relaxation following the 2018 Zakythos earthquake, VALI affected by the uplift of the Selinoutas river delta (see Section S4.3) as shown in Fig. S14, induced by the activity of the offshore fault north of Aigion.

**Figure S34.** The boundaries of the ten blocks of our model. The North Aegean block is the combination of two sub-blocks one for the north (Lemnos and Lesvos), one for the south (Chios, Samos). Here we use one single block (despite the large north-south extension and the localization zone between Lesvos and Chios) because we do not have yet, at the moment, enough GPS stations in the area and long enough time-series to divide the zone in two. The focal mechanisms are for events located between 0 and 35 km of depth, published by the National and Kapodistrian University of Athens (http://www.geophysics.geol.uoa.gr/) for the period 2003–2018.

**Table S1.** Coordinates and velocities in the ITRF2014 (Altamimi et al. 2016) of the 329 points.

**Table S2.** Glossary and coordinates of site names (Fig. 2 and Fig. S3) used in the text, figures and tables. In the column ‘Text’, M means main document and S means Supporting Information.

**Table S3.** Scatter of velocities for 110 pairs of stations separated by less than 15 km. The formal uncertainties provided by the calculations and the linear regression of the time-series are scaled to bet fit with the scatters of velocities. See Fig. 3 for graphic representation of the east, north and vertical scatters.

**Table S4.** Observed and modelled coseismic and post-seismic displacements for the earthquakes listed in Table 1. The faults parameters are in Table 2. We provide those numbers only for the events modelled in this paper. For the event modelled in previous works, the observed and modelled displacements can be found in Ilieva et al. (2016), Serpetsidaki et al. (2014), Elias (2013), Briole et al. (2015) and Ganas et al. (2016b, 2018, 2019, 2020). For the Methoni 2008 earthquake the vertical displacement cannot be assessed for the instantaneous coseismic but it can be assessed for the first 10 d.
**Table S5.** Observed and modelled displacements corresponding for some of the transients presented in Section S4.

**Table S6.** Comparison of solutions (velocities and uncertainties) for the stations AUT1, DYNG, PAT0, TUC2 (blue dots in Fig. 4). The solutions ROB-EUREF, LTK-EUREF and UGA-CNRS are published by the European Plate Observing System (EPOS, https://gnssproducts.epos.ubi.pt). The solution NGL is published by the Nevada Geodetic Laboratory (http://geodesy.unr.edu). Our initial solution is the one before correcting the velocity biases induced by the coseismic and post-seismic displacements.

**Table S7.** ITRF2014 coordinates of the 329 points at the epoch 2020.0 and component of the velocity due to the coseismic and post-seismic displacements (plotted in Fig. 4).

briole-et-al_GPS-points.kml

briole-et-al_analysed-earthquakes.kml

briole-et-al_blocks-model.kml

briole-et-al_tabS1.csv

briole-et-al_tabS3.csv

briole-et-al_tabS7.csv

briole-et-al_GPS-time-series

briole-et-al_helvel2020-model

Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the paper.