The Dynamic Crust of Northern Afar and Adjacent Rift Margins: New Evidence From Receiver Function Analysis in Eritrea and Ethiopia

Afar is undergoing the final stages of continental rifting and hosts the triple junction between the Red Sea, Gulf of Aden, and Main Ethiopian rifts. To better understand the nature of the crust and continental breakup in the region, we calculate teleseismic receiver functions across northeastern Afar and the Danakil microplate, using new data from a regional deployment in Eritrea. We estimate the Moho depth and bulk crustal VP/VS ratio using the H‐κ stacking method. The heterogeneity of our crustal thickness estimates (∼19–35 km) indicates that the Danakil microplate has undergone stretching and crustal thinning. By investigating the relationship between crustal thickness and topographic elevation, we estimate the regional crustal bulk density as ρc ≈ 2,850 ± 20 kg m−3, which is higher than expected, given the crustal thickness of the region. We show that topography is 1.5 ± 0.4 km higher than would be expected due to crustal isostasy alone. We propose that this topography is supported by the same hot mantle upwelling suggested to be responsible for the onset of rifting in East Africa. Uplift is generated due to the presence of a hot thermal anomaly beneath the plate and by thinning of the lithospheric mantle. Our results are consistent with a number of independent constraints on the thermal structure of the asthenospheric and lithospheric mantle. Evidence of melt within the crust is provided by anomalously high VP/VS ratios of >1.9, demonstrating that magma‐assisted extension continues to be important in the final stages of continental breakup.


Introduction
The Afro-Arabian Rift System is an ideal setting to investigate the temporal and spatial evolution of continental rifting.It consists of a network of rifts that developed asynchronously and are currently at different stages of the rifting process.Initial continental breakup is occurring in the south of the EARS in Mozambique, while further to the north, full seafloor spreading has been established in the Gulf of Aden and the Red Sea (e.g., Manighetti et al., 1998;Tesfaye et al., 2003).Elsewhere, the Main Ethiopian Rift displays the characteristics of a continental rift in its immediate stage of development, while Afar is undergoing the final stages of continental breakup as it develops into oceanic rifting in the transitional Red Sea zone (e.g., Kogan et al., 2012;Makris & Ginzburg, 1987).It is currently a matter of debate whether full seafloor spreading has begun on land in Afar.Magnetic anomalies have been identified in the region, analogous to those observed at oceanic spreading centers (e.g., Courtillot et al., 1980), which could indicate the onset of seafloor spreading, or could be explained as a result of magmatism during late stages of continental rifting (Bridges et al., 2012).
The Afar region is a diffuse extensional province which hosts a triple junction formed by the meeting of three diverging plate boundaries: the Main Ethiopian Rift, the Gulf of Aden, and the Red Sea, separating the Nubian, Somalian and Arabian tectonic plates respectively (Figure 1).The Gulf of Aden and the Red Sea do not directly connect in Afar, but are actively propagating through the region via the subaerial development of a series of disconnected, propagating rift segments (e.g., Manighetti et al., 1998Manighetti et al., , 2001)).Instead of directly connecting to the Gulf of Aden through the Bab el-Mandeb Strait, the Red Sea bifurcates into two branches south of 16°N (Barberi & Varet, 1977;McClusky et al., 2010).This process has isolated the Danakil microplate (or Danakil block), a sliver of continental crust rotating independently of Nubia (Figure 1a; Courtillot et al., 1984;Manighetti et al., 2001;Sichler, 1980).
The coherent rotation of the Danakil block has been established by several geodetic models and is associated with the formation of the Danakil Depression, a region of very low topography in northern Afar (e.g., Eagles et al., 2002;McClusky et al., 2010;Viltres et al., 2020), but there are fewer measurements around the proposed borders of the microplate, meaning the exact boundaries are not well-constrained by the models.Continental deformation models have two end-members: the block model, where deformation is confined to major block-bounding faults, with little internal deformation, and the continuum model, where deformation is diffused and quasi-continuous across a wider region (Thatcher, 2007).The transition between block and continuum models is gradational and has many intermediate stages.Thus, despite the observed coherent rotation of the Danakil block, it may have undergone internal deformation and stretching.Kinematic models have produced different estimates for the extent of the Danakil block.Viltres et al. (2020) highlight the diffuse character of the boundary between the Danakil and Nubia in the northeast, where inter-rifting deformation extends over a wide deformation zone.Varet (2017) place the southern boundary of the Danakil microplate near the Nabro Volcanic Range.This boundary choice is supported by the distribution of regional seismicity along this trend (Goitom, 2017).As a result of geographic remoteness, very few studies have been undertaken in northeastern Afar and the Danakil block, meaning there is a lack of information on the microplate and its role in partitioning strain between the two branches of the Red Sea Rift.
The nature of the crust beneath Afar has been widely debated.Some authors argue that new oceanic crust is being created in the region (e.g., Mohr, 1989), whilst others conclude that the area is entirely underlain by intruded, stretched continental crust (e.g., Makris & Ginzburg, 1987).An alternative proposal is that the overall crustal composition is transitional between continental and oceanic (e.g., Hammond, Kendall, et al., 2011;Redfield et al., 2003).These studies suggest that the region west of the current rift axis is predominantly made up of mafic material with oceanic affinities, whereas continental material exists to the east, potentially connected to the Danakil block.
Continental breakup involves extensional faulting, ductile plate stretching and intrusive and extrusive magmatism (e.g., Barton & White, 1997;Bastow & Keir, 2011;White & McKenzie, 1989).It is crucial to elucidate the extent to which mechanical-or magma-assisted extension processes dominate at each stage of continental breakup as it transitions into seafloor spreading.Indeed, this transition process may be more complex than a smooth evolution from mechanical-to magma-dominated extension into oceanic spreading.Previous studies in the Main Ethiopian Rift have shown that upper-crustal extension throughout much of Ethiopia is currently accommodated primarily by dyke intrusion along narrow, magmatic segments; strain has migrated away from the border faults and localized to rift-aligned magmatic intrusions since around 20 Ma (e.g., Ebinger & Casey, 2001;Wolfenden et al., 2004).However, the Danakil Depression to the north, which has been affected by magmatism since rift onset, now hosts significantly thinned continental crust, coupled with marked subsidence and thick evaporite deposition (Hutchinson & Engels, 1972;Keir et al., 2013;Wolfenden et al., 2005).The seismic reflection study of Hurman et al. (2023) found that upper crustal axial faulting in a magma-poor area of the depression accommodates 30% of the total extension.Thus it has been proposed that the final stages of continental breakup in northern Afar may be characterized by a return to ductile plate stretching and upper-crustal brittle faulting (Bastow et al., 2018).Hammond, Kendall, et al. (2011) comment that the lack of data from Eritrea limits overall regional interpretations of crustal structure.Eritrea hosts the Eritrean plateau (an extension of the western Ethiopian plateau), the Gulf of Zula, the Danakil Depression, the Danakil Alps and substantial off-rift volcanism in the Nabro Volcanic Range.The strike of the Nabro Volcanic Range is almost parallel to the trend of the Hanish-Zukur islands in the Red Sea and perpendicular to regional Red Sea rift trends.The presence of these volcanic centers across the Danakil block has led to proposals that they represent the surface expression of a "leaky" transform fault accommodating strain transfer from the Red Sea into Afar (e.g., Barberi & Varet, 1977;Viltres et al., 2020).
Crustal thickness (i.e., depth to the base of the Mohorovičić discontinuity) and the ratio of P-wave velocity to S-wave velocity (V P /V S ratio) provide key constraints on crustal composition and the amount of partial melt hosted by the crust (e.g., Dugda et al., 2005;Stuart et al., 2006;Zandt & Ammon, 1995).Previous studies of crustal structure in the region have shown that the crustal thickness varies from 40 to 45 km beneath the western Ethiopian plateau, decreasing to 20-26 km in central Afar, and becoming as thin as 16 km in the Danakil Depression (Ahmed et al., 2022;Berckhemer et al., 1975;Cornwell et al., 2010;Dugda et al., 2005;Hammond, Kendall, et al., 2011;Mackenzie et al., 2005;Maguire et al., 2006;Makris & Ginzburg, 1987;Reed et al., 2014;Stuart et al., 2006).The V P /V S ratio also shows spatial variation, from 1.7 to 1.9 beneath the plateau to >2.0 near magmatic segments where the crustal thickness is less than 26 km; the high values of V P /V S are interpreted as the result of aligned melt stored in interconnected stacked sills in the lower crust (Hammond, 2014).
Finally, it is widely accepted that a large low velocity province exists within the mantle beneath the African plate (e.g., Koelemeijer et al., 2017;Ritsema et al., 1999).It has been long-argued that this long-lived thermo-chemical anomaly has provided the impetus necessary to incite rifting processes in the region, as well as inducing mantle flow that dynamically supports the excess elevation of Africa (Ebinger & Sleep, 1998;Kendall & Lithgow-Bertelloni, 2016;Lithgow-Bertelloni & Silver, 1998).However, the exact form of mantle upwelling from the core-mantle boundary to the surface is still a matter of debate, particularly in northeast Africa (e.g., Civiero et al., 2015).Hence the regional effect of a sub-crustal thermal anomaly upon northeastern Afar and the Danakil block merits investigation.The primary controls on topographic elevation are the isostatic effects of variations in crustal thickness and density (e.g., Tiberi et al., 2005); however, the density structure of the mantle, particularly that of a hot thermal anomaly, can affect the pattern of topography observed at the surface (e.g., Bertelloni & Gurnis, 1997;Bird, 1978;Hager & Richards, 1989;Stephenson et al., 2021).By examining the relationship between topographic elevation and crustal structure in the study region and comparing it to the global relationship between these parameters, we can deduce information about mantle structure in northeastern Afar and the plateaus within the Danakil block.
In this paper, we exploit a new data set from a temporary regional seismic deployment, the Eritrea Seismic Project, to calculate teleseismic receiver functions.Using the H-κ method developed by Zhu and Kanamori (2000) and modified by Ogden et al. (2019), we calculate the crustal thickness (H) and bulk V P /V S ratio (κ) at each station.We supplement our analysis with data from another temporary deployment in central Afar, AFAR0911.The use of the AFAR0911 deployment enables us to update previous work by Hammond, Kendall, et al. (2011), Reed et al. (2014), and Ahmed et al. (2022), as well as validating the method of Ogden et al. (2019).Finally, we investigate the relationship between elevation and crustal thickness in order to understand the dominant control on topography within the region.

Seismic Network
To calculate receiver functions, we use seismic data from the Eritrea Seismic Project, a temporary regional seismic network deployed across the country over a period of 16 months from June 2011 until October 2012 (green inverted  (Hammond, Kendall, et al., 2011) and inverted orange triangles (Ahmed et al., 2022).Volcanoes from the Smithsonian Institution Global Volcanism Program are plotted as red triangles (Smithsonian, 2023).Inset shows the locations of the teleseismic earthquakes used in this study, with a yellow star marking the study area.triangles in Figure 1b; Hammond et al., 2011a).Six broadband seismic systems were provided by SEIS-UK (five CMGESP and one CMG3T).Stations were placed on the Eritrean plateau (EITE, CAYE), near the Gulf of Zula (DOLE), at the foot of the Danakil Alps (TIOE), near the Nabro Volcanic Range (FAME) and close to the Bab el-Mandeb strait (ASSE).
We supplement seismic data from this regional network with data collected by a temporary local seismic network deployed in the aftermath of Nabro volcano's 2011 eruption (plotted collectively as station NAB in Figure 1b; Hammond et al., 2011b).Eight 3-component broadband 30 s Güralp seismometers (five CMG-6TD and three CMG-40TD) were provided by SEIS-UK to monitor Nabro's post-eruptive state.The network was fully operational from 31 August 2011 until October 2012.The data were all initially recorded at 100 Hz sample frequency and then switched to 50 Hz sample frequency early in October 2011.We also calculate receiver functions for stations from the AFAR0911 deployment, which ran from 2009 until 2013 (yellow inverted triangles in Figure 1b; Keir & Hammond, 2009).These stations are three-component broadband 60 s CMG-ESP Güralp seismometers, recording at 50 Hz sample frequency.Due to the lack of cultural noise, all the data exhibit a generally high signalto-noise ratio, with excellent data quality (Hammond et al., 2014).

Teleseismic Data Set
We search the USGS earthquake catalog for seismic events with magnitude M b ≥ 5, occurring at an epicentral distance greater than 30°.Events are then sorted into two epicentral distance ranges, 30°-90°and >60°.The former range is standard for calculating P-to-S receiver functions, as it eliminates the effects of core phases and upper mantle triplications.To increase our azimuthal coverage, we also calculate PP-to-S receiver functions from events in the distance range >60°.
The data are filtered using a Butterworth bandpass filter with corner frequencies 0.04 and 3 Hz.Seismograms are inspected manually to remove poor quality traces, and then analyzed based on the signal-to-noise ratio.Sourcereceiver pairs are kept for further analysis if the P-wave signal-to-noise ratio is >3.This process results in the selection of 102 earthquakes for receiver function construction and processing (inset in Figure 1b).

Receiver Function Calculation
Receiver function analysis is a well-established technique for investigating velocity contrasts in the crust and upper mantle beneath three-component seismic stations (Langston, 1979).When a P-wave impinges on a velocity contrast, part of the energy is converted into an S-wave.In the absence of anisotropy and dipping boundaries, these P-to-S conversions are preferentially recorded on the radial components of a seismogram, whilst direct Pwave arrivals are recorded on the vertical component.For teleseismic waves, the incidence angle of the incoming wave will be close to vertical (27°or lower), meaning that we can use the vertical component of motion as a good approximation of the event source.Deconvolution of the vertical component from the radial components removes information common to both (i.e., the instrument response, path and source effects).This isolates the response from the local crustal structure along the incoming, near-vertical raypath from a teleseismic earthquake.The resulting receiver functions contain phase arrivals that represent the P-to-S conversion and subsequent reverberations within a layer (the positive polarity PpPs phase, and the negative polarity co-arriving PpSs + PsPs phases, see Figure 3 in Jenkins et al. (2020)).The amplitudes of the arrivals depend on the P-wave incidence angle (i.e., the ray parameter), the size of the velocity contrast, and whether the wave passes from a high velocity layer to lower velocity layer, or vice versa.The arrival times of the converted phase and reverberations depend on the depth of the velocity contrast, the P-wave and S-wave velocity between the contrast and the surface, and the ray parameter.
We calculate receiver functions using the iterative time-domain deconvolution method of Ligorría and Ammon (1999).This technique iteratively builds receiver functions in the time domain from Gaussian pulses of a set width.For our analysis, we use a range of Gaussian width factors at 0.2 intervals from 0.8 to 4.0 to produce a set of 17 individual receiver functions for each source-receiver pair.Since periods less than the Gaussian width are not resolvable, our calculated receiver functions have corresponding low-pass filters between 0.4 and 2.0 Hz.The radial receiver function is then convolved with the vertical component and cross-correlated with the original radial component to evaluate the iterative deconvolution variance.If this variance is ≤80%, the receiver function is rejected.We also undertake a final manual inspection of the receiver function traces, removing poor quality traces (e.g., those where the maximum peak arrival is not the direct P-wave arrival, those with anomalous long wavelength features, and those which are significantly different from others that sample the same region).This results in a range of 9-63 receiver functions per station (Table 1).

H-κ Stacking
We use the H-κ stacking procedure of Zhu and Kanamori (2000) to estimate bulk crustal properties beneath each seismic station (see Text S1 in Supporting Information S1 for further details on the technique).This method has produced accurate estimates of these properties in regions where the crustal structure is simple (e.g., Thompson et al., 2010).However, our study region is tectonically complex region, and we therefore anticipate some potential difficulties in receiver function calculation and interpretation.
The presence of thick sedimentary layers in the study area has been well-established (see Rime et al. (2023) for a comprehensive geological map).The Danakil Depression hosts a sedimentary basin likely containing up to ∼3-5 km of Pliocene-Pleistocene-Recent evaporites (Bastow & Keir, 2011;Bastow et al., 2018;Hutchinson & Engels, 1972;Rime et al., 2023).Extensional basins in southern and central Afar are filled with Pliocene-Holocene lacustrine, alluvial and volcaniclastic sediments, exceeding 200 m in places (e.g., Renne et al., 1999).Although parts of southern and central Afar are covered by the Stratoid basalts, this cover is lacking to the north in the Danakil Depression and on the Danakil microplate, which are dominated by sedimentary rocks (Barberi et al., 1980;Rime et al., 2023;Tortelli et al., 2022;Varet, 2017).Shallow sedimentary layers produce significant P-to-S converted energy that can mask the signal from the Moho and its crustal reverberations, increase the time delay of the direct P-arrival, and enhance the "ringing" nature of the receiver function signal (e.g., Ogden et al., 2019;Reed et al., 2014;Zelt & Ellis, 1999).
Due to the complicated tectonic history of Afar and the Danakil block, we expect significant local differences in the amount of extension and/or thickening due to underplating and intrusion.This can be seen in the variability in previous Moho depth estimates across Afar (e.g., Ahmed et al., 2022;Dugda et al., 2005;Hammond, Kendall, et al., 2011;Lavayssière et al., 2018;Reed et al., 2014).It is distinctly observable at single stations with backazimuthal variation in Moho depth, due to the receiver functions from different directions sampling crust of different thickness on either side of the station (e.g., Dugda et al., 2005;Hammond, Kendall, et al., 2011).This is particularly pronounced for stations in Afar that are close to the large border faults, such as ABAE (Hammond, Kendall, et al., 2011).Regions of lower crustal intrusions show a gradational transition between the crust and mantle, as observed by Lavayssière et al. (2018) beneath the western Ethiopian plateau.A gradational Moho lowers the amplitude of the Moho P-to-S conversions (Gallacher & Bastow, 2012).The frequency content of the receiver functions thus acquires greater significance for these low amplitude P-to-S conversions, as H-κ stacking will only be sensitive to the Moho using lower frequency receiver functions if the Moho is a more gradational boundary (e.g., Ogden et al., 2019).
For H-κ stacking, V P is held constant for the whole crust and has to be either known a priori or assumed.Previous receiver function studies in southern and central Afar use a wide range of values for average crustal P-wave velocity, from 4.65 to 6.5 km s 1 , based on criteria such as wide-angle seismic reflection studies or by examining the theoretical relationship between partial melt fraction and seismic velocities.Some authors allow for variability in V P between stations (e.g., Reed et al., 2014), whilst others assume fixed values for all stations in their study (e.g., Ahmed et al., 2022;Dugda et al., 2005;Hammond, Kendall, et al., 2011).
To produce robust estimates, we therefore employ the cluster analysis approach of Ogden et al. (2019).This methodology rigorously explores the H-κ parameter space, including frequency content of the receiver functions and assumed average crustal V P , to account for sensitivity of the method to complex crustal structure.The quality control criteria applied to the receiver functions and the cluster analysis provides quantitative evidence for the reliability of H-κ stacking at a station.The number of quality control criteria failed by a station gives an indication of why H-κ stacking has produced an unsatisfactory result.In regions where the Moho is gradational, the ability to adjust the frequency content of the receiver functions means that H-κ stacking remains reliable.The final result will successfully identify the center of the Moho depth range if the Moho layer thickness is in the range of 0-13 km (Ogden et al., 2019).See Supporting Information S1 (Text S2) for more detail on the procedure and the associated uncertainty estimates.

Results
Seven stations from the Eritrean Seismic Project, including the composite station NAB (made up of the six Nabro array stations), and seven stations from the AFAR0911 network produce reliable results using our H-κ stacking procedure (Table 1, Figure 2).Ten quality control criteria assess the quality of these results, and are described in the Supplementary Information (Text S3 in Supporting Information S1 for the Eritrean Seismic Project, Text S4 in Supporting Information S1 for AFAR0911).Diagnostic result figures for each station are plotted in Figures S1-S19 in Supporting Information S1.One of the stations on the Eritrean plateau (EITE) shows strong variation in the Ps arrival time with back-azimuth, a phenomenon previously observed in the region by Dugda et al. (2005) and Hammond, Kendall, et al. (2011).Following these examples, we split the receiver functions at EITE into those arriving from the east and those arriving from the west, carrying out H-κ stacking on each set separately.These results are plotted as EITE_E and EITE_W, respectively, at the first Moho bounce point of the PpPs multiple for the dominant back azimuth (Figure 2).We follow a similar procedure for the central Afar stations ABAE and MAYE, as the receiver functions are almost exclusively from the east for these stations (Figures S11 and S18 in Supporting Information S1).In total, we calculate 15 estimates of crustal thickness and Almost all of the receiver functions show an obvious Ps Moho arrival between 2.5 and 5 s, corresponding to a crustal thickness range of 18.9-34.8km across all stations.This arrival is not always the first positive arrival, which is evidence of a sedimentary layer at the surface.In general, the receiver functions show additional P-to-S arrivals as well as the Moho crustal reverberations, indicating the presence of complex intracrustal structure leading to these subsidiary conversions.The PpPs Moho arrival is observable in most of the receiver functions between 8.5 and 15.5 s, whereas the final arrival of PpSs + PsPs is more ambiguous, occurring between 12 and 20 s.The stations show a range of V P from 5.9 to 6.5 km s 1 , which we determine to be reasonable, given previous controlled source studies (Mackenzie et al., 2005;Maguire et al., 2006;Makris & Ginzburg, 1987).
Within error, the crustal thicknesses we observe are consistent with previous estimates by Hammond, Kendall, et al. (2011) although they tend to estimate a slightly higher V P /V S than our results.This consistency gives confidence to the H-κ stacking procedure we employ, as well as our selected values of V P .

Stretching of the Danakil Block
It is well-established that continental basement is exposed on the Danakil microplate, particularly within the Danakil Alps (e.g., Bosworth et al., 2005;Freund, 1970;Rime et al., 2023;Varet, 2017).Drury et al. (1994) propose that this crustal segment was shortened and thickened by transpression in late Neoproterozoic to Cambrian times.Our crustal thickness estimates for the stations located on the Danakil block are lower than the global average value for extended continental crust of 30.5 km (Christensen & Mooney, 1995).They are also 10-15 km thinner than estimates for stations found on the western Ethiopian plateau (Figure 2).Indeed, stations TIOE and FAME show thicknesses that are similar to those found in the rift valley of Afar.While we lack constraints from the interior of the Danakil Alps, we nevertheless propose that the Danakil is not "relatively unstretched" as has been previously claimed (Wolfenden et al., 2005), but has undergone a significant amount of deformation and crustal thinning.Two-dimensional area balancing in a previous plate kinematic reconstruction study supports this conclusion, suggesting that the Danakil block has undergone a minimum of 200% stretch since the onset of rifting (Figure 3; Redfield et al., 2003).A geometric estimate by Morton et al. (1975) also gave a stretching factor (i.e., ratio of stretched width to original width) of up to 2.5 in the Danakil Alps.
Station EITE (E) has a V P /V S ratio of 2.054, which is significantly higher than the average value for continental crust of 1.768 (Christensen, 1996).V P /V S ratios above 1.9 cannot be accounted for by differences in rock composition (Thompson et al., 2010).Instead, these elevated values are typically explained by the presence of fluid (e.g., melt) in the crust, which decreases the S-wave velocity without significantly affecting the P-wave velocity (e.g., Watanabe, 1993).This effect is enhanced in the presence of aligned melt due to the influence of anisotropy (Hammond, 2014).EITE (E) is located on the edge of the western plateau.It shows strong similarity to stations MAYE and ABAE, which are positioned on the rift shoulder in central Afar and also show high V P /V S ratios.We therefore propose that partial melt in the crust is present along the edge of the western plateau in central Afar and persists northward.
Another very elevated V P /V S value of 2.042 is observed at station TIOE.In general, crustal thickness in the region is negatively correlated with the bulk V P /V S ratio (see Figure S20 in Supporting Information S1), implying that the thinnest crust is also the most intruded.This indicates that the stretched, thinned nature of the crust beneath the Danakil block has resulted in the presence of partial melt and focused magmatism, similar to what is observed by Dixon et al. (1989) in the Red Sea.It is possible that this melting and volcanism was what initially caused the onset of crustal extension, as postulated by Drury et al. (1994).Our observations suggest that it is likely that magmatic extension continues to play a dominant role in the final stages of continental breakup in Afar.Stations ASSE, CAYE, and FAME have V P /V S values that are greater than 1.8, which are closer to the expected values for oceanic crust (Christensen, 1996).The most likely explanation for this is the presence of magmatic intrusions causing the continental crust to be more mafic in composition.

Crustal Sills Beneath Nabro Volcano
Our V P /V S result at Nabro volcano (station NAB, Figure 2) is relatively low (1.773).Nabro erupted in 2011 and therefore the crust beneath the volcano might be expected to have an elevated V P /V S ratio, reflecting melt storage.However, it is important to consider the pore geometry of melt storage and its influence on seismic velocities (Takei, 2002).Hammond (2014) demonstrated that the relative orientation and alignment of melt has a significant effect on H-κ estimates.If melt is vertically aligned, the shear wave conversion will be split as it passes through, with each of the split waves producing different κ estimates due to their different velocities.When there is good back-azimuthal coverage at a station, the overall stacking procedure will be dominated by the slow shear wave velocity, and thus produces a higher V P /V S estimate.This effect becomes less pronounced as melt becomes more horizontally aligned.The petrological study of Donovan et al. (2018) concluded that Nabro hosts distinct batches of magma in sills throughout the crust, which might go some way toward explaining the low bulk V P /V S ratio estimate produced by H-κ stacking.
Other geophysical studies have proposed crustal magma reservoirs taking the form of sub-horizontal lenses.
Ground deformation observed at the Rhenish Massif in Central Europe has been explained as melt accumulation in sub-horizontal sill-like structures (Silverii et al., 2023).A seismic reflection study of the Juan de Fuca Ridge concludes that sub-horizontal stacked magma lenses are present in multiple locations beneath the ridge, providing sources for eruption events (Carbotte et al., 2021).The volcanic plumbing system around Mayotte Island includes a subhorizontal magma reservoir lying between a new seamount and Mayotte, as suggested by the deformation study of Cesca et al. (2020).A local velocity model for Mayotte found an average V P /V S ratio of 1.6, an even lower value than our estimate for Nabro, which the authors attribute to the presence of low aspect-ratio cracks filled with gas-rich magma (Lavayssière et al., 2022).2014), Ahmed et al. (2022).Elevation estimates are from the Generic Mapping Tools Global Earth Relief grid, with the topography low-pass filtered for wavelengths >30 km (Wessel et al., 2019).The black line is the best-fitting regional isostatic relationship where (ρ a ρ c )/ρ a = 0.11; the dotted black line is the average global isostatic relationship (Stephenson et al., 2024).(b) Sub-crustally supported topography as a function of asthenospheric mantle temperature and lithospheric thickness, assuming a temperature anomaly is constrained to a 100 km thick channel and equilibrium lithospheric thickness is 120 km.Contours indicate regional elevation anomaly, ϵ, generated by lithospheric thickness and sub-plate temperature.The shaded area represents the range of regional sub-crustal support (Δϵ = 1.5 ± 0.4 km) estimated from the offset in intercept values between the regional crustal thickness-elevation relationship and the global one.The red square is plotted at the subplate temperature and lithospheric thickness given by the tomography model of Celli et al. (2020) after conversion into temperature following the methodology of Hoggard et al. (2020) and Richards et al. (2020).The yellow star is plotted at the sub-plate temperature and lithospheric thickness given by inverse modeling of rare earth element compositions (Ball et al., 2021).
However, the low V P /V S ratio cannot be solely explained by horizontally aligned melt.Alternatively, the overall crust at Nabro could be melt-depleted, with melt traveling directly up from the mantle to shallow storage regions just prior to eruption.Nabro's high crustal thickness could be due to these magmatic intrusions and additions occurring through time.There is petrological evidence for a significant amount of old eruptive products in the crust beneath Nabro, due to the presence of xenocrystic material in the erupted magmas Donovan et al. (2018).
Unfortunately, the back-azimuthal coverage is not comprehensive enough to carry out a crustal anisotropy study at Nabro to investigate these proposals further (e.g., Hammond, 2014;Liu & Niu, 2012).

Implications for Rifting Through the Gulf of Zula
Afar is a diffuse extensional province; Courtillot et al. (1980) note that it is difficult to identify plate boundaries in Afar due to the homogeneous stretching of the region.Manighetti et al. (2001) find that strain within Afar is localized in disconnected rift segments, currently opening and propagating through numerous active faults, open fissures and volcanism.The active rift axis passes from the Red Sea through the Gulf of Zula to the west of the Danakil block and becomes more distributed further south in Afar, eventually spanning distances of ∼175 km (Figure 26 in Manighetti et al., 2001;Kogan et al., 2012).The rotation of the Danakil block toward Arabia favors this propagation of rifting on land in Afar rather than a southward propagation of the Red Sea Rift beyond 16°N (Sani et al., 2017).
We note that station DOLE, despite its proximity to the rift axis through the Gulf of Zula, has a relatively thick crust of 27.2 km.However, it is important to notice that the receiver functions mainly sample the crust to the north-east of the station (and thus further from the rift), due to the back-azimuthal coverage.Therefore, the crustal thickness recorded at this station likely reflects a portion of unstretched continental crust on the Danakil microplate.
Viltres et al. ( 2020) conclude from their GPS study that inter-rifting deformation is not localized directly in the Gulf of Zula.Instead it is accommodated by volcanic vents and distributed faults 20 km west of the Gulf of Zula.Station EITE (E) lies to the west of this area of active deformation, with a crustal thickness of 25.5 km and a relatively high V P /V S ratio of 2.05.We postulate that extension may continue to shift westwards if the crust is rheologically weak and magma-rich at the edge of the plateau.

Topographic Elevation and Crustal Structure
Topographic elevation can primarily be attributed to the isostatic effects of crustal thickness and density variations (e.g., Airy, 1855;Pratt, 1855).However, a secondary control on continental topography is the density structure of the mantle (e.g., Hager et al., 1985;Hoggard et al., 2021;Lithgow-Bertelloni & Silver, 1998;Sembroni et al., 2016;Stephenson et al., 2021).Our results can therefore place indirect constraints on crustal and mantle density structure.We investigate the relationship between topographic elevation (ϵ) and crustal thickness for central and northeastern Afar and the Danakil block (Figure 3a).
There is a demonstrable linear relationship between H and ϵ, suggesting that 3 km of variation in the topographic elevation in the study area can be ascribed to variation in crustal thickness (Figure 3a).A similar relationship between topographic elevation and crustal thickness is widely accepted from global observations (see Figure 5 in Stephenson et al., 2024;Lamb et al., 2020 and plotted in Figure 3a for comparison).However, two characteristics of the distribution of H and ϵ in our regional study provide a rare opportunity to quantify the regional properties of the crust and uppermost mantle.First, the slope of the relationship between H and ϵ is on the lower end of globally averaged estimates (e.g., Lamb et al., 2020).This difference is not unexpected given the rift setting of the region since the value of the slope is largely controlled by the density of the crust, although it has seldom been quantified (e.g., Ebinger et al., 1989;Tiberi et al., 2005).Furthermore, the intercept crustal thickness value is shifted upwards, which implies significant levels of sub-plate topographic support which must arise from processes operating in the mantle beneath the region.In other words, topographic elevation is generally higher than would be expected for a given value of crustal thickness if crustal isostasy were the only control on topography (Airy, 1855;Pratt, 1855).Here, we use these observations to provide new constraints on lithospheric and subplate density structure.

Bulk Crustal Density
We can quantify the difference in average crustal density between our study region and the global average by using the gradient of the plot in Figure 3a to estimate crustal density beneath central and northeastern Afar and the Danakil block (e.g., Lamb et al., 2020;Stephenson et al., 2024).Assuming Airy isostatic equilibrium, the slope of the relationship between crustal thickness, H and topographic elevation, ϵ, is given by on the ratio (ρ a ρ c )/ρ a , where ρ c is crustal density and ρ a is asthenospheric mantle density.This relationship further assumes that topography is compensated in the asthenospheric mantle and that lithospheric mantle structure and crustal density are both consistent across the region of interest (cf.Lamb et al., 2020).Hence if we can estimate the local value of ρ a , it is possible to place a bound on the average density of the crust in our study region.
To estimate the value of ρ a , we must calculate the effects of temperature and pressure upon mantle rocks given local constraints on asthenospheric temperature structure.Mantle density as a function of pressure, P, and temperature, T, is given by where α = (3.28 ± 0.25) × 10 5 K 1 is the thermal expansivity, K = 140 ± 20 GPa is the bulk modulus of the mantle and ρ m • = 3330 kg m 3 is mantle density at surface pressure and temperature.In this way, we place constraints on T and P in the asthenosphere beneath our study region.Ball et al. (2021) use inverse modeling of rare earth element concentrations in mafic volcanic rocks to estimate that the sub-plate potential temperature beneath north-eastern Ethiopia is ∼1,442 ± 50°C and lithospheric thickness a 1 = 30 ± 20 km.Shear wave velocity can also be converted into temperature using the experimentbased parameterization of Yamauchi and Takei (2016), as carried out by Priestley and McKenzie (2013), Richards et al. (2020), andHoggard et al. (2020).This approach minimizes the weighted misfit between: predicted and observed values of seismic velocities as a function of oceanic plate age and depths ≤112.5 km; the predicted shear wave velocity as a function of the mantle isentrope between 225 and 400 km depth; the globally averaged shear-wave attenuation as a function of depth; and finally published and calculated values of mantle viscosity (Hoggard et al., 2020;Richards et al., 2020).When applied to the recent tomographic model of Celli et al. (2020), this approach yields an asthenospheric potential temperature of T P ≈ 1,440 ± 30°C.To extract a lithospheric thickness, this temperature model is contoured at the 1,175°C isotherm, which approximates the temperature at which convection begins.We find that the lithosphere beneath Afar and the Danakil microplate is between ∼40 and 50 km thick.To ensure that we place the depth of isostatic compensation within the asthenospheric mantle, we choose a depth to the base of the lithosphere of 50 ± 20 km.Assuming that the mantle is isentropic, we assume a mantle potential temperature beneath the study region of T P = 1,440 ± 50°C.Hence, using we obtain asthenospheric temperature at this depth of T a = 1,465 ± 52°C, given acceleration due to gravity, g = 9.81 m s 2 and heat capacity, C p = 1,187 ± 50 J°C 1 (Mckenzie & Bickle, 1988).
These values for asthenospheric temperature and lithospheric thickness are corroborated by several other studies.Petrogenic modeling by Ferguson et al. (2013) finds a mantle potential temperature of 1,450°C, corroborated by the geodynamic models of Armitage et al. (2015).Stephenson et al. (2023) use major element thermobarometry of mafic magmas to estimate a lithospheric thickness of 40-60 km beneath northern Afar (see methodology of McNab & Ball, 2023;Plank & Forsyth, 2016).Using a similar approach as we have followed above, Fishwick (2010) also obtains a lithospheric thickness of about 50 km.Finally, Ebinger and Hayward (1996) find that the elastic thickness of the plate in the Danakil Depression is between 3 and 8 km.Although they obtain a higher value on the Ethiopian plateau of 56 km, this value is still lower than the stable cratonic elastic thicknesses of 64-79 km observed elsewhere on the African plate by Ebinger et al. (1989).These results suggest that the lithosphere is anomalously weak beneath our study region and therefore consistent with thin lithosphere and elevated subplate potential temperatures.
Geochemistry, Geophysics, Geosystems 10.1029/2023GC011314 We use the relationship P 1 = ρga 1 to find the pressure at the base of the modern-day lithosphere, where ρ ≈ 3000 kg m 3 is overburden density, and a 1 is depth to the base of the lithosphere.We find that P 1 ≈ 1.3 ± 0.5 GPa.Therefore ρ a ≈ 3,210 kg m 3 .In Figure 3a we apply a non-biased linear regression (i.e., a total least squares regression; see Text S4 in Supporting Information S1) that assumes uncertainty in both topography and crustal thickness (±30 m and ±3 km, respectively) and find that (ρ a ρ c )/ρ a = 0.11, which is on the lower end of the global average value of around 0.11-0.13(Lamb et al., 2020;Stephenson et al., 2024).Using our estimated value of ρ a and by propagating uncertainties in T a and P 1 , we find that in our study region, ρ c ≈ 2,850 ± 20 kg m 3 .In our study region, crustal thickness is 26 ± 10 km (i.e., median and interquartile range).Stephenson et al. (2024) provide a function to predict average bulk crustal density as a function of crustal thickness, which for this range yields a predicted mean crustal density of 2,800 ± 10 kg m 3 , which indicates that the crust in our study region is somewhat denser than the global average.We compare the result to the results of Tiberi et al. (2005), who carry out inverse modeling of gravity anomalies in Afar.In the upper 10-20 km beneath northern Afar, they find a large positive lateral density contrast compared to other areas of their study region.This large change in density is interpreted to be a result of a combination of crustal thinning and the addition to the crust of dense mafic material due to magmatic intrusions.
In summary, this somewhat elevated crustal density likely reflects the substantial magmatic intrusions and additions to the crust by mafic volcanic rocks that have taken place in this region.It furthermore provides a new and independent constraint that adds weight to our earlier conclusion that the high V P /V S ratios we observe are likely to be due to magmatic activity.Indeed, it is probable that these density estimates represent a lower bound for the thinnest crystalline crust in Afar and the Danakil microplate, much of which is overlain by thick sedimentary layers.In the Danakil Depression, the layers of sediment reach 5 km in extent (Bastow & Keir, 2011;Bastow et al., 2018;Hutchinson & Engels, 1972;Rime et al., 2023).If crystalline crust of thickness 15 km is overlain by 5 km of sediment of density 2,700 kg m 3 , and if the bulk density of both layers combined is 2,850 kg m 3 , it follows that the density of the crystalline crust is ∼2,900 kg m 3 .This value is consistent with significant mafic intrusion.

Sub-Plate Topographic Support
Next, we turn our attention to the topographic effects of elevated mantle potential temperatures beneath Afar and the Danakil block.It is well-established that mantle upwelling from a large low velocity province has had a significant effect on the tectonics and topography of Africa (e.g., Ebinger & Sleep, 1998;Kendall & Lithgow-Bertelloni, 2016;Nyblade & Robinson, 1994;Ritsema et al., 1999).While the magnitude of sub-crustal support from the mantle has been estimated in the Main Ethiopian Rift, a paucity of data means that the extent of the hot mantle swell northwards toward Afar and the Danakil block has not been well constrained (Holdt et al., 2022;Sembroni et al., 2016;Stephenson et al., 2024).From the relationship between crustal thickness and elevation, it is clear that topography is anomalously high given the crustal thickness in our study region.It is notable that, for H ≈ 10-20 km, the topography is approximately at sea level.However, globally, the average crustal thickness is about 32.1 km for regions at sea level (Lamb et al., 2020;Stephenson et al., 2024).To illustrate this point, take station CAYE, which has crustal thickness of H = 32.7 km, which is very close to the thickness of crust that is globally expected to reside at sea level.Instead, CAYE is located at an elevation of about ϵ = 2.1 km.This marked offset implies a significant level of regional sub-crustal support that sustains the topography in Afar and the Danakil microplate.
Using our results, we can for the first time quantify this dynamic support regionally.Some proportion of this locally anomalous topography will be supported by the flexural rigidity of the lithosphere.However, the size of our study region exceeds the expected elastic half-wavelength given an elastic thickness of 5-49 km that would result in regional elastic support (Ebinger & Hayward, 1996).It is therefore possible to assess the extent of regional excess topography by using the regression of H against ϵ.The difference between the mean intercept values of the local and global regressions of H against ϵ indicates that the study region is on average 1.5 ± 0.4 km (i.e., two standard errors on the mean) higher than would be expected as a result of crustal thickness variations alone (Figure 3a).
A hot mantle upwelling can generate uplift in three interconnected ways (Hager & Richards, 1989;Hager et al., 1985;Hoggard et al., 2021;Molnar et al., 2015).First, the upward flow of hot mantle produces viscous forces that cause vertical tractions on the base of the plate (e.g., Forte et al., 1993;Molnar et al., 2015).Second, the buoyancy of the hot mantle material generates isostatic uplift at the surface (e.g., Richards et al., 2020;Rudge et al., 2008).Finally, erosion of the base of the lithosphere by heat and/or melt and fluid infiltration can thin the lithospheric mantle (e.g., Bird, 1979).Here, we explore the effects of the second and third processes in our study region by estimating the uplift that is generated by the combined effects of a hot thermal anomaly beneath the plate and by thinning of the lithospheric mantle.
Uplift above a hot asthenospheric anomaly of thickness h and excess temperature ΔT is given by where a 0 and a 1 are initial and final lithospheric thicknesses and T 1 is ambient mantle temperature.We assume an initial lithospheric thickness of a 0 = 120 km, and take an asthenospheric channel to be h = 100 km thick (McKenzie, 1978;Rudge et al., 2008).The effects of using different values of a 0 and h are reported in in Text S5 and Figure S21 in Supporting Information S1.
Next, we return to the independent estimates of mantle temperature and lithospheric thickness determined from tomographic and geochemical constraints.We carry out a grid search by sweeping through the range of reported values for a 1 and ΔT using Equation 3 (see, e.g., Ball et al., 2021;Ebinger & Hayward, 1996;Hoggard et al., 2020;Stephenson et al., 2024).The residual topographic estimates calculated from our simple regression between H and ϵ in Figure 3a are consistent with the magnitude of uplift generated above a temperature anomaly of around 110°C residing in an asthenospheric channel of 50-150 km thickness beneath lithosphere that has been thinned by around 80 km (i.e., from 120 to 40 km; compare hatched area with circles and star on Figure 3b).These values are corroborated by the aforementioned independent estimates for sub-plate temperature and lithospheric thickness.It is important to emphasize that this signal of sub-plate topographic support is regional, extending beneath regions of both thick and thin crust, necessitating the emplacement of a hot sub-plate thermal anomaly and removal of the lithospheric mantle beneath both rift and plateau.
These results are consistent with the presence of a major mantle upwelling beneath Afar and the Danakil microplate.We note that our results do not attempt to estimate the mantle flow component of sub-crustally supported topography, nor do they attempt to quantify the component arising from below the uppermost 200 km of the mantle.Furthermore, a potentially significant component of topographic support from the lithospheric mantle arises from lithospheric extension.However, the mantle swell is clearly a widespread feature that extends from the Main Ethiopian Rift into northern Afar and the Danakil microplate and sub-plate temperature anomalies appear to be required to explain the height of topography in the region.

Conclusion
We have produced the first estimates of crustal thickness and bulk crustal V P /V S ratio in Eritrea, across an area encompassing northeastern Afar and the Danakil microplate.We calculate receiver functions from a temporary regional seismic array and apply the modified H-κ stacking technique of Ogden et al. (2019) to constrain the bulk crustal properties at eight locations in Eritrea.To validate this procedure, we also carry out H-κ stacking at stations in central Afar, producing estimates that are consistent with previous studies in the region.
This study has revealed the heterogeneity of crust in Afar and the Danakil microplate.We calculate a variability in crustal thickness of ∼16 km throughout the study area, ranging from ∼19-35 km.This reveals that the Danakil microplate is not an undeformed, unified block, but rather has experienced extension and crustal thinning resulting in crust as thin as 18.9 km.We also see variation in V P /V S ratio across our study region.Some stations have values typical of continental crust (∼1.77).Others show higher values of 1.8-1.86,likely a consequence of alteration by magmatic intrusion or addition to the crust.Two stations have V P /V S greater than 2.0, which we associate with the presence of partial melt, suggesting that magmatic extension continues to be important in the final stages of continental rifting.
By examining the relationship between topographic elevation and crustal thickness, we observe a deviation from the global relationship between these parameters.We estimate a denser bulk crust (i.e., including both crystalline basement and sediment) across northern and central Afar and the Danakil microplate (ρ c ≈ 2,850 ± 20 kg m 3 ), which is evidence for intrusive magmatism altering the composition of the crust to be more mafic.Our results reveal that the topography is anomalously high for the observed crustal thicknesses in Afar and the Danakil microplate, implying a level of sub-crustal support from a hot thermal anomaly beneath the crust.
Our observations of crustal properties address the previously existing data gap in northeastern Afar and the Danakil microplate, provide a useful base for further investigations of the crust and mantle within the Danakil microplate and at its boundaries, and give insights into the processes responsible for topographic elevation in the region.

Figure 1 .
Figure 1.Topographic map.The boundaries of the Danakil microplate and other plate boundaries are outlined by a dashed green line (Viltres et al., 2020).(a) Countries are labeled with black text and solid black lines demarcate national borders.Other geographic domains are labeled with red text.The navy box encloses the region plotted in (b).RS: Red Sea, GOA: Gulf of Aden, MER: Main Ethiopian Rift, NP: Nubian Plate, AP: Arabian Plate, SP: Somalian Plate.(b) Seismic stations used in this study are plotted as inverted green and yellow triangles and labeled.Seismic stations from previous studies are plotted as inverted blue triangles(Hammond, Kendall, et al., 2011) and inverted orange triangles(Ahmed et al., 2022).Volcanoes from the Smithsonian Institution Global Volcanism Program are plotted as red triangles(Smithsonian, 2023).Inset shows the locations of the teleseismic earthquakes used in this study, with a yellow star marking the study area.

Figure 2 .
Figure 2. Regional map showing the results of this study (Table1) alongside previous estimates of crustal thickness and V P /V S ratio (see TableS1in Supporting Information S1 for a complete list).Seismic stations from this study are labeled and plotted as black circles.Seismic stations from previous studies are plotted as black inverted triangles.The boundaries of the Danakil microplate and other plate boundaries are outlined by a dashed green line(Viltres et al., 2020).Dashed black lines demarcate national borders.Volcanoes from the Smithsonian Institution Global Volcanism Program are plotted as red triangles(Smithsonian, 2023).(a) The crustal thickness estimates.Thick black lines mark the profiles for controlled source work byBerckhemer et al. (1975),Makris and Ginzburg (1987),Egloff et al. (1991),Mackenzie et al. (2005),Maguire et al. (2006).Other estimates of crustal thickness come fromHammond, Kendall, et al. (2011), Reed et al. (2014), and Ahmed et al. (2022).(b) V P /V S ratio estimates from this study and fromHammond, Kendall, et al. (2011), Reed et al. (2014), and Ahmed et al. (2022).
Figure 2. Regional map showing the results of this study (Table1) alongside previous estimates of crustal thickness and V P /V S ratio (see TableS1in Supporting Information S1 for a complete list).Seismic stations from this study are labeled and plotted as black circles.Seismic stations from previous studies are plotted as black inverted triangles.The boundaries of the Danakil microplate and other plate boundaries are outlined by a dashed green line(Viltres et al., 2020).Dashed black lines demarcate national borders.Volcanoes from the Smithsonian Institution Global Volcanism Program are plotted as red triangles(Smithsonian, 2023).(a) The crustal thickness estimates.Thick black lines mark the profiles for controlled source work byBerckhemer et al. (1975),Makris and Ginzburg (1987),Egloff et al. (1991),Mackenzie et al. (2005),Maguire et al. (2006).Other estimates of crustal thickness come fromHammond, Kendall, et al. (2011), Reed et al. (2014), and Ahmed et al. (2022).(b) V P /V S ratio estimates from this study and fromHammond, Kendall, et al. (2011), Reed et al. (2014), and Ahmed et al. (2022).

Figure 3 .
Figure 3. Crustal and mantle structure.(a) Topographic elevation plotted against crustal thickness estimates (colored circles)for the study region, from this study and the previous receiver function studiesof Hammond, Kendall, et al. (2011), Reed et al. (2014),Ahmed et al. (2022).Elevation estimates are from the Generic Mapping Tools Global Earth Relief grid, with the topography low-pass filtered for wavelengths >30 km(Wessel et al., 2019).The black line is the best-fitting regional isostatic relationship where (ρ a ρ c )/ρ a = 0.11; the dotted black line is the average global isostatic relationship(Stephenson et al., 2024).(b) Sub-crustally supported topography as a function of asthenospheric mantle temperature and lithospheric thickness, assuming a temperature anomaly is constrained to a 100 km thick channel and equilibrium lithospheric thickness is 120 km.Contours indicate regional elevation anomaly, ϵ, generated by lithospheric thickness and sub-plate temperature.The shaded area represents the range of regional sub-crustal support (Δϵ = 1.5 ± 0.4 km) estimated from the offset in intercept values between the regional crustal thickness-elevation relationship and the global one.The red square is plotted at the subplate temperature and lithospheric thickness given by the tomography model ofCelli et al. (2020) after conversion into temperature following the methodology ofHoggard et al. (2020) andRichards et al. (2020).The yellow star is plotted at the sub-plate temperature and lithospheric thickness given by inverse modeling of rare earth element compositions(Ball et al., 2021).

Table 1
(Smithsonian, 2023)20)stimates of crustal thickness and V P /V S ratio (see TableS1in Supporting Information S1 for a complete list).Seismic stations from this study are labeled and plotted as black circles.Seismic stations from previous studies are plotted as black inverted triangles.The boundaries of the Danakil microplate and other plate boundaries are outlined by a dashed green line(Viltres et al., 2020).Dashed black lines demarcate national borders.Volcanoes from the Smithsonian Institution Global Volcanism Program are plotted as red triangles(Smithsonian, 2023).