Seismograph networks
Joint inversion was conducted on 38 temporary broadband seismograph station deployments across the Turkana Depression and southern Ethiopia (Supplementary Data Table 1). Seismograph deployments include 3 stations from the YY Ethiopian Plateau network51 and 34 stations from the 6R and Y1 Turkana Rift Arrays Investigating Lithospheric Structure (TRAILS) seismograph networks52,53. Additional data were sourced from the permanent GEOFON seismic station LODK located to the west of Lake Turkana. Seismic data were acquired through the National Science Foundation (NSF) Seismological Facility for the Advancement of Geosciences (SAGE) data archive operated by EarthScope Consortium (NSF award 1724509) and GEOFON (GEOFOrschungsNetz) repositories.
Surface-wave dataset
Fundamental-mode Rayleigh-wave group velocity dispersion curves from a 4–60-second period for each station were extracted from the isotropic component of a local anisotropic Rayleigh-wave tomographic model27. Longer-period (70–100 seconds) group velocities were adopted from a global Rayleigh-wave group velocity dispersion model28. The local group velocity dispersion data27 were used for periods of 4–40 seconds or, resolution permitting, 4–60 seconds, with the final joint inversion result independent of the upper cut-off choice. To produce a final dispersion curve for each station, a three-point moving average was calculated to smooth the transition between short- and long-period group velocity data and to mitigate against sharp transitions between neighbouring periods that are otherwise not warranted given the surface-wave period–depth–sensitivity resolution kernels27. Collectively, these periods have peak sensitivities in a range of crustal and uppermost-mantle depths (about 5–150 km)27, with some partial sensitivity down to 400 km.
Receiver function dataset
Teleseismic earthquakes that yielded successful receiver functions were identified from a TRAILS receiver function study that used the same seismograph stations analysed here18. Radial P and PP receiver functions (Extended Data Fig. 1a,c) from 350 earthquakes with magnitude >4 and >5.5, respectively, were calculated using an time domain iterative deconvolution method54. The number of iterations are set at 250, but the deconvolution may terminate earlier if the improvement in fit is less than 0.1% for each station. For P receiver functions, earthquakes with epicentral distances of 30–90° were used, whereas for PP receiver functions, distances of ≥60° were used. Before calculating the receiver functions, seismograms were Butterworth bandpass filtered with corner frequencies of 0.04 Hz and 2 Hz. Subsequently, the seismograms were windowed from 20 seconds before to 100 seconds after the P-wave arrival, then rotated from the vertical–north–east to the vertical–radial–transverse coordinate system.
As per other receiver function18 and joint inversion studies performed in East Africa6,26, the iterative deconvolution variance—where the radial component receiver function is re-convolved with the vertical component seismogram and cross-correlated with the original radial component seismogram—was used to assess the quality of the receiver function: deconvolved traces accounting for less than 90% of the original signal were excluded from the analysis. To increase the signal-to-noise ratio of the accepted deconvolved traces, individual receiver functions for each station were stacked by taking a point-by-point average. However, variations in amplitude and timing of the phases, resulting from different incidence angles (ray parameter) of the incoming P wave can cause incoherence in stacking26,55. To avoid this, 4 different receiver function stacks were computed for each station by binning individual receiver functions around central ray parameter values of 0.045 s km−1, 0.055 s km−1, 0.065 s km−1 and 0.075 s km−1 (Extended Data Fig. 2c). To image details in lithospheric structure, 2 overlapping Gaussian filters (1.0 and 2.5; Extended Data Fig. 2c) were used within each of the ray parameter bins to allow for resolution of sharp versus gradational discontinuities26,56.
Receiver functions recorded by individual TRAILS stations were found to have high cross-correlation coefficients18, reflecting a lack of backazimuthal variation in structure associated with seismic anisotropy and/or heterogeneous structure (Extended Data Fig. 1a,c). Thus, additional subdivision in backazimuthal bins was not deemed necessary in this study. Tangential receiver functions were also computed to assess the potential influence from anisotropy and/or small-scale three-dimensional crustal heterogeneity26 below each station (Extended Data Fig. 1b,d). Although notable tangential energy is present in the melt-rich Afar Depression57, no notable tangential amplitudes were observed at any station below our study area (Extended Data Fig. 1b,d). An isotropic, laterally homogeneous layered structure beneath the TRAILS network and surrounding regions therefore suffices to explain the main features in the receiver function waveforms.
Of the 38 seismograph stations analysed, 26 produced reliable receiver functions. Seismograph stations located in thick sedimentary basins (for example, those nearest to Lake Turkana: BUBE, KALK, OMOE and TBIK) produced P-to-S converted energy that masked Moho arrivals, in most cases resulting in delayed receiver function P and PS phases18,58. Other stations situated directly on basalts which overlie 2–3-km-thick sediments (for example, in the Anza Basin: BASK, BOBE, KRGK and MAIK) also failed to produce receiver functions suitable for joint inversion18.
Joint inversion procedure
... continue reading