Further constraints on the shear wave velocity structure of Cameroon from joint inversion of receiver function, Rayleigh wave dispersion and ellipticity measurements

SUMMARY

To develop a new model of the Cameroon crust that would provide further constraints on its composition and help study regional tectonics, we simultaneously invert Rayleigh wave group and phase velocity dispersion, P-wave receiver function and Rayleigh wave ellipticity measurements using a joint inversion method based on the neighbourhood algorithm (NA). We obtain ensembles of 1-D isotropic Vs models along with layered Vp/Vs ratio and crustal thickness estimate beneath 32 broad-band seismic stations in the region. Due to the addition of ellipticity data, the resolutions of shallow structures are well improved in our new model, thereby showing a better spatial correlation with known tectonic features. By interpolating the ellipticity measurements and performing joint inversions with group and phase dispersion data at several grid nodes, we are able to image the geometry and constrain the locations of mafic bodies (Vs of 3.7–4 km s^–1) that intruded into the upper crust (2–8 km) during the formation of Gondwana. Investigating the relationship between crustal parameters, we found a positive correlation between crustal thickness, Vp/Vs ratio, the Bouguer anomaly and topography. Specifically, we found that the inverted Moho depth and surface topography within the Adamawa Plateau and the Garoua Rift agrees with Airy isostasy. The inverted Vp/Vs ratio shows variation that ranged from 1.67 to 1.85 with uncertainties that are generally less than 0.05 at crustal depths. A regional averaged Vp/Vs ratio of 1.73 and a standard deviation of ∼0.03 in the crust suggest a felsic to intermediate bulk crustal composition, indicating that the Cenozoic volcanism is most likely characterised by small volumes of mafic intrusion with limited alteration to the crust, thereby supporting plume-free hypothesis such as small-scale mantle convection for the origin of the Cameroon Volcanic Line.

Composition and structure of the continental crust, Africa, Joint inversion, Crustal structure

1 INTRODUCTION

The Cameroon region is located in the interior of the relatively stable African continent and contains a variety of tectonic features that span more than 3 billion years of the geological record. Examples of major tectonic features that characterize the region include the old Achaean Congo Craton, the Proterozoic Oubanguides Belts which contain several major shear zones such as the Central African Shear Zone (CASZ), Garoua Rift and Mamfe Basin, an eastward extension of the Benue Trough from Nigeria to Cameroon formed during the breakup of Gondwana, and most importantly the Cameroon Volcanic Line (CVL), a chain of active and dormant alkaline volcanoes which began ca. 30 Ma without age progression (Fitton & Dunlop 1985; Fitton 1987; Noel et al.2014). The detailed structure, composition and geodynamic evolution of the CVL and other tectonic features in Cameroon are yet to be fully understood and are still the subject of considerable debate in the literature (e.g. Reusch et al.2010; De Plaen et al.2014; Adams et al.2015; Guidarelli & Aoudia 2016). In order to provide new insights on the structure and composition of the Cameroon crust and contribute to the on-going debate on the origin of the CVL, it is crucial to take advantage of the new geophysical data as they become available for improving crustal model of the Cameroon (e.g. Lin et al.2012, 2014). This is especially important for understanding the CVL since the competing geodynamic models for its origin make different predictions for the crustal structure/composition. For example, if large quantity of melt exists today within the Cameroon crust, it will be expected to manifest as anomalously slow shear wave velocity and elevated bulk crustal Vp/Vs ratios, whereas cooled mafic addition to the crust will be expected to increase wave speeds (e.g. Darbyshire et al.1998; Stuart et al.2006; De Plaen et al.2014). Therefore, in this study, we attempt to improve on the work of Tokam et al. (2010) which applied the linearized least squares inversion method of Juliá et al. (2000, 2003) to obtain 1-D Vs models by performing joint inversions of teleseismic group velocity dispersion and receiver function beneath the 32 stations of the Cameroon Broadband Seismic Experiment (CBSE). Although the study obtained excellent results that helped advance the understanding of crustal structure in Cameroon, the linearized inversion method used in the previous study is less suitable for quantifying the uncertainties of the inverted models and requires an initial model that is close to the true solution to reduce the probability of converging to a local minimum. Additionally, the group velocity dispersion obtained from teleseismic earthquake data using the single station method of Dziewonski et al. (1969) as employed in the work of Tokam et al. (2010) limits the ability of obtaining accurate dispersion measurements at short periods (<20 s) needed for high-resolution imaging of the upper crust. Therefore, we adopt a two-fold improvement that involves the use of advanced inversion methodology and the inclusion of new seismic datasets.

In the past few decades, substantial efforts have been geared towards the joint inversion of multiple geophysical datasets to ameliorate the limited constraints by individual observations on Earth's structure. Many of such studies have involved the use of surface wave dispersion complemented by either receiver functions or Rayleigh wave ellipticity and applied in many regions of the world (e.g. Juliá et al.2000, 2003; Bodin et al.2012; Shen et al.2013; Akpan et al.2016). Through these joint inversion schemes, it has been possible to put a tighter constraint on the average S-wave velocity as well as the velocity contrast across intracrustal discontinuities by taking advantage of the strengths of different data sets while also reducing model non-uniqueness and the trade-offs in the inversion (e.g. Langston 1979; Ammon et al.1990; Lin et al.2012; Chong et al.2016; Li et al.2016). In the same manner, a new series of study that involves the simultaneous inversion of Rayleigh wave dispersion, receiver function and ellipticity data are emerging and have demonstrated a significant improvement in the resolution of the upper crustal structure. For example, Kang et al. (2016) used a joint Bayesian Monte Carlo inversion scheme to study the crustal structure in northeastern China and a similar study have also been published for the entire United States by Shen & Ritzwoller (2016). More recently, Zhang & Yao (2017) proposed a stepwise linearized joint inversion method for these multiple seismic data sets and demonstrated the robustness of their algorithm using both synthetic and field data in Northeast China (Kang et al.2016; Li et al.2016). The list of publications detailing this kind of study is rapidly increasing, as the joint inversion of these three data sets is being applied to many other regions. In this study, we develop a joint inversion scheme based on a different global inversion method referred to as the neigbourhood algorithm (NA, Sambridge 1999a, b) to achieve a more robust model of the crustal structure as well as error estimates in the Cameroon region.

The recent advancement in ambient noise tomography (ANT) method has paved a way for the extraction of reliable short-period group and phase velocity dispersion data from ambient seismic noise and has been routinely used to image the subsurface (e.g. Shapiro et al.2005; Ojo et al.2017). Similarly, various studies have demonstrated that the Rayleigh wave ellipticity data have a stronger sensitivity to the shear wave velocity structure at shallow depths as well as to small-scale structures (Boore & Toksöz 1969; Sexton et al.1977; Tanimoto & Alvizuri 2006; Poggi & Fah 2010; Maupin 2017). Lin et al. (2012) also showed that the Rayleigh wave ellipticity data can be used to constrain the Vp/Vs ratio and the density in the upper crust. In this study, we refer to the Rayleigh wave ellipticity data as 'ZH ratio’ and defined it as the amplitude ratio of the vertical and radial components of the fundamental mode Rayleigh wave at different periods following the nomenclature of Tanimoto & Rivera (2008). Being a single station measurement from three-component records, the ZH ratio is primarily sensitive to structures near the seismic station, thereby reducing the along-path contaminations and earthquake source complexities (e.g. Tanimoto & Rivera 2008; Yano et al.2009; Lin et al.2012, 2014; Chong et al.2015). In recent years, there are increasing numbers of studies employing the ZH ratio to obtain a more accurate image of the uppermost crustal structure in different regions of the world (e.g. Berbellini et al.2016, 2017; Kang et al.2016; Li et al.2016; Shen & Ritzwoller 2016). Therefore, the second improvement we offer in this study is the use of multiple seismic data sets which involve the addition of Rayleigh wave ellipticity data and short-period dispersions obtained from Ambient Seismic Noise data processing in order to put tighter constraints on the uppermost crustal structure in Cameroon (e.g. Shen & Ritzwoller 2016; Zhang & Yao 2017). This is helpful for understanding the surface geology and its relationship with deep crustal processes, ground motion simulations for geohazard assessments as well as constraining deeper structures (Lin et al.2012; Li et al.2016). Specifically, for the region of Cameroon, the inclusion of the ZH ratio data has the potential of improving the resolution of the uppermost crustal structure, enabling us to contribute significantly to the research on magmatism-induced crustal modification along the CVL and to investigate the existence of a high-velocity layer (HVL) in the upper crust in Cameroon (Tokam et al.2010).

One of the interesting features from the 1-D Vs models of Tokam et al. (2010) is a localized fast velocity (Vs of 3.6–3.8 km s^–1) layer in the upper crust beneath many of the CBSE seismic stations in Cameroon. The authors suggest that the HVL indicates a fair amount of heterogeneity within the upper crust both within and between regions and generally interpret them as mafic bodies that intruded into the upper crust during magmatic events of various ages. Although this may be true for the region of Cenozoic volcanism, a plausible geological interpretation still needs to be considered for the other regions. Besides, this anomalous feature needs to be investigated further since the thin-layer model parametrization for the seismic velocities in the joint inversion study of Tokam et al. (2010) may restrict the location of the discontinuities required by the data. Also, the inclusion of short-period dispersion data (<20 s) obtained from teleseismic one-station method may induce such HVL in the upper crust, as the dispersion is known to be less accurate due to intrinsic attenuation and scattering along teleseismic ray paths (Ojo et al.2018b). Based on this, the inverted models from the traditional least-squares inversions may generate a number of discontinuities in the Vs model in order to fit every small wiggle in the observed data (Shen et al.2013). For example, a satisfactory fit to the receiver function data in such inversions may result from the destructive interference of multiple phases, which are generated at many discontinuities over a range of depth rather than from the Moho or some specific intracrustal discontinuity. Therefore, we suspect that some of the localized HVL may be artifacts of inversion resulting from weak vertical smoothing constraints. Therefore, we examined the inverted 1-D models from this study for indications of this HVL at different tectonic regions and developed a new 3-D model of the upper crust in Cameroon by performing a joint inversion of the ZH ratio and short-period dispersion data.

In this study, we present improved estimates of crustal structure parameters (shear wave velocity, crustal thickness and Vp/Vs ratio) with meaningful uncertainties by performing joint inversion of multiple seismic datasets in Cameroon, West Africa. Our improvement over previous studies came from the inclusion of new geophysical dataset and the application of a new inversion methodology. As far as we know, there have been no reports of simultaneous inversion of Rayleigh wave group and phase dispersion from ambient seismic noise, P-receiver function and ZH ratio data using the Neighbourhood Algorithm (Sambridge 1999a,b) in the literature. Based on the newly developed model, we investigate the extent of magmatism-induced crustal deformation and the previously reported fast velocity layer in the upper crust in Cameroon. Finally, we consider plausible geological interpretations for the observed features in our new model and the tectonic implications for the study area.

1.1 Tectonic setting

The Cameroon region in West Africa lies in longitudes 9°–15.5°E and latitudes 2°–11°N and shares boarders with Nigeria in the west, the Central African Republic in the east, Lake Chad in the north and Equatorial Guinea, Gabon and the Republic of Congo in the south (Fig. 1). Interestingly, this region contains a number of tectonically diverse features whose geological record spans much of Earth's history (∼3 billion years). These include the Ntem Complex which represents the northern edge of the Archean Congo Craton, the Oubanguides Belt containing the Central African Shear Zone (CASZ) and other important fault/shear zones, the Garoua Rift and the Mamfe Basin which are eastern extensions of the Benue Trough from Nigeria to Cameroon, and most importantly, a chain of tertiary to recent alkaline volcanoes referred to as the Cameroon Volcanic Line (Fitton & Dunlop 1985; Fitton 1987; Noel et al.2014; Fig. 1).

Figure 1.

Map of the study area showing surface elevation, broad-band seismic stations and major tectonic features (e.g. Toteu et al.2004; Nzenti et al.2006; Reusch et al.2010; Tokam et al.2010), including the CVL (red patches), Benue Trough, Congo Craton, Manfe Basin (MB), Oubanguides Belt, Adamawa Plateau, Kumba Graben (KG), Mt Oku (MO) and Garoua Rift. The Foumban shear zone (FSZ), Tchollire-Banyo Shear Zone (TSZ), Sanaga Fault (S. Fault) and Central African shear zone (CASZ) are shown as blue lines in the figure (e.g. De Plaen et al.2014). Deep blue triangles indicate stations that were installed in 2005; light blue triangles show stations that were installed in 2006 and the dashed cyan lines show the location of profiles plotted in Fig. 15. The cyan and yellow stars along profile AA indicate the beginning and top of the observed high velocity layers.

The Archaean Congo Craton is one of the large Precambrian Cratons that make up the African continent. Its northern portion manifests in southern Cameroon and is often referred to as the Ntem Complex (Vicat et al.1996; Tchameni et al.2001; Tokam et al.2010). Although the Congo Craton has mostly been tectonically stable since it was formed, it contains some reworked and resedimented material known as protocrusts formed during the Palaeoproterozoic (Tchameni et al.2001). An intracrustal discontinuity delineated by an E–W striking dense body located around 4°N marked a possible tectonic boundary between the Congo Craton and the northern Pan-African block (Fig. 1; Toteu et al.2004).

Another important feature in Cameroon is the Oubanguides Mobile Belt which is characterized by several fault/shear zones such as the Sanaga Fault and the large Central African Shear Zone (CASZ) which trends NE–SW across the southern part of Cameroon to other neighbouring countries (Dorbath et al.1986; Castaing et al.1994; Ngako et al.2003; Toteu et al.2004 ). The Oubanguides Belt formed as a result of the collision between the São Francisco Craton, the Congo Craton and the West African Craton during the formation of Gondwana and forms part of the larger Neoproterozoic Pan-African-Brazilian Belt (Castaing et al.1994; Ngako et al.2003; Toteu et al.2004; Nzenti et al.2006; Fig. 1). The domain boundaries in this region are defined by large shear zones such as the Sanaga Shear Zone (SSZ) and the Tchollire-Banyo Shear Zone (TSZ), which are also part of the larger Precambrian Central African Shear Zone (CASZ). Likewise, the Precambrian Basement rocks that underlie most part of Cameroon are also part of the Pan-African Belt assembled during the Pan-African Orogeny (∼600 Ma) and consist mainly of schists and gneisses intruded by granites and diorites (Noel et al.2014; Dawaï et al.2017). The Cretaceous sediments are mostly sandstones with small amounts of limestone and shales in the coastal plain, and the basins are of variable thicknesses of 0–5 km across the study area (Fig. 1; Fitton 1987; Déruelle et al.1991, 2007; Laske et al.2013).

The Benue Trough is a NE–SW-trending Mesozoic Rift Basin that extends from Niger Delta (Gulf of Guinea) to Lake Chad (northern Nigeria) with eastward extension known as the Garoua Rift and Mamfe Basin into Cameroon (Fig. 1). The Benue Trough is the third arm of a failed triple junction and a continental rift zone that is thought to have formed during the opening of the South Atlantic Ocean in the Cretaceous times beginning at approximately 130 Ma and continued until 84 Ma, followed by a period of compression (Burke et al.1971; Guiraud & Maupin 1992; Gallacher & Bastow 2012). Researchers have suggested that the similar Y-shape of the Benue Trough and the CVL may imply a common geodynamic control which is attributed to the northeast-trending Pan-African dextral shear zones (Fitton 1987; Guiraud & Maurin 1992). Likewise, the alkali basalts of the Benue Trough that is geochemically and isotopically similar to those of the CVL corroborate this inference (Coulon et al.1996).

A unique feature that has fascinated scientists for decades is the Cenozoic Cameroon Volcanic Line (CVL) which consists of a chain of alkaline volcanoes that runs from the Pagalu Island in the southwest (oceanic sector) to Garoua Rift in the northernmost part of Cameroon (continental sector) with no clear age progression (Fig. 1, Fitton 1987; Lee et al.1994). The oceanic sector includes the islands of Pagalu (5 Ma), São Tomé (14 Ma), Príncipe (31 Ma) and Bioko (1 Ma) while the continental sector includes the Cameroon (<3 Ma), Manengouba (1 Ma), Bamboutos (21–14 Ma), Oku (31–32 Ma) and Mandara mountains (32 Ma) as well as the Adamawa (11–7 Ma) and Biu (<5 Ma) plateaus (Fig. 1; Toteu et al.2001). Along the CVL, magma generation and intrusion have been sustained for more than 70 million years over a 1600-km-long chain straddling the continent–ocean boundary. Mt Cameroon is located at the centre of the line and is regarded as the only active volcano with an observed eruption in the year 2000 (Tabod 1991; Tabod et al.1992; Suh et al.2003; Déruelle et al.2007; Ateba et al.2009). The lack of age progression among other factors has subjected the origin of the CVL to considerable debates (e.g. Fitton & Dunlop 1985; Fitton 1987; Poudjom Djomani et al.1997; Déruelle et al.2007), including plume and plume-free models which have been suggested for its origin (e.g. King & Anderson 1995, 1998; Meyers et al.1998; King & Ritsema 2000; Burke 2001; King 2007; Reusch et al.2010, 2011; Koch et al.2012; Gallacher & Bastow 2012; Milelli et al.2012; De Plaen et al.2014; Adams et al.2015). An example of a plume model involves the flow of plume material from the Afar region guided to the CVL by thinned lithosphere beneath the central Africa rift system (Ebinger & Sleep 1998). Plume-free models involve decompression melting beneath reactivated shear zones (Fairhead & Binks 1991) or small-scale mantle convection resulting from edge flow along the northern boundary of the Congo Craton (King & Ritsema 2000). Along the volcanic line, geologists have found several occurrences of mantle-derived (ultramafic) xenoliths in the basaltic lavas suggesting that the formation of the CVL is most likely related to a mantle process (Princivalle et al.2000; Deruelle et al.2007). Consequently, mostly alkaline basalts of CVL may have modified the overlying crust significantly if the present-day melt is retained in the crust and/or cooled mafic intrusions like other known hot spots in the world. On the other hand, we will expect a limited impact on the crust if the silica-poor magmas truly fractionate in the mantle rather than the crust prior to eruption (Fitton 1980; Suh et al.2003; Gallacher & Bastow 2012).

1.2 Previous crustal structure studies

Early constraints on the crustal structure in Cameroon came from the analysis of gravity data (e.g. Poudjom et al.1995; Nnange et al.2000; Toteu et al.2004; Tadjou et al.2009; Shandini et al.2010; Basseka et al.2011) and active/passive seismic datasets (e.g. Stuart et al.1985; Dorbath et al.1986; Tabod 1991; Tabod et al.1992; Plomerova et al.1993). Recently, several new and detailed seismological studies have emerged since the deployment of a 32-station broad-band seismic experiment across Cameroon (CBSE) from 2005 to 2007. These include the joint inversion of receiver function and teleseismic group velocity dispersion (Tokam et al.2010), the estimates of crustal thickness and Vp/Vs ratio (Gallacher & Bastow 2012), the Ambient Noise Tomography (Guidarelli & Aoudia 2016), and the more recent crustal radial anisotropy study (Ojo et al.2017). Previous non-seismological studies provided us with initial estimates of crustal thickness in Cameroon. For example, Stuart et al. (1985) revealed a crustal thickness of 33 km beneath the Adamawa Plateau and a thinner crust of 23 km beneath the Garoua Rift which are also in agreement with estimates from gravity studies (e.g. Poudjom et al.1995; Nnange et al.2000). Using the CBSE data, Tokam et al. (2010) investigated the crustal structure in Cameroon more broadly and reported a crustal thickness of 35–39 km beneath the CVL and the Oubanguides Belt and a relatively thicker crust of 43–48 km beneath the Congo Craton. The thickened crust was attributed to collisional tectonic activities that characterized the formation of Gondwana, whereas the relatively thinner crust featuring slow shear wave velocity along the CVL was used as evidence for crustal modification beneath the CVL. The authors also reported the existence of a thin mafic layer (Vs of 3.6–3.8 km s^–1) in the upper crust beneath most of the CBSE seismic stations (Fig. 1) and a relatively fast lower crustal shear wave velocity (Vs > 4.0 km s^–1). In addition to crustal thickness estimates that were generally in agreement with those of Tokam et al. (2010), Gallacher & Bastow (2012) observed markedly low Vp/Vs ratios (network average ∼1.74) across Cameroon. Surprisingly, this value is relatively low and atypical of hotspot regions. Therefore, the authors argued against the presence of either melt or cooled mafic crustal intrusions within the Cameroon crust due to CVL magmatism. Using the continuous ambient seismic noise data recorded by the CBSE array, Guidarelli & Aoudia (2016) presented a number of period-dependent group velocity maps that predominantly revealed a very slow velocity along the CVL and fast velocity beneath the Congo Craton. Based on their results, the authors favored processes combining small-scale upwelling at the edge of a thick lithosphere and reactivation of Precambrian basement structures to explain the distribution of Holocene to recent magmatism and plateau uplift in Cameroon. The 3-D radial anisotropy model of Ojo et al. (2017) revealed negative radial anisotropy (Vsv > Vsh) beneath Mt Oku, the Garuoa Rift and the CVL that is likely due to aligned microcrack networks related to rifting and magmatism. Likewise, positive radial anisotropy (Vsh > Vsv) is observed within the Oubanguides Belt and attributed to mineral alignment due to shearing in the fault zone. The authors also observed a widespread negative radial anisotropy in the middle crust that is likely caused by flow-induced alignment of anisotropic minerals that crystallized during magma intrusion. Although their results were in support of small-scale convection processes for the development of the CVL, the authors recommended further studies involving azimuthal anisotropy to obtain a complete picture of the stress field and geodynamic evolution in Cameroon (e.g. Ojo et al.2018b).

2 DATA AND METHODOLOGY

We used data from the 32 broad-band seismic stations belonging to temporary networks deployed and managed by the AfricaArray program and archived at IRIS-DMC (Wiens & Nyblade 2005; Fig. 1). The Cameroon Broadband Seismic Experiment (CBSE) is a 2-yr project from January 2005 to February, 2007 that primarily aimed at evaluating several geodynamic models proposed for the origin of Cameroon volcanic line (Fig. 1). These stations are almost evenly distributed across the entire country with interstation spacing of 50–150 km. Using ambient seismic noise and earthquake event data, we obtained local dispersion curves using the ambient noise tomography method, stacks of distance-moveout corrected P-wave receiver functions and Rayleigh wave ellipticity (ZH ratio) measurements at the 32 stations. Although our previous studies have shown the existence of radial and azimuthal anisotropy in the Cameroon crust, we assumed an isotropic case in this study and we did not investigate or remove the effect of anisotropy from the aforementioned data set (Ojo et al.2017, 2018b).

2.1 Rayleigh wave dispersion data

Following the standard ANT method (e.g. Bensen et al.2007, 2008; Lin et al.2008), we obtained inter-station Noise Correlation Functions (NCFs) by processing continuous records of ambient seismic noise data as detailed in Ojo et al. (2017). Using the NCFs, we employ the multiple phase-matched filter method of Herrmann (2013) to retrieve fundamental-mode Rayleigh wave group velocity dispersion in the period range of 5–30 s (Dziewonski et al.1969; Herrmann 1973; Levshin et al.1989; Liang & Langston 2008). Similarly, we used the automated frequency-time analysis method (AFTAN; Levshin et al.1972; Bensen et al.2007) to obtain fundamental-mode Rayleigh wave phase velocity dispersion in the period range of 4–40 s. To discriminate bad dispersion measurements, we impose a wavelength (1λ) criterion and specify a signal-to-noise ratio (SNR) greater than 10 in addition to manual inspection (e.g. Luo et al.2015). Subsequently, we used the fast-marching surface wave tomography (FMST) method of Rawlinson & Sambridge (2005) to estimate the Rayleigh wave group and phase velocity variations at different periods across the study area at 0.5° by 0.5° grid spacing (e.g. Ojo et al.2017). Then, we extracted pure path group and phase dispersion curves at the grid nodes closest to the seismic stations. Based on the tomography method employed in this study, the measurement uncertainties cannot be directly quantified. However, in our previous studies, we estimated the uncertainty in our phase velocity dispersion measurements to be about 10 m s^–1 and we expect significantly larger uncertainties at longer periods for group velocities (e.g. Luo et al.2013; Shen & Ritzwoller 2016; Ojo et al.2017).

2.2 P-Wave receiver functions

We downloaded seismic waveform data of 423 teleseismic earthquakes that occurred during the lifetime of the CBSE stations from IRIS-DMC (Trabant et al.2012). The earthquakes are chosen with epicentral distances in the range of 30°–90° and magnitude m_b ≥ 5.5 (Fig. 2). The seismograms are windowed between 20 s before and 100 s after the P onset time and we removed data with insufficient data points. Then, we rotate the seismograms into radial and transverse components, and apply a zero-phase Butterworth bandpass filter from 0.04 to 2 Hz. We then compute the radial receiver function using the time-domain iterative deconvolution method of Ligorria & Ammon (1999) for a Gaussian width of 2.5 (pulse length ∼1 s). Furthermore, we select high-quality and consistent receiver functions following an automatic procedure that stacks all the receiver functions at each station and measure the cross correlation (CC) between the stack and individual receiver functions. We discard receiver functions with CC values less than 0.8, those without positive amplitude around time 0.5–2 s or abnormal amplitude above 1.0. We correct the selected receiver functions for distance moveout based on the IASP91 earth model (Kennett & Engdahl 1991) aligning them to an epicentral distance of 60° (Liu & Niu 2011). The distance-moveout corrected receiver functions are then stacked linearly to obtain the mean and standard deviation which is used as proxy for the uncertainty in the receiver function data.

Figure 2.

Earthquakes used in this study for the computation of receiver function (blue dots with epicentral distance between 30° and 90°) and Rayleigh wave ellipticity measurements (red dots with epicentral distance between 45° and 165°). The blue triangle is the study area.

2.3 Rayleigh wave ellipticity (ZH RATIO)

Seismograms from teleseismic earthquakes of magnitude m_b ≥ 5.0 in the distance range 45°–165° were obtained from the IRIS DMC for the period of January 2005 to February 2007 (Fig. 2). We chose events that occurred with a depth of less than 100 km to ensure that the surface waves are well excited. We decimate the data and remove mean, linear trend, and instrument response. We rotate the horizontal components into the radial and transverse components. To measure the ZH ratio, we followed a procedure similar to that described by Tanimoto & Rivera (2008). We define the surface wave window (half length of 4T, where T is period) and cross correlate vertical and Hilbert-transformed radial components to measure their similarity and ascertain the presence of Rayleigh wave. We apply a three-wavelengths criterion and compute the Gaussian width using a subroutine from the Computer Programs in Seismology (Herrmann 2013). We then apply a series of narrowband Gaussian filters to the radial and vertical components to get the envelope of both components in the frequency domain (Herrmann 1973). Subsequently, we compute the ZH ratio by taking the ratio of the maximum amplitudes of the envelopes of the Gaussian filtered vertical component and Hilbert transformed radial component as defined in eq. (1).

\begin{eqnarray*} \eta \ \left( T \right) = \frac{{Z\left( T \right)}}{{H\left( T \right)}}, \end{eqnarray*}

(1)

where T is the period. According to Tanimoto & Rivera (2005), the inverse ratio is preferable to avoid situations where the ratio goes to infinity. The ZH ratio is measured for each earthquake and at each period based on the assumption of an isotropic medium, and we average the result over azimuths for each period and station. To select high-quality ZH ratios, we used waveforms with cross correlation value greater than 0.8 and with S/N ratio on both the radial and the vertical components greater than 10. We also made sure that the phase difference between the radial and vertical components is less than 22.5° (π/8 rad) and we limit the maximum and minimum values of ZH ratios to 3 and 0.2, respectively (Tanimoto & Rivera 2008; Lin et al.2012). Furthermore, we discard measurements that lay outside of the ±2σ bound (σ is the standard deviation of the measured ZH ratio data) and recomputed the mean and the standard deviation for periods with reliable ZH ratios measurements greater than 10.

2.4 Joint inversion scheme

2.4.1 Joint inversion for 1-D Vs model

To perform a robust joint inversion of the multiple datasets, we developed a new program based on the NA of Sambridge (1999a,b). NA is a global optimization algorithm that is well suited for solving highly non-linear geophysical inverse problems. Unlike the conventional linearized method, it avoids the calculation of partial derivatives and randomly samples the multidimensional parameter space in a self-adaptive and efficient manner using the Voronoi cells. The final result of the inversion process is an ensemble of optimal models that minimize the objective function with an acceptable level of fit to the measured data and uncertainties (e.g. Berbellini et al.2016, 2017; Ojo et al.2018a). By further Bayesian analysis, the NA method can also provide information about resolution and trade-offs between model parameters. Another advantage of this method is its independence of a good initial model, but rather a random initial model is generated based on the user-supplied range of model parameters. The method naturally avoids the problem of misfit function scaling thereby allowing any type of user-defined misfit measure to be employed (Sambridge 1999a,b; Sambridge & Mosegaard 2002). Compared to other global optimization methods, the NA is easily tuned (only two control parameters) and faster since it uses the rank of the misfit function to compare models rather than the misfit values themselves. This significantly reduces the number of forward computation and leads to an overall reduction in the computation time (Sambridge 1999a, b).

We compute theoretical surface wave data (i.e. Rayleigh wave group and phase velocity dispersion and ZH ratio) using the efficient normal-mode formalism of Herrmann (2013) for 1-D layered earth and synthetic receiver function using the method of Shibutani et al. (1996) for each model sample (Ferreira & Woodhouse 2007). Both forward modelling approaches are computationally efficient, making the NA a suitable inversion strategy. We define the objective function as a sum of the square of the misfit between the observed and computed data for each data type:

\begin{eqnarray*} {M_{joint}} = \mathop \sum \limits_{j = 1}^4 {w_j}\ \cdot \mathop \sum \limits_{i\ = \ 1}^{{N_j}} \frac{{{{\left[ {D_{ij}^{obs} - {G_{ij}}\left( m \right)} \right]}^2}}}{{\sigma _{ij}^2}},j \end{eqnarray*}

(2)

where j = 1…4 refers to the four types of data (receiver function, group velocity dispersion, phase velocity dispersion and Rayleigh wave ellipticity) used in the inversion, w_j is the weight of each data, |${G_{ij}}( m )$| and |$D_{ij}^{obs}$| are the synthetic and measured data, respectively, σ_ij is the measurement error, and N_j is the total number of measurement of a given data type. Since we do not have any preferential usage of the datasets, we set the weight to be equal for all the four types of data to allow equal contributions of the observables in constraining the inverted model (w₁ = w₂ = w₃ = w₄ = 0.5). After a series of synthetic tests, we set the two tunable NA parameters that control how explorative or exploitative the inversion proceeds to be 100 and 60, respectively. This implies that we generate 100 new models (Ns = 100) at each iteration and add it to the population of existing fitting models. Out of this, the best 60 models (Nr = 60) are chosen to produce a new distribution of Voronoi cells and the process continues iteratively until a maximum of 500 iterations and 50 500 model evaluations are performed. After several trials, we conclude that these values are good enough to allow a robust inversion that result in a well constrained model.

Since our goal is to study the crustal structure, we parametrize the 1-D model beneath each station with five crustal layers and one layer in the uppermost mantle. Each layer is assigned a range of thickness, shear wave velocity at the top, shear wave velocity at the bottom and a range of layer Vp/Vs ratio making a total of 24 free parameters to invert beneath each station (Table 1). We determine the density of each layer using the relation of Brocher (2005) in order to constrain the Vp/Vs ratio better without inverting for density. We also apply a physical dispersion correction (Kanamori & Anderson 1977) using the Q-model from AK135 in the crust (Kennett et al.1995). We neglect the effect of surface topography in our 1-D inversions, implying that zero depth has different elevations at different station locations.

Table 1.

Range of model parameters.

Layer	Thickness (km)	Top Vs (km s–1)	Bottom Vs (km s–1)	Vp/Vs ratio	Qp	Qs
Crustal Layer-1	0–2	1.50–4.00	1.50–4.00	2.00–3.00	100	25
Crustal Layer-2	0–3	1.50–4.00	1.50–4.00	1.65–2.00	675	300
Crustal Layer-3	1–15	2.60–4.50	2.60–4.50	1.65–1.80	1450	600
Crustal Layer-4	5–20	3.20–4.50	3.20–4.50	1.65–1.80	1450	600
Crustal Layer-5	5–20	3.20–4.50	3.20–4.50	1.65–1.80	1450	600
Mantle Layer-1	5–30	4.00–5.00	4.00–5.00	1.70–1.90	1450	600

The ranges of model parameters are carefully chosen according to previous studies and typical values from the global AK135 model (Kennett et al.1995). We invert for 1-D Vs models beneath the 32 CBSE stations constrained by a portion of the receiver function (time window 0–10 s), group velocity dispersion (period 5–30 s), phase velocity dispersion (period 4–40 s) and Rayleigh wave ellipticity (period 1062 s) measurements. The final 1-D Vs model is the mean of all accepted ensemble of models for each station and the standard deviation is used to quantify the model uncertainty. The crustal thickness is obtained by computing the sum of the five inverted crustal layer thicknesses for each model and the mean of all accepted ensemble of models for each station is taken as the final Moho depth within the bounds determined by the standard deviation. Although we refer to the inverted model as an isotropic Vs model, it is actually a vertically polarized shear wave speed (Vsv) model since Rayleigh wave data cannot constrain the horizontally polarized shear wave speed (Vsh).

2.4.2 Joint inversion for 3-D Vs model

To develop a new 3-D model of the uppermost crustal structure in Cameroon, we interpolate the ZH ratio measurements (Section 2.3) at spatial grid spacing of 0.5° x 0.5° across the study area and extract local variation of ZH ratio with period at these grid points. Our choice of grid size is based on previous resolution test results (e.g. Ojo et al.2017, 2018b). Interpolating the ZH ratio measurement across the study area may serve to minimize local differences and artifacts that may be due to the uncertainty of the ZH ratio estimation. However, the resolution of the map will be in the order of the inter-station spacing (∼50–150 km). Subsequently, we extracted pure path dispersion curves at the same grid points from our previous 2-D group and phase velocity tomographic maps in the period range of 6–30 s (e.g. Ojo et al.2017). We then performed a joint inversion of the ZH ratio, group and phase velocity dispersion data at each of these grid points for 1-D Vs model using the method described in Section 2.4.1. Using ordinary kriging method (Davis 2002), we combined the entire 1-D Vs model into a new 3-D Vs model from the surface down to 60-km depth.

3 RESULTS

Using NA program, we obtained a library of 1-D Vs, Vp/Vs ratio and crustal thickness estimates for stations in the CBSE network in Cameroon. According to previous studies, stations CM08 and CM14 experienced instrument failures during the deployment and may not record enough high quality waveform for data analysis (Tokam et al.2010). Gallacher & Bastow (2012) reported that they could not obtain receiver function results for stations CM08, CM14, CM15, CM18, CM19, CM28 and CM31; while Tokam et al. (2010) also failed to obtain a good data fit for stations CM03, CM05, CM09, CM15 and CM23 in their studies. Also in this study, no result was obtained for station CM08 and we did not obtain a satisfactory data fit at stations CM09 and CM31. To present the results obtained in this study, we divide the entire study area into six tectonic regions of interest, namely: Coastal Plain, Congo Craton, Oubanguides Belt, CVL, Adamawa Plateau and Garoua Rift (e.g. Tokam et al.2010). Table 2 provides a list of the broad-band seismic stations located in each of these tectonic subdivisions (Fig. 1).

Table 2.

CBSE seismic stations grouped by tectonic terrains.

S/N	Tectonic regions	Stations
1	Coastal Plain	CM01, CM05
2	Congo Craton	CM02, CM04, CM06, CM07, CM11
3	Oubanguides Belt	CM03, CM10, CM12, CM14, CM17
4	Cameroon Volcanic Line (CVL)	CM09, CM13, CM15, CM16, CM18, CM19, CM20, CM23
5	Adamawa Plateau (Upper CVL)	CM21, CM22, CM24, CM25, CM26, CM27
6	Garoua Rift	CM28, CM29, CM30, CM31, CM32

The results of the inverted crustal parameters (Vs, Vp/Vs ratio and Moho depth) are summarized in Table 3 based on the different lithospheric subdivisions (Table 2) and depth in the crust across Cameroon. We also present plots of the variation of average shear wave velocities and Vp/Vs ratios in the upper, middle and lower crust in different tectonic regions in Figs 3 and 4, respectively. Additionally, we present the joint inversion result for all the CBSE stations and a summary table of the data sampling in the Supporting Information.

Figure 3.

Mean shear wave velocity in the upper (red dots and dashed lines), middle (cyan dots and dashed lines) and lower crust (red dots and dashed lines) at different tectonic regions, namely Cameroon Volcanic Line, Adamawa Plateau, Congo Craton, Garoua Rift, Coastal Plain and Oubanguides Belt.

Figure 4.

Mean Vp/Vs Ratios in the upper (red dots and dashed lines), middle (cyan dots and dashed lines) and lower crust (red dots and dashed lines) at different tectonic regions, namely Cameroon Volcanic Line, Adamawa Plateau, Congo Craton, Garoua Rift, Coastal Plain and Oubanguides Belt.

Table 3.

Summary of inverted crustal parameters beneath the CBSE seismic stations.

		Shear wave velocity (km s–1)				Vp/Vs ratio				Vs ≥ 4.0 km s–1	Vs ≤ 3.0 km s–1	Vs ≥ 3.6 km s–1
Tectonic Terrain/Stations	Crustal thickness (km)	UC	MC	LC	CA	UC	MC	LC	CA	LC (km)	UC (km)	km
Coastal Plain
(1) CM01	40.53	3.52	3.95	4.09	3.85	1.75	1.71	1.73	1.73	9.5	0.5	7.5
(2) CM05	36.94	3.45	3.97	4.01	3.81	1.74	1.69	1.71	1.71	5.5	0.5	8
Terrain Average	38.735	3.485	3.96	4.05	3.83	1.745	1.7	1.72	1.72	7.5
Congo Craton
(1) CM02	40.56	3.62	3.86	4.01	3.83	1.74	1.72	1.72	1.73	8		2.5
(2) CM04	41.42	3.64	3.87	4.11	3.87	1.73	1.72	1.74	1.73	11.5		3
(3) CM06	39.03	3.6	3.85	4.08	3.84	1.74	1.71	1.72	1.73	9.5		4.5
(4) CM07	39.9	3.63	3.82	4.05	3.83	1.78	1.72	1.73	1.74	8.5		1
(5) CM11	42.43	3.6	3.85	4.05	3.83	1.78	1.72	1.73	1.74	8.5		2
Terrain Average	40.668	3.618	3.85	4.06	3.84	1.754	1.718	1.728	1.734	9.2
Oubanguides
(1) CM03	36.66	3.57	3.86	4.17	3.86	1.75	1.71	1.76	1.74	8.5		2.5
(2) CM10	35.74	3.5	3.94	3.94	3.79	1.74	1.7	1.73	1.72	1		9.5
(3) CM12	37.57	3.59	3.93	4.15	3.89	1.75	1.71	1.75	1.74	3		8.5
(4) CM14	40.43	3.53	3.83	4.07	3.81	1.79	1.72	1.76	1.76	5		0.5
(5) CM17	35.22	3.52	3.86	4.31	3.9	1.73	1.71	1.79	1.74	6		9
Terrain Average	37.124	3.542	3.884	4.128	3.85	1.752	1.71	1.758	1.74	4.7
CVL
(1) CM09 *MC	38.62	3.25	4.13	4.23	3.86	1.81	1.69	1.73	1.75	12.5	4	8
(2) CM13 *kG	30.4	3.43	3.76	3.95	3.71	1.73	1.67	1.7	1.7	2.5	0.5	8
(3) CM15	28.96	3.48	3.85	3.92	3.74	1.74	1.68	1.69	1.7	0.5	0.5	3
(4) CM16	35.52	3.45	3.86	4.01	3.77	1.74	1.71	1.74	1.73	4.5	0.5	9
(5) CM18 *MB	33.19	3.37	3.88	4.11	3.78	1.75	1.73	1.75	1.74	10.5	1	7.5
(6) CM19	37.58	3.52	3.8	4	3.77	1.73	1.71	1.74	1.73	5.5		8
(7) CM20	42.21	3.53	3.87	3.99	3.79	1.74	1.73	1.74	1.74	3.5		9
(8) CM23	35.78	3.45	3.86	3.92	3.74	1.74	1.72	1.76	1.74	1	0.5	8.5
Terrain Average	35.2825	3.435	3.8763	4.0163	3.77	1.7475	1.705	1.7313	1.7288	5.0625
Adamawa
(1) CM21	33.76	3.49	3.84	3.94	3.76	1.78	1.71	1.74	1.74	2	0.5	2.5
(2) CM22	37.03	3.46	3.83	4	3.76	1.74	1.71	1.73	1.72	4		9
(3) CM24	37	3.53	3.75	3.94	3.74	1.75	1.72	1.73	1.73	2.5		8.5
(4) CM25	41.12	3.53	3.8	4	3.77	1.75	1.73	1.74	1.74	4.5		11
(5) CM26	34.36	3.51	3.78	3.91	3.73	1.75	1.72	1.75	1.74	2		6.5
(6) CM27	33.86	3.46	3.79	4	3.75	1.74	1.7	1.73	1.72	4.5		9.5
Terrain Average	36.18833333	3.4967	3.7983	3.965	3.7517	1.7517	1.715	1.7367	1.7317	3.25
Garoua Rift
(1) CM28	28.44	3.48	3.63	3.93	3.67	1.75	1.7	1.71	1.72	2.5		11.5
(2) CM29	28.29	3.31	3.82	3.92	3.68	1.77	1.69	1.7	1.72	1.5	1	7
(3) CM30	32.65	3.49	3.73	4.11	3.77	1.75	1.7	1.74	1.73	6.5		7
(4) CM31	36.64	3.47	3.88	4.09	3.82	1.76	1.7	1.71	1.72	11	0.5	1
(5) CM32	35.97	3.53	3.82	3.95	3.76	1.78	1.73	1.74	1.75	3.5		9
Terrain Average	32.398	3.456	3.776	4	3.74	1.762	1.704	1.72	1.728	5
Network Average	36.38	3.4968	3.8384	4.009	3.7781	1.7535	1.7087	1.7281	1.7297	1.72968

CT: Crustal Thickness; UC: Upper Crust (1/3CT); MC: Middle Crust (1/3CT-2/3CT); LC: Lower Crust (2/3CT-CT); MB: Mamfe Basin; KG: Kumba Graben; MC: Mt. Cameroon. Grey row are stations with poor data fit (CM09, CM31) and Instrument Failure (CM14).

3.1 Coastal plain

The two CBSE stations located in this region have an average crustal thickness of 38.7 km (Table 3). In the uppermost crust (<1 km), the shear wave velocities at both stations are slow (≤3.0 km s^–1), but CM01 revealed a much lower shear wave velocity at depth when compared to CM05. Subsequently, the shear wave velocity increases monotonically with depth at a relatively fast rate exceeding 3.6 km s^–1 at around 8 km and ultimately reaching a velocity of 4.0 km s^–1 at CM01 and 3.79 km s^–1 at CM02 within the upper crust (Supplementary Material). In the middle crust, the shear wave velocity is similar between the two stations and relatively constant with a mean value of 3.96 km s^–1 (Fig. 3) The mean shear wave velocity in the lower crust is 4.0 km s^–1 (Table 3, Fig. 3). However, the thickness of crustal layers having shear wave velocity of 4.0 km s^–1 or higher in the lower crust is significantly greater beneath CM01 (∼9.5 km) compared to CM05 (∼5.5 km). The Vp/Vs ratio at both stations is high (>1.8) in the shallow crust (∼ 1 km) and has an average value of 1.7 and 1.72 in the middle and lower crust, respectively (Fig. 4). The Vp/Vs ratio is slightly higher (∼0.02) at CM01 than CM05 throughout the crust (Table 3, Fig. 4).

3.2 Congo craton

The crustal thickness at the five seismic stations within the Congo Craton ranged from ∼39–44 with an average value of ∼41 km for the region (Table 3). The 1-D model for the region is characterized by shear wave velocities greater than 3.6 km s^–1 within the top 5 km. However, a localized high velocity layer can be observed at stations CM07 and CM11 where a slight drop in the shear wave velocity at depths ranging from 3 to 7 km creates the appearance of a low velocity zone below the peak shear wave velocity (3.63–3.68 km s^–1) at around 2-km depth (Supporting Information). Subsequently, the shear wave velocity increases monotonically with depth with average values of 3.85 and 4.06 km s^–1 in the middle and lower crust, respectively (Fig. 3). The thickness of crustal layers with shear wave velocity exceeding 4.0 km s^–1 in the lower crust ranged from 8 km at station CM02 to 11.5 km at CM04 in this region (Table 3). In the upper crust, there are obvious differences in the Vp/Vs ratio between stations within this region, and the stations with relatively fast shear wave velocities are associated with lower Vp/Vs ratios. Stations CM07 and CM11 which are located more in the northern part of the region are characterized by an average Vp/Vs ratio of 1.78 while the remaining southernmost stations have slightly lower average Vp/Vs ratio of 1.73 (Fig. 1, Tables 2 and 3, Fig. 4). The Vp/Vs ratios seem constant in the middle crust at all the seismic stations within the region and slightly increase from an average value of 1.72 to 1.73 in the lower crust (Table 3, Fig. 4).

3.3 Oubanguides belt

The crustal thickness beneath the CBSE stations in the region ranged from 35 to 40 km with an average value of 37 km (Table 3). The shear wave velocity models for stations within the region exhibit slightly different characteristics in the crust. For example, at shallow depths (≤5 km), stations CM10 and CM17 have a slightly lower Vs (<3.5 km s^–1), while stations CM03, CM12 and CM14 show an indication of localized high velocity layer with shear wave velocities ≥3.6 km s^–1 (Supporting Information). The latter stations subsequently reveal a low-velocity zone beneath the HVL in the upper crust. This feature is more obvious at station CM14 with a localized low-velocity layer (<3.4 km s^–1) at the base of the upper crust. The remaining stations have shear wave velocity increasing with depths up to about 4.0 km s^–1 at station CM10 in the middle crust and increased to 4.31 km s^–1 at station CM17 in the lower crust (Table 3, Fig. 3). At station CM10, the shear wave velocity is relatively lower in the lower crust with average value of 3.9 km s^–1 (Table 3, Fig. 3). Beside station CM03 which is located more in the southern part of the region (Fig. 1), the thickness of the crustal layers with shear wave velocity exceeding 4.0 km s^–1 is ≤6 km and comparatively thinner than other regions like the Congo Craton and the Coastal Plain. The Vp/Vs ratio is similar among the stations in the region with average values ranging from 1.75 in the upper crust to 1.71 in the middle crust and subsequently increases to 1.76 in the lower crust (Table 3, Fig. 4).

3.4 Cameroon volcanic line (CVL)

The crustal thickness in the southwestern part of the CVL ranges from ∼29 km at station CM15 to 42 km at CM20, with an average value of 35 km within the region. In addition to station CM15, the crustal thickness is relatively small (30–33 km) in the Mamfe Basin (CM18) and Kumba Graben (CM13). The observed data beneath Mt Cameroon (CM09) are quite complex and the inverted model has relatively poor data fit. The shear wave velocity in the uppermost crust (≤5 km) for the stations in this region is anomalously low (<2.8 km s^–1) especially at stations CM15, CM18 (Mamfe Basin) and CM09 (Mt Cameroon). The thickness of upper crustal layers with shear wave velocities ≤3.0 km s^–1 reduces from 4 km at CM09 to ≤0.5 km as we move towards the northernmost stations (Fig. 1, Table 3). Subsequently, the shear wave velocity increased at a relatively fast rate with depth, reaching about 3.8 km s^–1 at station CM15 in the upper crust (Fig. 3). The shear wave velocity in the middle crust is similar for most stations with an average value of 3.87 km s^–1, but there are obvious differences in the lower crust (Fig. 3). At stations CM13 (Kumba Graben) and CM15, there is a sharp increase in shear wave velocity (∼3.8–4.4 km s^–1) at relatively small depth interval (27–35 km); while at other stations, we observe a gradual increase in velocity with depth (Supporting Information). Generally, the average shear wave velocity increased from 3.87 km s^–1 in the middle crust to ∼4.0 km s^–1 in the lower crust (Fig. 3). Apart from stations CM09 and CM18, the thickness of the crustal layers with shear wave velocity exceeding 4.0 km s^–1 is thin (<5 km) compared to other regions. The average Vp/Vs ratio ranged from 1.75 in the upper crust to 1.70 in the middle crust and subsequently increases to 1.73 in the lower crust (Table 3, Fig. 4). When compared to the average values at these crustal depths, we found that the Vp/Vs ratio is relatively higher for station CM09 (Mt Cameroon) in the upper crust and lower for stations CM13 (Kumba Graben) and CM15 in the middle crust (Table 3, Fig. 4).

3.5 Adamawa plateau (Northern CVL)

Stations in the northern part of the CVL which contains the Adamawa Plateau have crustal thicknesses ranging from 33.7 km at CM21 to ∼41 km at CM25, with a regional average of ∼36 km which is close to that found in the southern part of the CVL. In the top 6 km, the shear wave velocity is less than 3.6 km s^–1 at all the stations except at CM21 which reveals an anomalously low velocity (<2.8 km s^–1) with a thickness of ∼0.5 km (Table 3). Also, we observe a localized high-velocity layer with maximum velocity of 3.7 km s^–1 at depths above 10 km. The shear wave velocity increases with depth into the middle crust at all stations with an average value of 3.8 km s^–1 (Table 3, Fig. 3). The shear wave velocities at most of the stations reveal the presence of a low-velocity zone (≤4.0 km s^–1) and fast velocity gradient (δVs ≥ 0.2 over a depth of ∼3 km or more) at stations CM21, CM22, CM24 and CM26 in the lower crust (Supporting Information). The average thickness of lower crustal layers with shear wave velocity exceeding 4.0 km s^–1 is significantly reduced (≤3.3 km) compared to all other regions. The Vp/Vs ratio is similar for all stations in the region in the upper and middle crust with average values of 1.75 and 1.71, respectively, but differs in the lower crust (Table 3, Fig. 4). Stations CM21, CM26 and CM27 revealed a much higher Vp/Vs ratio (∼4.5 per cent) compared to other stations at the base of the crust (Fig. 4).

3.6 Garoua rift

The crustal thickness at station CM28 located to the south of the Rift and station CM29 near the Rift is about 28 km and comparatively thinner than stations to the northern part of the Rift (CM30, CM31 and CM32) whose average crustal thickness is 35 km. Also, the average crustal thickness in this region (∼32 km) is smaller than other regions in Cameroon. There are obvious differences between the inverted shear wave velocity models for stations located in this region. The Vs model at Stations CM29 and CM31 reveals an anomalously low velocity (<3.0 km s^–1) with a thickness of ∼1 km in the uppermost crust and the shear wave velocity subsequently exceeds 3.6 km s^–1 at depth of ∼10 km for most stations within the region (Table 3, Supporting Information). In the middle crust, the average shear wave velocity beneath station CM28 is relatively low (∼3.63 km s^–1), while the average Vs beneath the Garoua Rift (CM29, CM31 and CM32) is relatively fast (3.82–3.88 km s^–1, Table 3, Fig. 3). In the lower crust, the average Vs beneath CM30 and CM31 is above 4.1 km s^–1, while the average Vs beneath other stations in the region are relatively low (∼3.93 km s^–1, Fig. 3). Likewise, the average thickness of layers with shear wave velocity exceeding 4.0 km s^–1 in the lower crust is significantly thicker (6.5–11 km) beneath stations CM30 and CM31 and relatively thin (≤2.5 km) beneath stations CM28, CM29 and CM32 (Table 3). The Vp/Vs ratio in the upper crust has slightly higher average values (1.77–1.78) beneath stations CM29 and CM32; while in the middle crust, CM29 exhibits the lowest Vp/Vs ratio within the region (∼1.69) and station CM32 reveals a slightly higher Vp/Vs ratio of 1.73 (Table 3, Fig. 4). Other stations are characterized by average Vp/Vs ratios of 1.75, 1.70 and 1.71 in the upper, middle and lower crust, respectively (Table 3, Fig. 4).

3.7 Variation of Vp/Vs ratio with depth and tectonic domain

Since this is the first depth-dependent Vp/Vs ratio estimates in Cameroon, we examine it more closely to establish how it varies as a function of depth and tectonic domain (Fig. 4). The mean Vp/Vs ratio in the upper crust ranged from 1.762 in the Garoua Rift to 1.745 in the Coastal Plain with a mean difference of ∼0.02 across the entire study area (Table 3, Fig. 4). The mean Vp/Vs ratios along the CVL are similar to those in the Coastal Plain (∼1.74) and so are the mean Vp/Vs ratios in the Adamawa Plateau comparable to those in the Oubanguides Belt and Congo Craton (∼1.75). Deeper into the middle crust, the mean Vp/Vs ratios are characteristically lower across the study area with the highest value of 1.718 in the Congo Craton and lowest value of 1.704 in the Garoua Rift (Table 3, Fig. 4). There is a slight difference (∼0.01) in the mean Vp/Vs ratio at this depth range between the northern (i.e. Adamawa Plateau) and southern CVL. Other regions such as the Oubanguides Belt and the Coastal Plain revealed a mean Vp/Vs ratio of ∼1.7. The mean Vp/Vs ratio ranged from 1.758 in the Oubanguides Belt to 1.72 in the Garoua Rift and Coastal Plain in the lower crust. These values are slightly higher than the values in the middle crust at different tectonic regions but they are still lower than the mean Vp/Vs ratios in the upper crust (Fig. 4). The mean Vp/Vs ratios are slightly higher (∼0.01) in the Adamawa Plateau than the CVL and are about 1.73 in the Congo Craton. The regional crustal Vp/Vs ratio averages indicate a slightly higher value of 1.74 in the Oubanguides Belt, median value of 1.73 in the Congo Craton, Adamawa Plateau, CVL and Garoua Rift, and the lowest value of 1.72 in the Coastal Plain.

Furthermore, we compared the whole-crustal averages to the result of earlier receiver function HK-stacking study in Cameroon (Gallacher & Bastow 2012) and likewise the upper crustal averages to the network average of 1.74 found by De Plaen et al. (2014) during their seismicity study in Cameroon (Fig. 5). We found that the average Vp/Vs ratios in the upper crust is close (±0.01) to the network average of 1.74 found by De Plaen et al. (2014) at most of the seismic stations with the exception of a few (e.g. CM07, CM11, CM14, CM09, CM21, CM29 and CM32). We note that some of these stations are either located close to a fault line (e.g. CM07, CM11 and CM21) in the Oubanguides Belt or at locations where shallow sediment are found (e.g. CM29). Such near-surface structural heterogeneities may have complicated the retrieved receiver functions and affected the accuracy of the previously estimated Vp/Vs ratios by Gallacher & Bastow (2012).

Figure 5.

Comparison of upper-crustal and whole-crustal Vp/Vs averages from this study with bulk Vp/Vs ratios from receiver function Hk-stacking (Gallacher & Bastow 2012) and previously reported Vp/Vs network average from seismic anisotropy studies in Cameroon (De Plaen et al.2014).

A comparison of the crustal-averaged Vp/Vs ratios to the bulk crustal Vp/Vs ratios from receiver function (Gallacher & Bastow 2012) revealed similarities to within ± 0.02 for all the CBSE stations in the Adamawa Plateau region and some stations in other regions (e.g. CM02, CM07, CM10, CM12, CM13, CM30 and CM32, Fig. 5, Table 2). However, large discrepancies (∼0.1) exist for the two stations in the Coastal Plain and other stations (CM11 and CM23) where the Vp/Vs ratios were reported to have been estimated using approximation methods due to the lack of clear reverberation phases on the estimated receiver functions beneath the stations (Gallacher & Bastow 2012). At other stations such as CM04, CM06, CM20, CM23 and CM29, the bulk Vp/Vs ratios from receiver functions are higher (δ_Vp/Vs∼0.1); while they are lower (δ_Vp/Vs∼0.1) at stations CM11, CM03, CM17 and CM16 when compared to the average crustal estimate from this study (Fig. 5). We attribute the large discrepancies to methodological differences and possible error in the observed data. For example, the Vp/Vs ratios from the HK-stacking method may be biased by some inherent assumptions and assumed parameters (e.g. planar horizontal interface, initial P-wave velocity, and phase weighting values) and data noise; while our inverted Vp/Vs ratios may also be affected by our choice of model parametrization and errors in Moho depth estimate.

4 DISCUSSION

4.1 Evaluating the effect of including rayleigh wave ellipticity data

Since this study is the first to compute ZH ratio in Cameroon, we present maps of the mean ZH ratio and associated uncertainties (standard deviation from the mean) for periods in the range 10–62 s beneath the CBSE stations in Figs 6 and 7 and evaluate the effect of including the ZH ratio data in the joint inversion. From Fig. 6, we observe that the heterogeneities in the upper crust are more pronounced at short periods (<30 s). In this period range, high values of ZH ratio coincides with the Congo Craton and the Oubanguides Belt; while relatively low ZH ratio values are observed along the CVL, Adamawa Plateau and the Coastal Plain (Fig. 6). The differences in the values of ZH ratio are indicative of the different seismic properties of the rock types beneath the stations and show a positive correlation with well-known features in study area. The relatively low ZH ratio values observed in the sediment-filled Coastal Plain may be related to the characteristic slow shear wave velocity at shallow depth in the plain (e.g.., Ojo et al.2017). The high ZH ratio values in the Congo Craton and Oubanguides Belt may indicate a weak velocity gradient in the uppermost crust suggesting the presence of high-velocity mafic materials or near-surface outcrops (Lin et al.2014). The stations closer to the Garoua Rift zone filled with shallow sediments (e.g. CM29) have a relatively low ZH ratio value than other nearby stations, indicating a strong positive velocity gradient in the uppermost crust due to the large impedance contrast between the shallow sediments and the underlying crystalline bedrock (Tokam et al.2010). The lowest ZH ratio value is observed in the area around the active Mt. Cameroon volcano (CM09) and may be related to small-scale heterogeneities due to magmatism. At longer periods (Figs 6g-i), the ZH ratio values appear to be much more homogeneous across the different tectonic regions. The uncertainties of the ZH ratio measurements at periods in the range 10 - 62 s are generally smaller than 3 per cent of the estimated values except for stations CM14, CM27 and CM13 (Fig. 7). This is comparable to the estimates from studies in other regions (e.g. Lin et al.2012, 2014; Kang et al.2016; Shen & Ritzwoller 2016). However, we note that the standard deviation from the mean used as a proxy for the measurement error can also be biased based on the number of measurements used to compute it at each period. Hence, it is likely not an ideal estimate of the measurement error but provides useful information about the relative uncertainty between the stations.

Figure 6.

Rayleigh wave ellipticity measurements for the CBSE stations at periods of (a) 10 s, (b) 14 s, (c) 18 s, (d) 22 s, (e) 26 s, (f) 30 s, (g) 42 s, (e) 50 s, and (f) 62 s. The station names are displayed above the color dots. Geological structures are as per Fig. 1.

Figure 7.

Uncertainty estimates for the Rayleigh wave ellipticity measurements at (a) 10 s, (b) 14 s, (c) 18 s, (d) 22 s, (e) 26 s, (f) 30 s, (g) 42 s, (e) 50 s, and (f) 62 s. Geological structures are as per Fig. 1.

Since the use of ellipticity is still considerably less of a standard technique than the joint inversion of receiver function and surface wave dispersion, we compute the depth sensitivity kernels of ZH ratio in order to understand the aspect of the model that is constrained by the observed data and evaluate the exact influence of using this additional constraint by directly comparing two sets of inversions with and without incorporating the ZH ratio data. In Fig. 8, we display the depth sensitivity kernels of ZH ratio to shear wave velocity (dZH/dVs), P-wave velocity (dZH/dVp) and density (dZH/dRho) computed for a reference 1-D model using subroutines modified from the program of Herrmann (2013) (e.g. Yuan et al.2016; Zhang & Yao 2017). It is obvious that the ZH ratio is most sensitive to the shear wave velocity (Vs), a little sensitive to density and no appreciable sensitivity to the Vp, in agreement with other studies (Berbellini et al.2016, 2017; Maupin 2017). We also note that the sensitivity of ZH Ratio to Vs is confined to the shallow crust (≤10 km) even at long periods. However, the sensitivity changes sign from positive to negative with increasing depth and period (Fig. 8b). This indicates that the inclusion of the ZH ratio data in the joint inversion scheme primarily helps constrain the uppermost crust (Chong et al.2015, 2016).

Figure 8.

Depth sensitivity of Rayleigh wave ellipticity (ZH ratio) for a reference 1-D velocity model shown in (a) to shear wave velocity (b), to P-wave velocity (c), and to density (d) at different periods. Blue and red colors indicate positive and negative sensitivities, respectively.

Further examination using two sets of real data in inversions with and without ZH ratio in the study area revealed that the inverted Vs model using three types of data (including the ZH ratio) has lower shear wave velocities (δ_Vs ≤ 0.8 km s^–1) in the top 10 km compared to the model resulting from the inversion using two data types (Fig. 9a). The shear wave velocity difference between the inverted models using two and three types of data at most of the CBSE stations is higher than the mean model error (∼ 0.2 km s^–1) at depths less than 10 km (Fig. 9a). Likewise, the computed model error for the inversion using the ZH ratio data is slightly lower (δ_e ≤ 0.4) than the models resulting from the inversion of only two types of data (Fig. 9b). This is especially true at the uppermost crust where we expect some improvement based on the ZH ratio sensitivity kernel (Fig. 8). Therefore, we think that the inclusion of the ZH ratio data is effective in better constraining the absolute shear wave velocity in the uppermost crust in Cameroon, in agreement with synthetic tests (e.g. Chong et al.2015, 2016). There are also notable differences between the inverted Vp/Vs ratio from the two set of inversions in the top 5 km (Fig. 9c). The inverted model that includes the ZH ratio revealed lower Vp/Vs ratios and smaller model error values (∼0.2) when compared to those resulting from two data types (Figs 9c and d). However, the computed differences in the Vp/Vs ratio seem to be in the same order of magnitude with the differences in the Vp/Vs ratio model error (Figs 9c and d). The computed crustal thickness from both sets of inversions showed similar differences (∼3 km) for most of the CBSE stations (Fig. 9e). However, we obtain smaller error range for the inversion including the ZH ratio (Fig. 9f). The reduced Moho depth errors may be a result of the more improved uppermost crustal structure.

Figure 9.

Differences between inverted crustal parameters with and without ZH ratio data (A) Vs differences, (B) Vs error differences, (C) Vp/Vs ratio differences, (D) Vp/Vs ratio error differences, (E) Moho depth differences, and (F) Moho depth error differences.

4.2 1-D velocity model from joint inversion of receiver function, rayleigh wave dispersion and ZH ratio measurements

4.2.1 Data fit for the joint inversion result

To evaluate the inversion result, we compute the data fit (DF) defined as |$1 - {M_{joint}}$| (eq. (2)) at each station for the different data types (Fig. 10). A value of unity means perfect fit, while lesser values indicate an increasing misfit between the synthetic data and the observed data. Figs 10(a–d) show that the inversion achieved fit with 65 per cent, 93 per cent, 69 per cent and 93 per cent of the stations having DF values greater than 0.95 for group velocity dispersion, receiver function, ZH ratio and phase velocity dispersion, respectively. Overall, the fit to the phase velocity dispersion and receiver function is better than those for the ZH ratio and the group velocity dispersion data. This may be due to errors in the manual picking process adopted for the group velocity dispersion measurements, low-quality data for ZH ratio measurements at some stations (e.g. CM14 and CM31), and small-scale heterogeneities as expected for stations CM09 and CM13 near the active volcano (Fig. 10d). The near-unity DF values obtained in this study (0.87–0.99) demonstrate the robustness of the newly developed inversion scheme and generally indicate that our current model parameterization (Fig. 10; Table 1) is good enough to allow satisfactory data fits for most of the CBSE stations. However, problematic stations like CM09 and CM31 where we could not achieve a good data fit may need further efforts to investigate the appropriate range of model parameters needed to fit all the data types.

Figure 10.

Station-dependent variation in the fit between the synthetic data from the mean model to the observed (a) group velocity dispersion data, (b) receiver function data, (c) phase velocity dispersion data, and (d) Rayleigh wave ellipticity measurements (ZH ratio). The labelled stations are denoted by the color-filled circles representing the value of the misfit according to the legend. The data fit is defined as 1 – M_joint in eq. (2), and unity means perfect fit. Blue lines are major faults (Fig. 1).

4.2.2 Fit improvement and comparison to previous model

To demonstrate the model fitness achieved in this study, we compare the DF values between the synthetic data computed from the model of Tokam et al. (2010) and those from this study using a fit improvement (FI) index defined as

\begin{eqnarray*} &&Fit\ Improvement\ \left( {FI} \right)\nonumber \\ &&\quad = \left[ {\frac{{{{\left( {DF} \right)}_{ThisStudy}} - {{\left( {DF} \right)}_{Tokametal.\left( {2010} \right)}}}}{{{{\left( {DF} \right)}_{Tokametal.\left( {2010} \right)}}}}} \right]\ \times 100 \end{eqnarray*}

(3)

We obtained a Fit Improvement (FI) value ranging from –6 to 37 per cent as shown in Fig. 11. The FI values for the receiver function data ranged from 2.2 to 12.6 per cent, suggesting that our new 1-D Vs model provides better fit to the receiver function data at all stations than the model of Tokam et al. (2010) (Fig. 11b). For the ZH ratio, the FI values ranged from –2.1 to 37.9 per cent (Fig. 11d). Since the models from Tokam et al. (2010) were obtained from joint inversion of receiver function and group dispersion data only, the highest FI values were obtained for the newly added ZH ratio dataset as expected. The FI values for the group and phase velocity dispersion data ranged from -6.1 per cent-18.5 per cent (Figs 11a and c). The negative FI values imply that the model of Tokam et al. (2010) could explain the observations better than our Vs model at some stations and for some data types. However, our new model generally provides better fit to all the types of data at most stations than the model of Tokam et al. (2010) (Fig. 11). In Fig. 12, we show examples of the joint inversion result overlain by the model of Tokam et al. (2010) for four representative stations (CM02, CM12, CM26 and CM29) located in different tectonic settings (Table 2; Fig. 1). In comparison to the model of Tokam et al. (2010), the ZH ratio provides a tighter constraint on the uppermost crustal structure without tradeoff with the deeper structure as shown by the mismatch between the observed and predicted data at short periods (gray line in Fig. 12). Hence, significant differences are observed in the near-surface structure (e.g. CM02, Fig. 12) while deeper structures are similar in most cases (e.g. CM12, CM26, and CM29, Fig. 12). For example, the localized high-velocity layer observed in the model of Tokam et al. (2010) in CM02 (Fig. 12) is absent in our model, and we find this smooth model with gradual increase of velocity with depth geologically more plausible for a cold stable cratonic region. In the Supplementary Material, we plot our joint inversion result overlain by the model of Tokam et al. (2010) for all the CBSE stations.

Figure 11.

Fit improvement (eq. (3)) of the new models over the existing models of Tokam et al. (2010) for (a) group velocity dispersion data, (b) receiver function data, (c) phase velocity dispersion data, and (d) Rayleigh wave ellipticity measurements (ZH ratio) at each station. The texts on top of the color dots are the station names and the colors are related to the computed fit improvement values at that station and for specific data type. Positive percentages imply that the new model have better fit to the observed data while negative ones imply that the model of Tokam et al. (2010) have a better fit to the observed data. The blue lines inside each of the images represent the major faults in the study area as denoted in Fig. 1.

Figure 12.

Joint inversion results for four representative stations located in the Congo Craton (CM02), Oubanguides Belt (CM12), Adamawa Plateau (CM26) and Garoua Rift (CM29) (Fig. 1). In each of the figures, the upper-left box shows the histogram distribution of the inverted Moho depth for the best 300 models with the mean Moho depth and standard deviation stated above it. The bottom-left box shows the distribution of 1-D models and their probability. Overlain on it is the corresponding model of Tokam et al. (2010) for the station. The upper-right box shows two black lines representing the upper and lower bound of the standard deviation of the measured receiver function data. Red lines are synthetic receiver functions for the best 300 models and the thin black line is the synthetic receiver function from the mean model. The data fit for the mean model and Tokam's model is written above it together with the station name. The middle-right box shows the observed ZH ratio data in black dots, the synthetic ZH ratio from Tokam's model in gray dots, the synthetic ZH ratio from the mean model in black line, and the synthetic ZH ratio from the best 300 models in red lines. Above it are the fit values for the mean model and Tokam's model. The bottom-right box shows the observed dispersion data in black dots and the synthetic dispersion curves from Tokam's model plotted in gray dots; the synthetic dispersion curves from the mean model is plotted in black lines and the synthetic dispersion curve from the best 300 models in red lines. Above and below the plotted data are the fit values for the mean model and Tokam's model for phase and group velocity dispersion data.

4.2.3 Differences in crustal architecture across cameroon

Since the inclusion of ZH ratio primarily helps to constrain the uppermost crustal structure (Fig. 8), we therefore examine the features of our new models at crustal depths less than 10 km. In summary, we observe a localized high-velocity layer (HVL) beneath stations CM21 and CM07 which may be related to their proximity to known fault zones in the study area (Fig. 1). In the cratonic region, the shear wave velocity (Vs) is faster than other regions (3.6–3.8 km s^–1) even at shallow depth and increases monotonically with depth into the crust. In locations filled with near-surface sediment like the Coastal Plain, Mamfe Basin and Garoua Rift (e.g. CM01, CM05, CM15, CM18 and CM29), the shear wave velocity (Vs) at the surface is typically very low (<3 km s^–1) and increases at a fast rate with increasing sensitivity to the basement that underlies the sediment. Generally, the Vs beneath the sediments at these stations have values exceeding 3.6–3.8 km s^–1 at shallow depths. In the Oubanguides Belt, stations CM03 and CM14 also reveal indications of localized high-velocity layer with Vs exceeding 3.6 km s^–1 below 10-km depth. However, for all other stations, the shear wave velocity is found to be less than 3.6 km s^–1 at depths above 10 km.

Furthermore, in Figs 13 and 14, we present depth slices of the inverted Vs models and associated uncertainties across the study area. From these maps, it is obvious that the shear wave velocity generally increases across the study area from the upper crust to the lower crust, in agreement with our previous studies (e.g. Ojo et al.2017, 2018b). Except for some discrepancies in the uppermost crustal structure detailed above, our shear wave velocity estimates are in agreement with the results of Tokam et al. (2010). Fig. 13(a) reveals slow shear wave velocity along the Coastal Plain, around the active Mt. Cameroon Volcano and the Mamfe Basin (Fig. 1). The slow velocities in this location are likely due to the presence of shallow sediment in the basin and along the coastline but may also indicate the presence of melt in the shallow crust beneath the active Mt. Cameroon region that erupted recently (e.g. Suh et al.2003). In Fig. 13(g), we observe that the shear wave velocities in the Garuoa Rift and southeastern part of CVL are higher than the threshold for crustal material (∼4.2–4.3 km s^–1), suggesting an increasing sensitivity to the uppermost mantle velocities beneath the Moho at a relatively shallower depth than other locations (e.g. Christensen & Mooney 1995; Christensen 1996; Tokam et al.2010). The lower crust in Cameroon is characterized by high shear wave velocities greater than 4.0 km s^–1, which could result from the process of metamorphism or possible evidence of mafic addition to the lower-crust during the Pan-African orogeny (Figs 13h and i; Tokam et al.2010; Gallacher & Bastow 2012). However, the observed Vp/Vs ratios in the lower crust are atypical for large volumes of cooled mafic intrusion due to mantle plume-crust interaction (e.g. Stuart et al.2006; Heit et al.2015).The estimated uncertainty of our inverted Vs model is generally less than 0.15 km s^–1 (∼3 per cent of the measured value) in most part of the study area except for Station CM14 where limited quality data were available due to instrument failure during the deployment (Fig. 14).

Figure 13.

Shear wave velocity distribution beneath each seismic station at different crustal depths (a) 2 km, (b) 6 km, (c) 10 km, (d) 15 km, (e) 20 km, (f) 25 km, (g) 30 km, (e) 35 km, and (f) 40 km. The blue lines inside each of the images represent the major faults in the study area and all other labels are same as denoted in Fig. 1.

Figure 14.

Uncertainty estimates for the inverted shear wave velocity at (a) 2 km, (b) 6 km, (c) 10 km, (d) 15 km, (e) 20 km, (f) 25 km, (g) 30 km, (e) 35 km, and (f) 40 km. The blue lines inside each of the images represent the major faults in the study area and all other labels are same as denoted in Fig. 1.

4.2.4 Localized High-velocity layer (HVL) in the upper crust

To take further advantage of the new ZH ratio dataset and broadly investigate the existence of HVL in the upper crust, we developed a new 3-D model of the Cameroon crust (Section 2.4.2) and examined several profiles cross-cutting the study area in different directions and tectonic regions. Although the new 3-D model has a good resolution in the entire crust, we focus on the depth range of 0–20 km (Fig. 15) where we have achieved improved resolution and is more relevant to the question of the existence of HVLs being addressed here. The cross sections shown in Figs 15(a) and (b) reveal the presence of localized HVL (3.7–4 km s^–1) at depths of 2–8 km in the Oubanguides Belt and Adamawa Plateau concealed beneath a thin layer (1–2 km) characterized by relatively slow Vs ranging 2.5–3 km s^–1 (cyan and yellow stars in Fig. 1). Below the HVL is an obvious low-velocity zone characterized by shear wave speed of about 3.4–3.5 km s^–1 (Fig. 15a). The high shear wave velocity suggests that the rock type is of mafic composition and we hypothesize that these are images of cooled mafic intrusions in the uppermost crust probably emplaced when the Sao Francisco Craton, the Congo Craton and the West African Craton collided during the formation of Gondwana (Castaing et al.1994; Toteu et al.2004). The HVL can also be observed on the 1-D Vs models obtained from the joint inversion results for stations CM03, CM07, CM21 and CM14 (Supplementary Material). We took several profiles across the study area in order to constrain the locations of the HVL and found that they are localized structure within the Oubanguides Belt and the Adamawa Plateau. We did not find such localized high velocity features along the CVL as shown in Fig. 15(c). We think this is reasonable and in support of the observation of low Vp/Vs ratios across the study area (Fig. 17) as large volumes of mafic materials would have been expected to change the bulk crustal composition from being predominantly felsic to intermediate compositions. However, the shear wave velocity in the Congo Craton is faster than those in other tectonic regions as expected and the fast velocity is consistent across the entire Craton and not a localized structure as we observed in the Oubanguides Belt and the Adamawa Plateau (Fig. 15d). The extreme right portion of Fig. 15(a) reveals the fine structure of the sediment filled in the Garoua Rift with shear wave velocity of 2.7–2.8 km s^–1 and thickness of about 2 km. The shear wave velocities of 3.6–3.8 km s^–1 can be observed in the upper crust in Fig. 15 as indicated earlier by Tokam et al. (2010).

Figure 15.

Cross-sections of the 3-D shear wave velocity model from joint inversion of ZH ratio, group and phase velocity dispersion data using the Neigbourhood Algorithm. (a) Profile AA’, (b) Profile BB’, (c) Profile CC’, and (d) Profile DD’ denoted in Fig. 1 by the cyan colored dash lines. Plotted on top of the cross-sections are the topographic profiles depicted in red. Approximate locations of geologic features of interest are labelled as CC: Congo Craton; OB: Oubanguides Belt; AP: Adamawa Plateau; GR: Garoua Rift; MC: Mt. Cameroon; CVL: Cameroon Volcanic Line; and CP: Coastal Plain (Fig. 1).

4.3 Comparison of crustal thickness estimates and variation with topography

The Moho depth estimate in this study came from the summation of the five inverted crustal thicknesses as detailed in the model parameterization in Table 1, and the topographic data were obtained from the ETOPO1 global relief model (Anante & Eakins 2009). Fig. 16 presents a plot of the surface topography and crustal thickness estimates from this study and others for comparison purposes. The inverted Moho depth and associated uncertainties from our study is plotted with error bars while estimates from the traditional HK-stacking method (Zhu & Kanamori 2000) and two other studies (Tokam et al.2010; Gallacher & Bastow 2012) are plotted on it. The inverted Moho depth ranges from 28 km in the Garoua Rift to 45 km in the Congo Craton and the estimated uncertainties are in the range of 2.2–7 km. Along the CVL, the Coastal Plain and the Oubanguides Belt, the crustal thickness ranges from 33–38 km (Fig. 16). For most stations, the Moho depth from other studies fell within the range of uncertainties of our crustal thickness estimates. However, discrepancies exist for the estimates of Tokam et al. (2010) for stations CM01, CM05 and CM09 along the coast and near the active Mt. Cameroon volcano. The new Moho estimate also seems to be consistently higher than other estimates for stations CM29, CM30, CM31 and CM32 in the Garoua Rift and slightly lower for some stations in the Congo Craton (CM02 and CM03). We attribute the minor discrepancies to different methodology and data quality (Fig. 16).

Figure 16.

Surface topography and Moho depth estimates from this study and previous studies beneath the Cameroon Broad-band Seismic Experiment (CBSE) stations (Fig. 1).

We investigated the relationship between crustal thickness and topography in Cameroon and found that the variations are clearly different for the northern and the southern part of the study area (Fig. 16). While the Moho depth is positively correlated with the topography in Adamawa Plateau and Garoua Rift, there are no obvious trends in southern CVL, the Congo Craton and the Oubanguides Belt. This implies that the variation of crustal thickness does not have a first-order influence on the topography in southern Cameroon; while a similar Airy isostasy behavior as we observed in the Adamawa Plateau and the Garoua Rift has also been seen in the eastern margin of the Tibetan Plateau and other regions of the world (Clark & Royden 2000; Clark et al.2005; Royden et al.2008).The obvious anti-correlation of the Bouguer anomaly and topography in the Adamawa Plateau and Garoua Rift is consistent with Airy isostasy and corroborates our interpretation (Figs 16 and 17). Based on the results of experimental crustal modeling, the relationship between the Moho undulation and topography observed in the Adamawa Plateau may indicate structures formed during a compressional regime while the opposite relationship observed in the Garoua Rift may depict structures formed during an extensional regime (Martinod & Davy 1992; Eyike et al.2018). A comparison of the Moho estimate from gravity data (Figs 17a and b; Eyike et al.2018) and those from this study (Fig. 17c) reveals a maximum difference of about 14 km in the Garoua rift and southwestern part of the CVL (Fig. 17d). Since the Moho estimates from the gravity study (Fig. 17b) is computed relative to the geoid and not to the topography, the anomalies in the Garoua Rift and southwestern CVL could provide additional evidence for crustal thinning in these locations, in agreement with previous studies (e.g. Tokam et al.2010; Ojo et al.2017).

Figure 17.

Comparison of crustal thickness estimates from gravity and seismic data: (a) Bouguer gravity anomaly, (b) Moho depth from gravity (Eyike et al.2018), (c) Moho depth from joint inversion results and (d) absolute difference between (b) and (c).

4.4 Possible explanations for the observed Vp/Vs ratio in cameroon

The crustal Vp/Vs ratio (or the Poisson ratio) is important for constraining the bulk crustal composition and has been shown to offer better constraints than individual Vp or Vs parameters because of its reduced sensitivity to variations in pressure and temperature (Christensen 1996). Vp/Vs ratios are often related to rock types and mineral contents through comparisons with published laboratory data. For common rock types, the Poisson ratio (σ), which is directly related to the Vp/Vs ratio (κ), varies from 0.20 (κ = 1.63) to 0.35 (κ = 2.1) and is particularly sensitive to composition. An increase in the silica content lowers σ, whereas the higher abundance of mafic components will increase the value (Christensen 1996). For lower-crustal rocks, low (<0.26), intermediate (0.26–0.28) and high (>0.28) σ values are characteristics of felsic, intermediate, and mafic compositions, respectively (Zandt & Ammon 1995). To improve the estimate of this geophysical parameter and allow for multidimensional imaging that permits better interpretation, we directly invert for layered Vp/Vs ratio and associated uncertainties in our joint inversion scheme (Figs 18 and 19). The inverted Vp/Vs ratio shows variation in the range of 1.67–1.85 (Poisson ratio from 0.17 to 0.24) with uncertainties generally less than 0.05 (<3 per cent of measured values) in the entire crust (Figs 18 and 19). As detailed in Section 3.7, the distribution of average crustal Vp/Vs ratios in all the six tectonic regions in the upper crust suggests an intermediate composition, while a predominantly felsic composition characterized the middle crust, and a low to intermediate composition in the lower crust (Fig. 4, Tables 2 and 3). This is surprising as it is generally expected that the upper crustal rocks will be more felsic than the lower crustal rocks (e.g. Christensen 1996, Fig. 4, Table 3). This anomalous compositional trend may suggest lower-crustal partial delamination as a result of eclogitization or small-scale modification of the original mafic lower crustal rocks beneath the study area (e.g. Thompson et al.2010; Gilligan et al.2016). We also observe a positive correlation between the average crustal Vp/Vs ratios and the estimated crustal thicknesses described by the linear relationship Vp/Vs ratio (κ) = 0.0018*crustal thickness (H) + 1.6651 in Cameroon. Although the coefficient of determination (R² = 0.2972) does not explain all of the observed variations, it generally indicates a decreasing abundance of felsic composition in the crust with increasing crustal thickness. The slightly elevated regional average Vp/Vs ratio in the Oubanguides Belt may be related to the presence of mafic intrusions in the upper crust as observed in this study (Figs 4 and 15). The average regional Vp/Vs ratio value of 1.73 in the Congo Craton is comparable to values usually encountered in cratonic regions due to the Tonalite–Trondhjemite–Granodiorite (TTG) geology of typical Archean domain (e.g. Thompson et al.2010). However, the range of Vp/Vs ratios is comparable to other tectonic regions in the study area. Alternatively, the observation of a thick crust with an overall felsic composition as we observed in the Congo Craton is typical of a palaeo-collisional zone, which may suggest that the Cratonic crustal root has deflected the rising magmas from the upper mantle sideways without any significant modification over time to form the CVL, in line with the model proposed by King & Anderson (1995, 1998).

Figure 18.

Inverted Vp/Vs ratio beneath the CBSE stations at crustal depths of (a) 2 km, (b) 6 km, (c) 10 km, (d) 15 km, (e) 20 km, (f) 25 km, (g) 30 km, (e) 35 km and (f) 40 km. The blue lines inside each of the images represent the major faults in the study area and all other labels are same as denoted in Fig. 1.

Figure 19.

Uncertainty estimates for inverted Vp/Vs ratio at (a) 2 km, (b) 6 km, (c) 10 km, (d) 15 km, (e) 20 km, (f) 25 km, (g) 30 km, (e) 35 km and (f) 40 km. The blue lines inside each of the images represent the major faults in the study area and all other labels are same as denoted in Fig. 1.

The surprisingly low range of Vp/Vs ratios along the CVL with mean value of 1.73 is less than the global average value of 1.765 and atypical of hotspot volcanic regions (e.g. Christensen 1996; Stuart et al.2006; Gallacher & Bastow 2012). Similar observations have been made by Gallacher & Bastow (2012) who reported a markedly low network average Vp/Vs ratio of 1.74 and ruled out the presence of large amounts of melt, mafic intrusions or underplating beneath Cameroon. There are several possibilities that can explain the observed low Vp/Vs ratios in Cameroon. The first and obvious one is the presence of large volumes of quartz-rich rocks (e.g. sandstone, granodiorite and granite-gneiss) or increasing contents of silica in the crust, which in turn makes the bulk crustal composition to be more felsic when compared to the average continental crust (Zandt & Ammon 1995; Christensen 1996). Secondly, since the composition of the lower crust plays an important role in the determination of the bulk Vp/Vs ratio, it is possible that there is a significant compositional modification of the mid-lower crust or a removal of the mafic part of the crust by processes such as basal erosion or lithospheric delamination such that we now have an overall felsic crust (Griffin & O'Reilly 1987; Mengel & Kern 1992; Elsheikh et al.2014). Although the range of observed Vp/Vs ratios is broadly in favor of delamination of the lower crust, the relatively fast shear wave velocity observed in the mid-lower crust (Fig. 13) does not lend support to a large-scale removal of the mafic component of the mid-lower crust in the study area (Table 3, Figs 4, 5 and 13). However, the thickness of layer with Vs ≥ 4.0 km s^–1 in the lower crust is significantly thin beneath the Adamawa Plateau and relatively thin along the CVL, the Oubanguides Belt and the Garoua Rift (Table 3). Therefore, we cannot rule out the possibility of incomplete delamination with partial preservation of the original lower crust.

Additionally, low Vp/Vs ratio values could arise as a result of the misinterpretation of the high-velocity lower crust for the uppermost mantle, such that the Moho is misplaced to the top of the lower crust like the case of Greenland (e.g. Dahl-Jensen et al.2003; Artemieva & Thybo 2008; Thybo & Artemieva & Thybo 2013). Although we observe high velocity in the mid-lower crust in Cameroon, which typically results from metamorphism (e.g. eclogite facies) or mafic additions to the base of the crust, we think this is not likely the case for Cameroon based on the similarity of Moho depth estimates from different methods and datasets (Figs 16 and 17). Based on the available evidence, we think that the first case mentioned above offers a more plausible explanation for the observed low Vp/Vs ratio in Cameroon, although other possibilities also remain. The existence of large-scale formations with low densities as observed on the Bouguer gravity anomaly map (Fig. 17a) further corroborates our interpretation. However, a strong discrimination against the different possibilities for the entire study area or a specific region may require further studies that will be the focus of our future research.

5 CONCLUSIONS

In this study, we developed a joint inversion scheme based on the NA (Sambridge 1999a,b) to simultaneously invert receiver function, Rayleigh wave dispersion and ellipticity data for improved 1-D shear wave velocity, crustal thickness and Vp/Vs ratio models in Cameroon, West Africa. We obtained results that provide better constraint primarily to the uppermost crustal structure and offer misfit improvements over previous models. The new Rayleigh wave ellipticity (ZH ratio) data included in this study showed a good correlation with the local geology. High ZH ratio coincides with the Congo Craton and the Oubanguides Belt while relatively low ZH ratio is found along the CVL, Adamawa Plateau and the Coastal Plain, with uncertainties generally smaller than 3 per cent of the measured ZH ratio value. We obtain good data fit from the joint inversions with 65, 93, 69 and 93 per cent of the stations having fit values greater than 0.95 for group velocity dispersion, receiver function, ZH ratio and phase velocity data, respectively. Since the improvement in our model is in the uppermost crust, we focus on the features in these depth ranges and investigate the reported fast velocity layer in Cameroon. The Vs beneath the Congo Craton is generally faster than other regions ranging between 3.6 and 3.8 km s^–1 at shallow depth. Likewise, stations CM01, CM05, CM15, CM18 and CM29 with significant near-surface sediment exhibit a very low Vs at shallow depth that increases rapidly above 3.6 km s^–1 with increasing sensitivity to the underlying basement. For stations along the CVL, the shear wave velocity is consistently slower than other regions with a value of about 3.6 km s^–1 at shallow depths. At greater depths, our shear wave velocity estimates are in agreement with the results of Tokam et al. (2010). The uncertainty of our Vs models is generally less than 0.2 km s^–1 (i.e. ∼4 per cent of the measured value) in most part of the study area. The crustal thickness estimates ranged from 28 km in the Garoua Rift to 45 km in the Congo Craton with uncertainties in the range of 2.2–7 km, in agreement with previous studies (Tokam et al.2010; Gallacher & Bastow 2012). Similar to observations in other regions of the world, the Moho depth is directly correlated to the topography within the Adamawa Plateau and the Garoua Rift in Northern Cameroon, but no discernible relationship exists for stations in the southern part of the study area. In agreement with Gallacher & Bastow (2012), the inverted Vp/Vs ratio ranged between 1.67 and 1.85 (Poisson ratio from 0.17 to 0.24) with uncertainties generally less than 0.05 (<3 per cent of measured values) at different depths within the crust, reflecting the heterogeneity of crustal composition. Assuming isotropy, the Vp/Vs ratio is atypical of hotspot volcanic regions; while the network mean Vp/Vs ratio of 1.73 is less than the global average (∼1.765), thereby suggesting a felsic to intermediate crustal bulk composition in Cameroon (Zandt & Ammon 1995; Christensen 1996). We attribute this to large volumes of quartz rich rocks and/or increasing contents of silica in the crust and possible lower-crustal delamination. However, we note that there are several other possible scenarios that may be responsible for the unexpectedly low Vp/Vs ratio, which demands further research efforts. Based on the results obtained in this study, we conclude that the newly estimated crustal parameters (1-D Vs, Vp/Vs ratio and Moho depth) are well constrained by the sensitivities of the different data types, especially in the uppermost crust, and therefore place a tighter constraint on the structure and composition of the crust in Cameroon than previous studies.

SUPPORTING INFORMATION

Figure S1. Joint inversion results for all the CBSE stations (31 in total). In each of the 31 plots, the upper-left box shows the histogram distribution of the inverted Moho depth for the best 300 models with the mean Moho depth and standard deviation written above it. The bottom-left box shows the distribution of 1-D models and their respective probability. Overlain on it is the corresponding model of Tokam et al. (2010) in gray thick line for each station. The upper-right box shows two black lines representing the upper and lower bounds of the standard deviation of the measured receiver function data. Red lines are synthetic receiver functions for the best 300 models and the thin black line is the synthetic receiver function from the mean model. The fitness values for the mean model and Tokam's model are written above it together with the station name. The middle-right box shows the observed ZH ratio data in black dots, the synthetic ZH ratio from Tokam's model in gray dots, the synthetic ZH ratio from the mean model in black line, and the synthetic ZH ratio from the best 300 models in red lines. Above it are the fitness values for the mean model and Tokam's model. The bottom-right box shows the observed dispersion data in black dots, the synthetic dispersion curve from Tokam's model is plotted in gray dots; the synthetic dispersion curve from the mean model is plotted in black line; and the synthetic dispersion curve from the best 300 models is plotted in red lines. Above and below it are the fitness values from the mean model and Tokam's model for phase and group velocity dispersion data.

Figure S2. 2-D Maps showing group velocity dispersion (a, d, g), phase velocity dispersion (b, e, h) and ZH ratio (c, f, i) in Cameroon at periods of 10, 20 and 30 s.

Table S1. Number of data points from receiver function, ZH Ratio, phase velocity and group velocity used in the joint inversion.

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

ACKNOWLEDGEMENTS

We are grateful for the detailed and constructive comments received from two anonymous reviewers and the Editor, Professor Michael Ritzwoller, which greatly improved the quality of the submitted manuscript. The seismic dataset used in this study was retrieved from the IRIS Data Management Center (http://www.iris.edu) using the Standing Order for Data (SOD) program (http://seis.sc.edu/SOD). We appreciate the time and efforts of the participants of the Cameroon Broadband Seismic Experiment in deploying and maintaining the temporary network. This research was supported by the Major Program of the National Natural Science Foundation of China (41590854). We thank Malcolm Sambridge of the Australia National University, Australia for making his Neighbourhood Algorithm code openly available. Most of the figures in this paper are plotted using the Generic Mapping Tools (http://gmt.soest.hawaii.edu). We thank Huajian Yao of the University of Science and Technology of China for insightful discussions. We are grateful to Dr Eyike of the University of Douala, Douala, Cameroon for providing us his recent Moho estimates and acknowledge the Cameroon gravity dataset made available to him by the authorities of the Bureau Gravimetrique International (BGI).