## 2 Data, Methodology and Results
Nevertheless, Figure 6 shows that the azimuthal variation in phase velocities is not due to anisotropy but to lateral heterogeneity (Figure 9a; [10]), because Equation (1) does not fit the phase velocities satisfactorily within their errors in general. The best fit of Equation (1) is achieved at the period of 40 s when the coefficients \(\mathrm{Ai}\) (with i varying from 2 to 5) are in general negligible within their error (Figures 7 and 8a), and the error \(\varepsilon\) is greater for the period range in which those are not negligible within their error (e.g., 5-40 s; Figures 7 and 8a). This lateral heterogeneity is compatible with anisotropic models of azimuthal isotropy (i.e., hexagonal symmetry or transversely isotropic [11]), which can be calculated if the study area is divided into regions (Figures 1a and 9b). These regions are established considering the areas with similar Rayleigh-wave phase velocities. The Love-wave phase velocities are not considered because they show bigger scatter than the Rayleigh-wave ones (Figures 7 and 8a). In Figure 10b, the Love (black circles) and Rayleigh (black triangles) phase velocities corresponding to the region 1 are shown. They are calculated with the mean of the values determined for each period for the paths plotted in the upper central part of Figure 10 and the small vertical bars present for the Rayleigh-wave values (Figure 10b), which confirms that the dispersion curves considered for this region have similar values. Then, for each region, an isotropic model is defined to calculate the anisotropic model as a small perturbation of this isotropic model ([5,6]). For region 1 (Figure 9b), Table 1 shows this isotropic model, which is calculated from the isotropic inversion of the Rayleigh-wave dispersion, which is determined for the paths that are within this region (Figure 10b [7]). The S-velocity of this model (Figure 10a) is plotted only for depths above 260 km due to the poor resolution obtained for greater depths in this inversion process (Figure 10c). The Rayleigh-wave dispersion is properly fitted by this model, but the Love-wave dispersion is not fitted (Figure 10b), i.e., the Love- and Rayleigh-wave dispersions are not compatible with a unique isotropic model. This Love-Rayleigh discrepancy makes it necessary to consider the existence of seismic anisotropy. This anisotropy is considered with hexagonal symmetry [11], in the depth range from 10 to 100 km (Figure 11a [5,6]), to fit satisfactorily both dispersion curves (Figure 11b). For the other regions (Figure 9b), the same process is followed, selecting an isotropic model for each region (Supplements T1, T2, T3, T4 and T5) and calculating its perturbation in terms of S-velocities (Supplements F1, F2, F3, F4 and F5  
Figure 2: Love-wave group velocity determination: (**a**) the observed seismogram (red line) corresponding to the event 70 (Supplement 1) recorded at station: RAYN (Supplement 2, transverse component), instrument corrected. The time-variable filtered seismogram (blue line) calculated from the observed seismogram (red line) and the initial group velocity (**b**,**c**). (**b**) Contour map of relative energy (normalized to 99 decibels) as a function of the period and the group time, calculated from the observed seismogram ((**a**), red line) with the MFT. The red points denote the group times inferred from the energy map. (**c**) The initial (red points, (**b**)) and the final (magenta line, (**d**)) group velocities calculated from the group times and the epicentral distance. (**d**) Contour map of relative energy calculated from the time-variable filtered seismogram ((**a**), blue line) with the MFT. The color scale is the same as in (**b**). Figure 4: [PERSON] and Rayleigh group velocities (blue line) determined for the traces of the event 70 (Figures 2c and 3c, magenta line), and their corresponding phase velocities (red line) determined from their respective group velocities by inversion.