Abstract
Earth’s only active carbonatite volcano Oldoinyo Lengai, Tanzania, is a peculiar endmember of volcanism in a young rift segment. Targeted by many petrological studies due to its effusive, cold carbonatite eruptions, its subsurface structure has only been recently explored using geophysical approaches. Here we use data from a short-term seismic deployment, to exploit coherent seismic signals in the frequency domain to detect, characterize and locate volcanic tremor. We show that narrow-band tremor is related to degassing and the ascent of carbonatite melt across the crust close to the Natron border fault, which was previously linked to the ascent of carbonatites. We further show how narrow-band and quasi-harmonic tremor partly alternate and interpret this as the interplay of deeper magmatic injections into a shallow reservoir which causes resonance. Our findings yield important implications on how physical properties, structures and dynamics of volcanic plumbing systems are linked to tremor signals.
Introduction
Seismology has been at the forefront of monitoring volcanoes and imaging magmatic systems at depth. Yet, only few systems on Earth are well monitored and even fewer well understood. In particular, we lack constraints on how the structure and dynamics of a magmatic system relate to the seismic signals we observe. Accordingly, studying magmatic systems with different properties and consequently seismic signals may illuminate aspects of magmatic plumbing systems we do not understand or even know. Here, we focus on Oldoinyo Lengai volcano which is located in the Natron Basin of the North Tanzanian Divergence, which is part of the longest continental rift worldwide: the East African Rift System (Fig. 1a). Oldoinyo Lengai is situated at the southern end of the Natron basin adjacent to thee Natron border Fault in the west. It is located next to the ~1 Ma old Gelai volcano which hosts the Naibor Soito volcanic field, a monogenetic cone field, on the southern flanks of1,2. Eruptions at Oldoinyo Lengai are known for alternating between large explosive eruptions and small-scale effusive eruptions and initiated ~0.37 Ma ago3,4,5 and. Recent explosive volcanism, generally involving nephelinitic silicate magma, occurred in 1917, 1940–1941, 1966–1967 and 2007–20083,4,6,7. Oldoinyo Lengai volcano is unique due to its production of natrocarbonatite lava with extremely low viscosity (10−1 to 102 Pa s) and low temperature (495–590 °C)8,9.
a Tectonic setting and geophysical network of the SEISVOL project. Red triangles are station locations, dashed lines indicate major faults, yellow lines depict dike intrusions in 2007. Inset shows research area within the EARS. b Schematic of Oldoinyo Lengai’s plumbing system altered from ref. 11, using results by Reiss et al.12,18 for the melt reservoir at the top of Oldoinyo Lengai and Petibon et al.72 for the melt reservoir at the foot of the volcano. NBF and GF are the Natron Border Fault and Gelai Fault.
Recent geophysical studies of the North Tanzanian Divergence and the Natron Basin have revealed a lower crustal magmatic body10 as well as a shallower melt body in the upper crust11 (Fig. 1b). This is thought to supply the 95% of eruption products that are non-carbonatite, while the production of carbonatite melts may be enhanced by a fluid-filled border fault11 which also degasses deep carbon which is sourced from the base of the Tanzania craton2. Recent geodetic observations point to a subsidence of the crater rim of 3.6 mm/year in the last ten years12 while remote sensing shows significant changes in the crater morphology including significant filling of the crater since the last explosive eruption in 2007/0813.
An enigmatic signal commonly observed near volcanoes is volcanic tremor14,15,16,17. It is generally defined as a long-lasting signal (days to months) without clear onset or end, originally thought to occur in a frequency range below 5 Hz14,15,16,17. At Oldoinyo Lengai, seismo-acoustic tremor has only been observed recently and used to infer effusive eruptions and the shallow dynamics of the crater lid18. This study, as well as others19, are increasingly challenging the assumption that tremor is a low frequency signal (<5 Hz) by reporting tremor frequencies up to 30 Hz. Many studies assume that the occurrence of tremor is a sign of an impending eruption20,21.
Generally, there is a systematic distinction between broad-band and harmonic tremor16,22,23. Broad-band tremor is characterized by distributed energy in a large frequency band and is mostly associated with eruptive/explosive activity and processes in the volcanic edifice, e.g., self-sustained fluid oscillations, magma–hydrothermal interactions, magmatic degassing, brittle fracture of melt, and solid extrusion dynamics and plug stick–slip melt14,19,24,25,26. Harmonic tremor differs from broad-band tremor in its power spectrum which usually consist of multiple peaks16, but definitions range from having only one fundamental frequency (‘monochromatic tremor’) to a fundamental with one or more overtones23. Harmonic tremor is mostly associated with fluid (magma and/or gas) flow through fractures causing resonance27,28,29, while non-linear excitation30 and pressure transients have also been suggested31,32.
Here, we use the network covariance matrix approach33, which allows tremor location without the need to pick tremor onsets34,35 as tremors have a mostly emergent onset and traditional earthquake location methods cannot be used. We characterize, detect and locate tremor at Oldoinyo Lengai volcano during a period of strong tremor activity in 2020 which falls in the period of elevated carbonatite eruptive activity starting late 201813,18. Studying tremor signals of this peculiar volcano sheds light on the dynamics of carbonatite melt and how melt transport and storage are linked to the observation of volcanic tremor.
Results
Tremor detection and location
We use seismic data from the SEISVOL network36 to detect, characterize and locate volcanic tremors during a 7-week long intense period of tremor activity from end of February to beginning of April in 2020. Within the 15-month deployment, this period only becomes outstanding considering smaller time windows, while a stronger phase of effusive eruptions and formation of a lava pond lead to strong tremor signals in 201918 (Fig. 2a). It is characterized by strong, minutes- to hours-long, recurring tremors (Methods, Fig. 2b). Tremor occurs at a weak fundamental frequency of ~1.9 Hz with overtones at ~4.8 Hz and ~7.6 Hz (these are also gliding over time). Therefore, this tremor should be called quasi-harmonic, as the overtones do not appear as integers of the fundamental. Quasi-harmonic tremors alternate with other tremor signals of different frequencies, predominantly in a band between ~2 and 4.5 Hz (Fig. 3b) which we thus refer to as ‘narrow-band’. Short-lived episodes of broad-band (2–18 Hz) tremor are also evident (see Methods for tremor nomenclature; Supplementary Figs. S1–S3). Strong frequency gliding occurs on April 2nd after nearly a week of relative tremor quiescence.
Results are displayed by frequency between February 21st and April 9th 2020 for several narrow frequency bands shown as filled circles. The size of the symbol corresponds to the number of located tremors at this location (cyan for 1–1.6 Hz, dark purple for 1.4–2.4 Hz, berry for 2.2–4.4 Hz, orange for 4.4–5.5 Hz, dark blue for 7–8 Hz). Gray dots are earthquakes from ref. 37. The dashed gray lines indicate the cross sections for the topography in the right and lower subpanels. Dashed black lines are major faults. Brown solid lines are dike intrusion from 2007/08.
We first locate tremor using a frequency band of 1–8 Hz which incorporates the most dominant tremor energy (Methods, Supplementary Figs. S4, S5). Using this frequency band, we locate ~677 h of tremor (~40% of the time period). Because these locations alternate depending on frequency content and time, which suggests that a wide tremor band consists of a mix of different tremor subtypes, we further calculate locations for narrow frequency bands to capture the differing frequency content of the narrow band and quasi-harmonic tremor. Tremor locations are clearly separated in frequency (Fig. 4) and time (S6).
Sketch of a transcrustal magmatic system with subtypes and mechanisms of volcanic tremor and cube of magmatic processes and their associated manifestations in the time, space and frequency. Observed seismic signals are: VLP Very Long Period, LP Long Period, DLP Deep Long Period, sometimes referred to as deep low frequency (DLP), VT EQs volcano-tectonic earthquakes. Stars and curly braces indicate tremor depths while the color of the star refers to the subtype of tremor signal. Mechanisms are explained in black boxes.
The best located tremors stem from the narrow-band tremor (Fig. 3a; ~547 h of tremor) and cluster dominantly underneath the northern flank of Oldoinyo Lengai at ~4–7 km depth (Fig. 4; all depths bsl). Tremor locations form tube-like structures from this location as well as from underneath the Natron border fault to shallower depth (Fig. 4). Tremor locations of the first quasi-harmonic overtone (4.4–5.5 Hz, ~330 h of tremor) are dominantly located in the shallow subsurface of the southeastern foot of Oldoinyo Lengai except during the week of 6th to 12th March, where it also appears on the northern flank (Supplementary Fig. S6). Tremor locations encompassing the fundamental frequency of the quasi-harmonic tremor (1.4–2.4 Hz) show a high variability which likely corresponds to an overall weak signal (Fig. 3a, b).
This drastically changes during the appearance of a low frequency gliding tremor on April 2nd which is superimposed on a tremor in the 1.4–2.4 Hz frequency band that produces the strongest measured spectral width (Fig. 3a) we observe. Its locations are roughly split between southwest and northeast of Oldoinyo Lengai between depths of 3 km bsl and within the edifice (Fig. 4, Supplementary Fig. S6). The locations of the gliding tremor are split in one deep cluster underneath the Natron border fault from 15 km and upwards and underneath the northeastern part of the volcano at depths of 4 km and the surface (Supplementary Fig. S3). The highest frequency range we use for the analysis (7–8 Hz) results in much fewer located tremors (~110 h) which corresponds to only short and weak signals in the spectral width (Fig. 3). Locations do not cluster very well but appear mostly scattered with few tremors appearing in the edifice.
Discussion
Tremor is ubiquitous at Oldoinyo Lengai volcano. In comparison to the seismicity patterns37,38, tremor locations cluster in the previously assumed aseismic part of the upper 6 km of the plumbing system (Fig. 4) bsl. Evidence for melt in this part of the plumbing system is manifold10,11,39,40,41. While melt bodies inferred by the different geophysical methods as well as petrological observations only partially agree, they all point towards a heavily altered subsurface, which may reduce the yield stress of rocks below the limit of earthquake generation42. This suggests that the generation of tremor in this part of the plumbing system may not be explained by an earthquake-like failure mechanism.
Attenuation imaging11 showed the main magmatic body is located away from the volcano towards Gelai, wedged between two seismicity clusters (Fig. 4; between 36° and 36.05° longitude). This was interpreted as a silicic sill structure11 which may provide the 95% of melt that is non-carbonatitic (phonolite and nephelinite3,43). Interestingly, very few tremors are located in the vicinity of this body and may yield the interpretation that the melt that was erupting in 2020 was not sourced from this reservoir. The imaging study further suggested that the Natron border fault may play an important role in the genesis and transport of carbonatites11. This seems to be corroborated by our results, as these clearly show tremor location from underneath the border fault to the edifice which invoke pipe-like structures.
The near-absence of tremor locations within the edifice and the location of the quasi-harmonic tremor at either side of the edifice seem puzzling at first. As comparison, we applied our methodology to a week of tremor data in 2019 for which we also have infrasound and thermal observations which all implied eruptive activity and shallow melt storage within the edifice18. Tremor locates in a similar location as the quasi-harmonic tremor and not on top of the edifice (Supplementary Fig. S7). Likely, the uncertainty of our velocity model in the edifice, as well as strong scattering and attenuation effects (local earthquakes were seldom observed at the crater station of the SEISVOL network) and the strong topography could factor towards a shift in locations. We further considered particle motion analysis, which showed that the crater station, which was working in May 2019 but not March 2020, has significant movement in the vertical plane (Supplementary Fig. S8) which is not evident for stations at lower elevations. This suggests that the source of the eruptive tremor signal in May 2019 may, at least partially, be closer to the foot of the volcano edifice than the top. Eruptive activity for March 2020 is unknown, as satellite images were mostly obstructed by clouds and the summit station was not working. However, particle motion analysis of signals in March 2020 supports a partially deep tremor source (Supplementary Fig. S9).
Narrow-band tremor coalesces in a cluster at around 5 km depth underneath the northern flank of Oldoinyo Lengai (Fig. 4), but also occurs at more distributed locations along pipe-like features crossing the Natron border fault. The two instances of intense broad-band tremor correspond to locations at 3–4 km depth beneath the edifice and underneath the rift flank (Supplementary Fig. S6 i-l). Other instances of broad-band tremor are too weak to be investigated in detail. Commonly, broad-band tremor is observed during eruptive periods24.
The frequency content of the narrow-band tremor is within the band for LP events (0.5–5 Hz, Fig. 3a). Instances where individual LPs and tremor have merged into and out of each other have lead to the idea that they share a common source and that tremor may simply represent a high rate of LPs19,44,45. Gliding tremors have been explained by a change of the inter-event time between individual LP events or a change in the depth46 or by changes in the parameters such as bubble concentrations in the magma which will change the depth of the tremor trigger29. Generally, fluids are thought to play an important role in the generation of LPs which are interpreted as the impulse response of a tremor-generating system47. The latter study envisioned the trigger mechanism as the rapid exsolution of gases from the fluid phase during magma ascent, where an LP event corresponded to a single triggering of the system, while continuous tremor would result from continuous triggering22. A study48 estimates that CO2 exsolution takes place at 700 °C and 100 MPa, which would correspond to a depth of ~4 km. Therefore, the cluster of tremors at a similar depth could be linked to CO2 exsolution (Fig. 4).
Regardless, ascribing the observations of LPs to any one trigger mechanism has proven to be difficult and previous studies invoke a myriad a models, e.g. self-sustained oscillations, magmatic-hydrothermal interactions, magmatic degassing, brittle failure of melt14,19,24,25,26. In Kamchatka35, tremor was ascribed to fast fluid pressure transients propagating throughout the crustal depths. Such tremor would emerge from a valving behavior. Deep tremor (40 km) has also been found in Kilauea49 and was modeled using a stationary crack model related to magma transport. At Hakone50, deep harmonic tremor (44 km) has been ascribed to magmatic fluid flow within a channel radiating seismic waves close to the upper mantle/lower crust. In both cases, the frequency content was broader than observed here. For the 2021 Taiogaite (Cumbre Vieja/La Palma) eruption, narrow-band tremor was interpreted to be caused by the fragmentation process close to the surface51. This is further supported by a modeling study of conduit flow52.
The lack of broad-band frequencies or harmonics in our narrow-band tremor observations implicate a point source process, potentially with fluids attenuating part of the frequency spectrum53. We surmise the trigger for these tremors could either be exsolution of gases or constricted flow (i.e. high rate of discrete LPs) or fluid perturbations connected to the ascent of carbonatite melt, which is surmised to be rapid and occur in pulses54 and may provide unsteady melt supply to the system. Thus, at this stage, the exact nature of narrow-band tremor generation is still unclear and requires further study (Fig. 4).
Much of the observed tremor is quasi-harmonic, a signal type which has been observed during the Cordón Caulle 2011–2012 eruption and was related to effusion rates55. Overtones at non-integer values have also been observed for (V)LP signals56,57,58. For quasi-harmonic tremor, the mechanisms for producing harmonics, such as resonance of a simple pipe-like structure, do not apply. Mechanisms of quasi-harmonics are also found in music theory of the timpani instrument as modes58, also known as Chladni pattern59,60. Depending on how the timpani is struck, different modes can be excited leading to overtones.
Different modeling studies have investigated resonance features with non-integer overtones, which are better described as ‘eigenmodes’56,60. Fluid flow in cracks with stratified fluid properties can lead to non-integer overtones of the oscillatory process60. Further, fluid-filled cracks can sustain guided waves whose properties are controlled by crack geometry and properties of the fluid56. In particular, two intersecting cracks (essentially a dike-sill structure) can produce a complex set of eigenmodes which were used to explain observed (very) long period (VLP) signals at Mayotte submarine volcano56. At Kilauea, VLP signals were inverted to constrain the conduit-reservoir system via its eigenmodes61. Water transport in glaciers through conduit-crack structures can produce complex overtones which also depend on the source function62.
The quasi-harmonic tremor we observe can thus either be generated by a surface, i.e. a penny-shaped crack, at the bottom of the edifice or be the result of the oscillation of a stratified fluid in a fluid-filled crack or the combination of a complex sill/dike geometry56,60,61,62 (Fig. 4). We observe a weak fundamental (Fig. 3a, b), which may be a result of how the modes are activated by a source. The gliding of the quasi-harmonics over time (Fig. 2a) suggests an evolution of the vibrating surface, perhaps due to rapid changes in crack stiffness.
The partially alternating behavior between quasi-harmonic and narrow-band tremor that we observe during some episodes (e.g., Fig. 3b) suggests that the system is connected: when a new pulse of carbonatite melt traverses the system at depth, it activates the vibration of the shallow sill-type structure (Fig. 4). Considering a model for the ascent of carbonatite melt, which is thought to occur in rapid pulses fracturing surrounding rock54, unsteady melt transport exciting and sustaining resonance in a complex fracture network may be the trigger for the tremors we observe. Future modeling of tremor generation in the frequency bands we observe will have to explore how pulses of carbonatite ascent can trigger resonance in dike-sill structures, fluid-filled cracks or if our tremor observations are better explained by oscillations of stratified fluids.
Our observations (and ref. 18) clearly show that tremor at Oldoinyo Lengai is not necessarily a low frequency (<5 Hz), shallow signal, and that one volcano can produce a wide variety of different tremor signals which reflect different processes at various depths. We can only speculate to what extent the tremor signals are controlled by the ultra-low viscosity of the melt, or the complexity of the plumbing system. A recent study of tremor of the 2021 La Palma eruption shows that the tremor amplitude is controlled by the viscosity51. Considering fluid dynamics in conduits63, they conclude that lower viscosity fluids produce higher frequency, long-lived tremor while higher viscosity magma damps high frequencies. Similar findings are reported by Balmforth et al.64 who find that oscillations of the fluid are faster when the viscosity is lower.
Compared to studies at different volcanoes, we do not locate tremor deep into the crust (not below 20 km; Fig. 4) and have not detected deep-long period events (DLPs) yet. Deep harmonic tremor underneath Hawaii49, and Hakone50 are linked to fluid flow while deep narrow-band tremor observations in Kamchatka are attributed to pressure transients in a fluid-saturated permeable medium35. In both systems, deep tremor appears in spatial and/or temporal connection to DLPs, which may be connected to deep degassing in some cases.
While narrow-band tremor appears through the crust in Kamchatka35 and at Oldoinyo Lengai (Fig. 4), where we invoke a point source or pressure excitation mechanism, eruptions may often produce broad-band tremor24 and at Oldoinyo Lengai, extreme high-frequency tremor18. However, detailed studies of different tremor types and their locations are missing. Therefore, we can only speculate how they exactly relate to proposed magmatic processes at depth (Fig. 4). Our study shows the need for careful re-analysis and location of tremor subtypes at active volcanoes which may elucidate their specific plumbing systems. Comparisons across systems may lead to crucial insights into commonalities and differences of tremor observations and their connections to pre- and syn-eruptive processes. This may help to better understand this particular seismic signal in the future.
Online methods
Estimating the covariance matrix
Based on the work of ref. 33, this study65 was the first to use the network covariance matrix approach to study volcanic tremor across a seismic network. Our approach is based on their work as well as refs. 34,35. We use data from a subset of the SEISVOL stations shown in Fig. 1. First, we read hour-long files of seismic data. We detrend the data and remove the instrument response. We filter the data between 0.01 and 20 Hz and downsample it by a factor of 2. We then divide the data into averaging windows of 10 min. Each averaging window of length \(\Delta t\) (10 min) is further subdivided into M subwindows which overlapped r = 50%, so that \(\Delta t={Mr}\delta t\). We applied spectral whitening on these subwindows. We did not apply amplitude clipping as this may lead to loss of important signal information outlined in ref. 65. We then calculate the cross spectra for each subwindow and station combination. The covariance matrix C(f) is then defined as the cross spectra for each subwindow and station combination33.
We used a sub window length of 48 s and M = 20 and stacked all covariance matrices in a given averaging window. This allows for a reasonably small resolution in time (compare refs. 34,35).
The covariance matrix is then decomposed into eigenvectors, which we can use to define the spectral width33,34,35,65 which provides a measure of the number of independent sources. It displays high values if the wavefield is composed of many different sources and its eigenvalue distribution is flat, which, in sum, are equal to noise, and is minimal when only one source coherently traverses the seismic network34.
Tremor nomenclature
One of the main problems of tremor descriptions is the lack of fixed terminology. We define three different tremor categories based on their appearance in the spectral width (Fig. 2), i.e. the coherent part of the tremor energy across the array: quasi-harmonic, narrow-band and broadband. Quasi-harmonic tremor occurs at a weak fundamental frequency of ~1.9 Hz with overtones at ~4.8 Hz and ~7.6 Hz (Supplementary Fig. S1). We call this tremor quasi-harmonic, as the overtones do not appear as integers of the fundamental and the term ‘harmonics’ are often used by seismologists to describe peaked spectra with overtones at integer values. Quasi-harmonic tremor is interrupted by a tremor with frequencies between 2 and 7 Hz, but most dominantly its energy is ~2–4.5 Hz which we thus refer to as ‘narrow-band’. This tremor also seems to have quasi-harmonics in some time windows, but for better distinction from the peaked quasi-harmonic spectra we call it ‘narrow-band’ (Supplementary Fig. S1). We define broad-band tremor as a signal between 2 and 18 Hz.
These characteristics are of course also present in single station spectrograms (Supplementary Fig. S2) and in their power spectra for single time-windows (Supplementary Fig. S3) albeit slightly different from station to station. This is likely the result of the radiation pattern and site effects. Quasi-harmonics are clearly visible, and the ‘narrow-band’ tremor seems to be slightly more broadband than in the spectral width plots and sometimes also has overtones. Over time, there is more variability than we can account for here, the frequency bands chosen for the location procedure are those that have the most coherent energy over the time period.
Automatic localization of tremors
After visual inspection of the spectral width, we choose several frequency bands which showed tremor sources (see Fig. 2 for 10-min intervals for March 2020). We use these frequency bands to create the wavefield as a cross-correlation function in the time domain by an inverse Fourier transform of the ‘filtered’ covariance matrix computed from the complex outer product of its first eigenvector and its Hermitian transpose34. This approach is used to denoise the data and enhance the dominant sources. We then calculate the envelope of the cross-correlation signal. In contrast to refs. 34,35, we apply minimal smoothing e.g. (three orders of magnitude less than ref. 35). Smoothing was originally applied to dampen the effects of scattering and imprecision of the velocity model34, but we find it reduces the sharpness of the localization.
We create a 3D-grid and use a 1-D velocity model based on ref. 66 to calculate theoretical rays and travel times from every grid point to every station location using a simple ray tracing approach. We assume the tremor field consists of S-waves34. For every grid point \(r\) of the 3D-grid and every cross-correlation envelope we look-up theoretical travel times between this point and the stations, calculate the differential travel time, and sum the value of the cross-correlation envelope at this lag time which results in the spatial likelihood function. We further normalize by the total grid volume and the maximum likelihood value over a given time period which gives us the network response function.
We automatically compute a solution for every time step in our spectral width. To choose which correspond to a tremor signal, we impose three criteria: the mean value of the spectral width has to be below 0.80, the number of elements of the grid which contain 95% of the maximum grid value has to be less than 1% of the grid and the location must not be equal to the boundaries of the grid. See Supplementary Fig. S4 for a single location example.
Location errors
Thus far, location errors have been estimated from the distribution of the network response34,35 considering the volume of values that constitute 95% of the maximum network response (not to be confused with a 95% confidence level). We calculate statistical errors for the derived tremor localization based on the traditional error calculations of an earthquake location67,68. This is usually defined as the misfit or residual between calculated and observed arrival times \(r\) and given by the variance-covariance matrix \({\sigma }_{x}^{2}\) of the hypocentral parameters x,y,z and origin time. It can be calculated by
where \({\sigma }^{2}\) is the variance \({\sigma }^{2}\) of the arrival times multiplied with the identity matrix. \(G\) is the matrix of partial derivatives of all hypocentral parameters which relates changes in location to changes in location via
where \({r}_{i}/r\) is the residual and \({t}_{i}^{{tra}}\) the travel times.
\({G}^{T}\) is \(G\) transposed. Every partial derivative can be expressed by
Where v is the velocity, x, y, and z hypocentral parameters and \({x}_{i}\), \({y}_{i},{z}_{i}\) station parameters.
The standard deviation of the hypocentral parameters is given by the square root of the diagonal elements of \({\sigma }_{x}^{2}\).
Given that we locate using cross-correlation functions based on two stations, our spatial derivatives, \(G\) do not only depend on one partial derivative, but two for each station contributing to the cross-correlation function
and \({\sigma }^{2}\) is the variance between the maximum of the cross-correlation function and the differential travel time.
Locations often display small statistical errors in latitude and longitude of 2–3 km, while the depth error is often very large, with smallest values of 10 km to infinity (Supplementary Table S1). Errors derived from the area containing the 95% of the maximum network response value yield smaller errors in depth, but considering the network response overall, which is often smeared across a wide depth range, assuming larger uncertainties in the vertical location is appropriate. Tests have demonstrated that the error depends on the station distribution. Adding more distant stations may decrease the error but it also decreases the signal in the spectral width plot so the overall location is worse. Shallow locations display extremely large vertical errors, often to infinity while deeper locations have depth errors around 10 km. Errors are also dependent on the used frequency band, with very narrow frequency bands generally producing larger errors. This is likely caused by the fact that the cross-correlation functions are very mono-frequent and other maxima, which are used in the error calculations, can easily be dominant. Another problem arises from the quasi-harmonic nature of the tremor. As described in the discussion, quasi harmonics are caused by surfaces, along which any location is feasible. Hence, the error may reflect that.
Tremor density
To better show which areas produce tremor more often, we simply sum the number of tremors per grid point. Tremor locations are highly frequency and time dependent. For the broad frequency band of 1–8 Hz, we locate ~677 hours of tremor (Supplementary Fig. S5). Based on the maxima of the spectral width plots, we use frequency bands of 1–1.6 Hz, 1.4–2.4 Hz, 2.2–4.4 Hz. 4.4–5.5 Hz and 7–8 Hz, and we locate ~94, ~214, ~547, ~330 and ~110 h of tremor. Tremors of the same frequency often appear in distinct locations; particularly north and south of the edifice in different depths down to ~15 km (Fig. 3) but there is also variability. To look at this in more detail, Supplementary Fig. S3 shows a weekly break-up of the data set.
The first week of tremor we consider (21st–27th of February, Supplementary Fig. S6a, c) shows locations for three dominant tremor bands. The broad frequency (1–8 Hz) locates mostly south of the edifice and invokes the shape of tubes from ~18 km depth to the surface. Considering individual locations plots, this may be due to some smearing in the network response function. When broken up into smaller frequency bands, only the narrow-band tremor (2.2–4.4 Hz) and the 4.4–5.5 Hz tremor bands produce tremor locations with non-infinite depth errors. While the narrow-band tremor locates north of the edifice and towards the border fault at depths down to 5 km, the 4.4–5.5 Hz tremor locates very shallowly at the southern flank of Oldoinyo Lengai. This split in location between the two frequencies remains quite stable over the next weeks (Supplementary Fig. S6b-l), with the exception of the week of 6th to 12th March (Supplementary Fig. S6e, g), where it also appears on the northern flank. Depending on which of the two signals is more dominant, the 1–8 Hz frequency analysis band seems to shift to either locations or appear somewhat in between.
The 7–8 Hz frequency band only shows few locations and tremor locations encompassing the fundamental frequency of the quasi-harmonic tremor (1.4–2.4 Hz) show a high variability which likely corresponds to an overall weak signal (Fig. 2, Supplementary Fig. S6) within this frequency band. This drastically changes during the appearance of a superimposed gliding tremor in a lower frequency band on April 2nd. Then, the 1.4–2.4 Hz frequency band produces the strongest measured spectral width meaning it is the most coherent signal we observe. Its locations are roughly split between southwest and northeast of Oldoinyo Lengai between depths of 3 km of the crust and within the edifice (Supplementary Fig. S6j, l). The gliding tremor is well located using the lowest frequency band we employ for the analysis. It is also split in locations, one deep cluster underneath the Natron border fault from 15 km and upwards and underneath the northeastern part of the volcano at depths of 4 km and the surface.
Comparisons to an eruptive period: locations May 2019
For comparison, we locate tremors for an eruptive period from 26th May 2019 to 1st of June 2019 (see ref. 18). Unfortunately, only 4 stations around Oldoinyo Lengai were operating at this time, including the crater station. We use two frequency bands for the analysis, the broader band 1–8 Hz for straightforward comparison to our other data and a more narrow-band 3.5–6 Hz frequency band which incorporates the most dominant energy during this tremor episode (for spectral width, see Supplementary Fig. S7a). For the broader frequency band, tremor locates partially in the edifice but much more dominantly surrounding the foot of the edifice as well as northeast and south west of the edifice forming pipe-like structures down to depths of ~7 km (Supplementary Fig. S7b). For the narrower frequency band, tremor mostly locates in the very shallow surface south of the edifice.
Particle motion analysis
Particle motion analysis is often used to infer the polarization and angle of incidence of an incoming wavefield69. While it is tricky and erroneous to infer source locations, it can be used to check whether the wavefield originates close to the surface or from depth and can be used as an independent measurement to check whether our locations are in line with the observed angle of incidence. We show a particle motion analysis of a strong tremor episode associated with eruptive activity in May 2019 in S8 and compare this to and a non-eruptive tremor episode in March 2020 in S9, which clearly shows the difference in signals due to different source processes. We use a bandpass filter between 3 and 6 Hz, after removing the mean and applying a linear detrending. The polarization features are computed using obspy70 which implements the formulation of Vidale71 with a 20-s window length and a sliding window fraction of 0.2.
Data availability
Tremor locations results have been uploaded to figshare https://doi.org/10.6084/m9.figshare.30053218. The dataset is archived at GEOFON Data Center, https://doi.org/10.14470/4W7564850022.
Code availability
Codes used in this work can be made available upon request.
References
Mana, S., Furman, T., Turrin, B. D., Feigenson, M. D. & Swisher, C. C. Magmatic activity across the East African North Tanzanian Divergence Zone. J. Geol. Soc. 172, 368 (2015).
Muirhead, J. D. et al. Displaced cratonic mantle concen-trates deep carbon during continental rifting. Nature 582, 67–72 (2020).
Dawson, J. Sodium carbonate lavas from Oldoinyo Lengai, Tanganyika. Nature 195, 1075 (1962).
Dawson, J. B., Keller, J., Nyamweru, C. Historic and recent eruptive activity of Oldoinyo Lengai. In Carbonatite volcanism: Oldoinyo Lengai and the petrogenesis of natrocarbonatites, IAVCEI Proc Volcanol, 4 (eds, Bell K., Keller J.) 4–22 (Springer, 1995).
Sherrod, D. R., Magigita, M. M. & Kwelwa, S. Geologic map of Oldoinyo Lengai (Oldoinyo Lengai) and surroundings, Arusha Region, United Republic of Tanzania. U.S. Geol. Surv. Rep. 65, 2013-1306 (2013).
Keller, J., Klaudius, J., Kervyn, M., Ernst, G. G. J. & Mattsson, H. B. Fundamental changes in the activity of the natrocarbonatite volcano Oldoinyo Lengai, Tanzania. Bull. Volcanol. 72, 893–912 (2010).
Klaudius, J. & Keller, J. Peralkaline silicate lavas at Oldoinyo Lengai, Tanzania. Lithos 91, 173–190 (2006).
Krafft, M. & Keller, J. Temperature measurements in carbonatite lava lakes andflows, Oldoinyo Lengai, Tanzania. Science 245, 168–170 (1989).
Fischer, T. P. et al. Upper-mantle volatile chemistry at Oldoinyo Lengai volcano and the origin of carbonatites. Nature 459, 77–80 (2009).
Roecker, S. et al. Subsurface images of the Eastern Rift, Africa, from the joint inversion of body waves, surface waves and gravity: Investigating the role of fluids in early-stage continental rifting. Geophys. J. Int. 210, 931–950 (2017).
Reiss, M. C., De Siena, L. & Muirhead, J. D. The Interconnected magmatic plumbing system of the Natron rift. Geophys. Res. Lett. 49, e2022GL098922 (2022).
Wauthier, C. & Ho, C. Satellite geodesy unveils a decade of summit subsidence at Ol Doinyo Lengai Volcano, Tanzania. Geophys. Res. Lett. 51, e2023GL107673 (2024).
Tournigand, P.-Y. et al. Remote volcano monitoring using crowd-sourced imagery and Structure-from-Motion photogrammetry: a case study of Oldoinyo Lengai’s active pit crater since the 2007–08 paroxysm. J. Volcanol. Geothermal Res. 443, https://doi.org/10.1016/j.jvolgeores.2023.107918 (2023).
McNutt S. R. Volcanic tremor. In Encyclopedia of Earth System Science (ed. Nierenberg W. A.) 4 (Academic Press, 1992).
McNutt, S. R. Seismic monitoring and eruption forecasting of volcanoes: a review of the state-of-the-art and case histories. In Monitoring and Mitigation of Volcanic Hazards (ed. Scarpa, T.) 100–146 (Springer, 1996).
Konstantinou, K. I. & Schlindwein, V. Nature, wavefield properties and source mechanism of volcanic tremor: a review. J. Volcanol. Geotherm. Res. 119, 161–187 (2002).
McNutt, S. R. & Roman, D. C. Volcanic seismicity. In: The Encyclopedia of Volcanoes, 1011–1034 (Academic Press, 2015).
Reiss, M. C. et al. Overview of seismo-acoustic tremor at Oldoinyo Lengai, Tanzania: Shallow storage and eruptions of carbonatite melt. J. Volcanol. Geothermal Res. https://doi.org/10.1016/j.jvolgeores.2023.107898 (2023).
Hotovec, A. J., Prejean, S. G., Vidale, J. E. & Gomberg, J. Strongly gliding harmonic tremor during the 2009 eruption of Redoubt volcano. J. Volcanol. Geotherm. Res. 259, 89–99 (2013).
D’Agostino, M. et al. Volcano monitoring and early warning on Mt. Etna, Sicily, based on volcanic tremor: Methods and technical aspects. In Complex Monitoring of Volcanic Activity (ed, Zobin, V. M.) 53–92 (Nova Science Publishers, Inc, 2013).
Chouet, B. A. Long-period volcano seismicity: its source and use in eruption forecasting. Nature 380, 309–316 (1996a).
Matoza, R. & Roman, D. One hundred years of advances in volcano seismology and acoustics. Bull. Volcanol. 84, 86 (2022).
Roman, D. C. Automated detection and characterization of harmonic tremor in continuous seismic data. Geophys. Res. Lett. 44, https://doi.org/10.1002/2017GL073715 (2017).
Chouet, B. A. & Matoza, R. S. A multi-decadal view of seismic methods for detecting precursors of magma movement and eruption. J. Volcano. Geotherm. Res 252, 108–175 (2013).
Chouet B. A. (1996b) New methods and future trends in seismological volcano monitoring. In Monitoring and Mitigation of Volcano Hazards (eds, Scarpa, R. & Tilling R. I.) 23–97 (Springer-Verlag, 1996).
Aki, K., Fehler, M. & Das, S. Source mechanism of volcanic tremor: fluid driven crack models and their application to the 1963 Kilauea eruption. J. Volcano. Geotherm. Res. 2, 259–287 (1977).
Chouet, B. A. Resonance of a fluid-driven crack: radiation properties and implications for the source of long-period events and harmonic tremor. J. Geophys. Res. 93, 4375–4400 (1988).
Lesage, P., Mora, M. M., Alvarado, G. E., Pacheco, J. & Metaxian, J. P. Complex behavior and source model of the tremor at Arenal volcano, Costa Rica. J. Volcano. Geotherm. Res. 157, 49–59 (2006).
Benoit, J. P. & McNutt, S. R. New constraints on source processes of volcanic tremor at Arenal Volcano, Costa Rica, using broadband seismic data. Geophys. Res. Lett. 24, 449–452 (1997).
Julian, B. R. Volcanic tremor: nonlinear excitation by fluid flow. J. Geophys. Res. 99, 11859–11877 (1994).
Schlindwein, V., Wassermann, J. & Scherbaum, F. Spectral analysis of harmonic tremor signals at Mt. Semeru volcano, Indonesia. Geophys. Res. Lett. 22, 1685–1688 (1995).
Hagerty, M. T., Schwartz, S. Y., Garces, M. A. & Protti, M. Analysis of seismic and acoustic observations at Arenal Volcano, Costa Rica, 1995–1997. J. Volcano. Geotherm. Res. 101, 27–65 (2000).
Seydoux, L., Shapiro, N. M., De Rosny, J., Brenguier, F. & Landes, M. Detecting seismic activity with a covariance matrix analysis of data recorded on seismic arrays. Geophys. J. Int. 204, 1430–1442 (2016).
Soubestre, J. et al. Depth migration of seismovolcanic tremor sources below the Klyuchevskoy Volcanic group (Kamchatka) determined from a network-based analysis. Geophys. Res. Lett. 46, 8018–8030 (2019).
Journeau, C. et al. Seismic tremor reveals active trans-crustal magmatic system beneath Kamchatka volcanoes. Sci. Adv. 8, https://doi.org/10.1126/sciadv.abj1571 (2022).
Reiss, M. C., Rümpker, G. SEISVOL - Seismic and Infrasound Networks to study the volcano Oldoinyo Lengai. GFZ Data Services. [Dataset]. https://doi.org/10.14470/4W7564850022 (2020).
Reiss, M. C. et al. The impact of complex volcanic plumbing on the nature of Seismicity in the developing Magmatic Natron Rift, Tanzania. Front. Earth Sci. 8, 609805 (2021).
Weinstein, A., Oliva, S. J., Ebinger, C. J., Roecker, S. & Tiberi, C. Fault magma interactions during early continental rifting: seismicity of the Magadi-Natron-Manyara basins, Africa. Geochem. Geophys. Geosyst. 18, 3662–3686 (2017).
Plasman, M. et al. Lithospheric low-velocity zones associated with a magmatic segment of the Tanzanian Rift, East Africa. Geophys. J. Int. 210, 465–481 (2017).
Biggs, J., Chivers, M. & Hutchinson, M. C. Surface deformation and stress interactions during the 2007–2010 sequence of earthquake, dyke intrusion and eruption in northern Tanzania. Geophys. J. Int. 195, 16–26 (2013).
Daud, N. et al. Elucidating the magma plumbing system of Ol Doinyo Lengai (Natron Rift, Tanzania) Using satellite geodesy and numerical modeling. J. Volcanol. Geothermal Res. 438, https://doi.org/10.1016/j.jvolgeores.2023.107821 (2023).
Liu, K. & Zhao, J. Progressive damage behaviours of triaxially confined rocks under multiple dynamic loads. Rock. Mech. Rock. Eng. 54, 3327–3358 (2021).
Berkesi, M., Bali, E., Bodnar, R. J., Szabó, Á., & Guzmics, T. Carbonatite and highly peralkaline nephelinite melts from Oldoinyo Lengai Volcano, Tanzania: the role of natrite-normative fluid degassing. Gondwana Res. 85, 76–83 (2020).
Fehler, M. Observations of volcanic tremor at Mount St. Helens Volcano. J. Geophys. Res. 88, 3476–3484 (1983).
Neuberg, J., Baptie, B., Luckett, R. & Stewart, R. Results from the broadband seismic network on Montserrat. Geophys. Res. Lett. 25, 3661–3664 (1998).
Neuberg, J., Luckett, R., Baptie, B. & Olsen, K. Models of tremor and low-frequency earthquake swarms on Montserrat. J. Volcano. Geotherm. Res. 101, 83–104 (2000).
Chouet, B. Excitation of a buried magmatic pipe: a seismic source model for volcanic tremor. J. Geophys. Res. 90, 1881–1893 (1985).
de Moor, M. et al. Volatile-rich silicate melts from Oldoinyo Lengai volcano (Tanzania): implications for carbonatite genesis and eruptive behavior. Earth Planet. Sci. Lett. 361, 379–390 (2012).
Aki, K. & Koyanagi, R. Deep volcanic tremor and magma ascent mechanism under Kilauea, Hawaii. J. Geophys. Res. 86, 7095–7109 (1981).
Yukutake, Y. et al. Harmonic tremor from the deep part of Hakone volcano. Earth Planets Space 74, 144 (2022).
Longpré, M. A. et al. Shifting melt composition linked to volcanic tremor at Cumbre Vieja volcano. Nat. Geosci. 1–9, https://doi.org/10.1038/s41561-024-01623-x (2025).
Coppess, K., Lam, F., & Dunham, E. Volcanic eruption tremor from particle impacts and turbulence using conduit flow models. Seismica, 4, https://doi.org/10.26443/seismica.v4i1.1285 (2025).
Kumagai, H. & Chouet, B. A. The dependence of acoustic properties of a crack on the resonance mode and geometry. Geophys. Res. Lett. 28, 3325–3328 (2001).
Walter, B. F. et al. Fluids associated with carbonatitic magmatism: A critical review and implications for carbonatite magma ascent. Earth Sci. Rev. 215, 103509 (2021).
Bertin, D. et al. High effusion rates of the Cordón Caulle 2011–2012 eruption (Southern Andes) and their relation with the quasi-harmonic tremor. Geophys. Res. Lett. 42, 7054–7063 (2015).
Liang, C., Peng, J., Ampuero, J. P., Shauer, N. & Dai, K. Resonances of fluid-filled cracks with complex geometry and application to very long period (VLP) seismic signals at Mayotte submarine volcano. J. Geophys. Res.Solid Earth 129, e2023JB027844 (2024).
Kumagai, H. Temporal evolution of a magmatic dike system inferred from the complex frequencies of very long period seismic signals. J. Geophys. Res. Solid Earth 111, https://doi.org/10.1029/2005JB003881 (2006).
Jones, R. K. The well-tempered timpani, in search of the missing fundamental. https://wtt.pauken.org/preface. Last accessed on 25th September 2025.
Ramsey, G. P. Percussion Instrument Group. In: The Physics of Music. Undergraduate Lecture Notes in Physics. Springer, https://doi.org/10.1007/978-3-031-53507-9_8 (2024).
Karlstrom, L. & Dunham, E. M. Excitation and resonance of acoustic-gravity waves in a column of stratified, bubbly magma. J. Fluid Mech. 797, 431–470 (2016).
Crozier, J. & Karlstrom, L. Evolving magma temperature and volatile contents over the 2008–2018 summit eruption of Kīlauea Volcano. Sci. Adv. 8, eabm4310 (2022).
McQuillan, M. & Karlstrom, L. Fluid resonance in elastic-walled englacial transport networks. J. Glaciol. 67, 999–1012 (2021).
Ferrick, M. G., Qamar, A. & St. Lawrence, W. F. Source mechanism of volcanic tremor. J. Geophys. Res. Solid Earth 87, 8675–8683 (1982).
Balmforth, N. J., Craster, R. V. & Rust, A. C. Instability in flow through elastic conduits and volcanic tremor. J. Fluid Mech. 527, 353–377 (2005).
Soubestre, J. et al. Network-based detection and classification of seismovolcanic tremors: Example from the Klyuchevskoy volcanic group in Kamchatka. J. Geophys. Res.Solid Earth 123, 564–582 (2018).
Albaric, J. et al. Contrasted seismogenic and rheological behaviours from shallow and deep earthquake sequences in the North Tanzanian Divergence, East Africa. J. Afr. Earth Sci. 58, 799–811 (2010).
Stein, S & Wysession, M. An introduction to seismology, earthquakes, and earth structure. (John Wiley & Sons, 2009).
Shearer, P. M. Introduction to seismology. (Cambridge university press, 2019).
Caudron, C. et al. Anatomy of phreatic eruptions. Earth Planets Space 70, 1–14 (2018).
Krischer, L. et al. ObsPy: A bridge for seismology into the scientific Python ecosystem. Comput. Sci. Discov. 8, 014003 (2015).
Vidale, J. E. Complex polarization analysis of particle motion. Bull. Seismol. Soc. Am. 76, 1393–1405 (1986).
Petibon, C. M., Kjarsgaard, B. A., Jenner, G. A. & Jackson, S. E. Phase relationships of a silicate-bearing natrocarbonatite from Oldoinyo Lengai at 20 and 100 MPa. J. Petrol. 39, 2137–2151 (1998).
Acknowledgements
We would like to extend our gratitude to everyone involved with the fieldwork, especially members of Engaresero Village, and Amani Laizer and Emmanuel Kazimoto. We thank the Tanzania Commission for Science and Technology (COSTECH) for supporting us and providing the research permits. The equipment for the seismic station was kindly provided by the Geophysical Instrument Pool Potsdam (GIPP) and the Goethe-University Frankfurt. Data were archived by GEOFON (GEOFON Data Centre, 1993). MCR acknowledges funding for the SEISVOL project by the DFG (German Research Council), grant number RE 4321/1-1.
Funding
Open Access funding enabled and organized by Projekt DEAL.
Author information
Authors and Affiliations
Contributions
M.C.R. is the PI of the SEISVOL project, designed this study, coded the analyses tools for tremor location, processed the seismic data, composed all figures. C.C. contributed the particle motion analyses. P.H. assisted in coding the tremor location tools. D.R. and C.C. contributed to the data interpretation. M.C.R. wrote the manuscript with contributions from C.C. and D.M.
Corresponding author
Ethics declarations
Competing interests
The authors declare no competing interests.
Peer review
Peer review information
Communications Earth & Environment thanks Leif Karlstrom and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Primary Handling Editors: Lucia Pappalardo, Joseph Aslin and Carolina Ortiz Guerrero. A peer review file is available.
Additional information
Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Supplementary information
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
About this article
Cite this article
Reiss, M.C., Caudron, C., Hering, P. et al. Tremor signals reveal the structure and dynamics of the Oldoinyo Lengai magmatic plumbing system. Commun Earth Environ 6, 781 (2025). https://doi.org/10.1038/s43247-025-02804-1
Received:
Accepted:
Published:
DOI: https://doi.org/10.1038/s43247-025-02804-1