Reappraisal of effective elastic thickness in the south-west Indian Ocean , and its possible implications

The bathymetry of the south-west Indian Ocean is dominated by three mid-oceanic ridge systems: the Chagos–Laccadives Ridge, the Central Indian Ridge, and the Mascarene Plateau. Although there have been a number of geophysical and geological investigations over the region, the genesis of these morphological features is still contradictory. Most of the estimations of effective elastic thickness in this region have been carried out in the spectral domain, either by transfer function analysis or by freeair admittance analysis. As these investigations were along some one-dimensional profiles or discrete blocks, spatial variation of the effective elastic thickness was not achieved. Here, we reappraise the estimation of effective elastic thickness in the south-west Indian Ocean by performing the computation in the spatial domain using flexure inversion. During this process, we also estimate the Moho depth throughout the region by two independent processes: gravity inversion, and flexural inversion. The Te values (effective elastic thickness) are estimated in the spatial domain, which match well with the results in the spectral domain obtained with the free-air admittance method. In addition, there is spatial variation of the Te values over the area analyzed. Our estimated Te values are low (1-6 km) along the Chagos-Laccadives ridge, implying its proximity to a spreading ridge at the time of creation. The Te values along the Mascarene Plateau show spatial variation with a seafloor age from north (Te, ca. 4 km) to south (Te, ca. 20 km). These findings substantiate earlier data and suggest that Réunion was created due to intraplate volcanism.


Introduction
The south-west part of the Indian Ocean is an amalgam of complex geological features that include, for example, oceanic ridges, seamounts, rifted margins, and hotspots.Among these different tectonic elements, three active midoceanic ridge systems dominate the topography: the Chagos-Laccadives Ridge (CLR), the Central Indian Ridge (CIR), and the Mascarene Plateau (Figure 1).The CLR is one of the most prominent aseismic ridges of the south-west Indian Ocean, and it lies within the north-east sector.The 2,500-km-long CLR extends from 15˚N to 12˚S along the western continental margin of India.The CIR is a divergent plate boundary that separates the Indian plate from the African plate.The CIR is a part of the global Mid-Oceanic Ridge system that extends southwards from the Carlsberg Ridge at 2˚N, to the Rodriguez Triple Junction at 25˚S, in the Indian Ocean.The third of these prominent aseismic ridges, the Mascarene Plateau, is one of the most prominent shallow submarine plateaus in the south-west Indian Ocean.The Mascarene Plateau lies to the north and north-east of Madagascar, starting from the Seychelles in the north, and extending through the Saya de Malha Bank and Nazareth Bank down to Réunion in the south.
Although several geophysical and geological investigations have been carried out to understand the crustal and the tectonic evolution of these morphological features, to date there is no unanimous theory regarding their genesis.For example, the CLR was described as a micro-continent by Krishnan [1968], whereas according to Narain et al. [1968] it is a transition zone between the oceanic and continental crusts.Fischer et al. [1971] and Sclater and Fischer [1974] proposed that the CLR consists of volcanic edifices that originated from leaky transform faults; this was later supported by Avraham and Bunce [1977].On the contrary, according to Naini [1980] and Naini and Talwani [1982], the CLR is part of a stretched continental crust that is somewhat similar to the Laxmi Ridge in the Arabian Sea.Chaubey et al. [1995] suggested that the CLR is formed of volcanic edifices emplaced over an oceanic crust.In contrast, Chaubey et al. [1998] and Chaubey et al. [2002] suggested that the CLR is part of a continental crust that is mingled with volcanic edifices.Several other studies and analyses have favored the possibility of a hotspot trace concept for the origin of these ridges [Francis and Shor 1966, Dietz and Holden 1970, Whitmarsh 1974, Duncan 1981, Morgan 1981, Tiwari et al. 2007].
The Mascarene Plateau is considered to be a small continental outlier that is similar to Madagascar, and is separated from the continent.The Seychelles Islands are considered anomalous as they are granitic, as opposed to the basaltic nature of most of the oceanic rock.In the northern part of Réunion, the Saya de Malha bank is considered either to be an extension of the Seychelles micro-continent [McKenzie and Sclater 1971] or a volcanic platform formed as a consequence of Deccan volcanism [Courtillot et al. 1986].Several contradictory hypotheses exist to explain the genesis of the Mascarene Plateau: e.g., according to Kamen-Kaye and Meyerhoff [1980], the Mascarene Ridge is a relict island, whereas the hotspot theory is equally popular.An equally important suggestion is that it was created due to active volcanism in a transform zone.
As the effective elastic thickness of an oceanic plate increases with age, and as observations at long distances from active plate boundaries give values that approximately correspond to the depth of the 450˚C isotherm [Watts andTorné 1992, Wessel 1992], the spatial variation of the elastic thickness mapped on a regional scale can explore the interactions between the tectonic processes and the longterm lithospheric strength.However, most of the earlier estimations of elastic thickness in this region used either transfer function analysis along some one-dimensional profiles, or free-air admittance over some discrete blocks [Ashalata et al. 1991, Tiwari et al. 2007].Hence, these have failed to capture the spatial variation of the elastic thickness.
In the present study, we adopt a method that necessarily operates in the spatial domain [Braitenberg et al. 2006].Analysis in the spectral domain has certain limitations; e.g., the estimated values have poor stability in the case of a small topographic load.Moreover the two-dimensional spectral analysis requires a square map of the topography and gravity data, which is not always available, for geographical reasons.In addition, the resolution of the Te map is much higher in the case of the spatial domain.
In the present study, we aim to reappraise the variations in the effective elastic thickness at the mid-oceanic ridges and the adjoining oceanic regions in the south-west Indian Ocean.These display considerable variations in their sediment thickness, crustal thickness, and gravity signatures.

The data
Our analysis requires bathymetry data, Bouguer gravity anomaly data, sediment thickness data, seafloor age data, and crustal thickness data.The bathymetry, Bouguer anomaly, and sediment thickness data are shown in Figure 2. We used the bathymetry grid extracted from the Digital Atlas of the General Bathymetric Chart of the Oceans, published by the British Oceanographic Data Centre on behalf of the International Oceanographic Commission of UNESCO and the International Hydrographic Organization.This grid has a resolution of 30 × 30 arc.s.The free-air anomaly data come from the Global Marine Gravity Field from ERS-1 and the GEOSAT geodetic mission altimetry of Andersen and Knudsen [1998] and Andersen et al. [2008].Both the bathymetry and the free-air gravity data have been corrected for the effects of the cooling lithosphere, using the Parsons and Sclater [1977] subsidence formula.The seafloor age data come from Müller et al. [1997].We restricted our analysis to well within the oceans, to avoid intricacies involved in the generation of satellite measurements from continental regions.We have converted the marine free-air gravity anomaly data ( ) to the Bouguer gravity anomaly ( ) using the slab formula: (1) where = 1,670 kg.m -3 is the density contrast between the surface rock and water, H is the bathymetry, and G is the gravitational constant.The sediment thickness model is taken from the Total sediment thickness of the world's oceans and marginal seas (World Data Center for Marine Geology and Geophysics, Boulder).The digital total sediment thickness database was compiled by the National Geophysical Data Center of the National Oceanic and Atmospheric Administration, with a resolution of 5 arc.min×5arc.min.These data estimate the total sediment thicknesses from the seafloor to the top of the basement.The crustal thickness data used here were obtained from the CRUST 2 model [Bassin et al. 2000].

Methodology
The analysis consists of a sequence of steps that include gravity forward and inverse modeling, and the inverse modeling of flexural rigidity [Braitenberg et al. 2006].In more detail, the steps of the methodology applied are shown Figure 3.

Gravity forward and inverse modeling
The gravity forward calculation is performed for a single boundary that defines a density discontinuity at the Moho.The gravity effect of an undulating discontinuity that divides two bodies with a constant density contrast is calculated using the Parker algorithm [Parker 1972], which makes use of a series expansion up to order 5 of the gravity field generated by an oscillating boundary.As the gravity potential field is proportional to 1/r, r being the distance from the observer to the mass point, the short wavelength variations of the gravity field are filtered out with increasing distance from the mass source [Zadro 1986, Blakely 1995].Therefore the Moho undulations contribute to the long wavelength part of the observed gravity field, whereas the short wavelength part is generated by the superficial masses due to sediment layers or intracrustal density inhomogenities.Hence, to determine the Moho undulations, the gravity values are corrected for the bathymetry (the Bouguer correction) and the sediment thickness variations.The forward modeling of the gravity field produced by a sedimentary basin for which the sediment thickness is known can be carried out by allowing different density variations with depth, as either a sediment compaction model [Sclater and Christie 1980] or a model of linear sediment density increase with depth.The final corrected gravity field is then inverted by applying an iterative algorithm that alternates downward continuation with direct forward modeling [Braitenberg and Zadro 1999].Thus this method is somewhat similar to the Oldenburg-Parker inversion approach [Oldenburg 1974].This analysis requires some initial parameters to be set (as shown in Table 1); e.g. the reference depth of the density interface, the density contrast across the interface, and the cut-off wavelength.

Flexure modeling: the Te inversion:
The next step is the isostatic modeling that follows the isostatic lithosphere flexure model [Watts 2001].This is an independent process to determine the Moho undulations, and it allows us to compare these with the Moho undulations obtained from gravity inversion.This inverse flexural analysis is superior to the spectral methods, in the sense that it has a relatively high spatial resolution.The flexure is calculated by the convolution approach with the numerical flexural response functions [Braitenberg et al. 2002[Braitenberg et al. , 2003].In the flexure inversion process, the crustal load is required, which is obtained from the density model and the Moho undulations and other initial parameters shown in Table 1.The crustal flexural rigidity is then inverted to match the known load for the given crustal thickness model.This results in a model of the Moho undulation based on the theory of flexural isostasy.

Parameter Symbol Value
Mean crustal density Table 1.Parameters used in the analysis.
We discuss briefly how the uncertainties in the input data and those coming from modeling assumptions affect the estimations of the elastic thickness.The main objective is the recovery of the bathymetry, which should be as exact as measured by direct observations.The procedure used here is to null the gravity signal in the presence of sediment cover.To achieve this goal, we used the whole gravity signal, and therefore integrated the gravity inversion with the flexural inversion.Thus, we followed a procedure as shown in Figure 3, and similar to that used by Ebbing et al. [2007] in the South China Sea.The inversion procedure integrates the available crustal structure, the crustal thickness, the sediment thickness model, and flexural models of isostasy.The accuracy of the derived Te estimates/flexural Moho depends on the gravity anomaly data.The gravity data that is analyzed in the study region have an accuracy that ranges from 2 mGal to 3 mGal [Anderson and Knudsen 1998].Thus the long wavelength part is dependent on the long wavelength part of the gravity field, the sediment model and the gravity effects of the crust mantle interface.In the present case, the sediment model is seen to have little effect on the derived Te inversions.However, neglecting the sediment thickness can result in underestimation of the basement depth.The density contrast captures the undulations of the crust mantle.In certain cases, the real density contrast used is either larger or smaller with respect to the actual observations, which results in over estimation or under estimation of the crust mantle interface.

Results
The sediment thickness (1,000 m) in most of the region is not very significant (Figure 2d).However, to obtain improved accuracy, we made the necessary corrections for the Bouguer gravity field considering the sediment density as ca.2,250 kg/m 3 .The corrected Bouguer anomaly shows a variation of around 50 mGal to 100 mGal in most of the ridge portions, whereas in the oceanic part, its value is high and positive, reaching almost ca.300 mGal to 400 mGal.This corrected gravity field is used as input in the constrained gravity inversion procedure.The cut-off wavelength is taken as 40 km, to filter out the shortwavelength anomalies.We repeated the inversion process of the Bouguer gravity anomaly by varying the reference depth and the density contrast within standard values of 20 km to 30 km and 3.5 kg/m 3 to 6.5 kg/m 3 .The final result was obtained for a reference depth of ca.23 km and a density contrast of ca.0.55.The crustal thickness values obtained (Figure 4a) are in good agreement with the values of the CRUST2 model.Using the gravity-inversion-derived Moho as input and also the equivalent topography as input, we proceed with the flexure analysis (Te inversion).We have inverted the Te values within a range of 1 km to 20 km in windows of size 110 km 2 ×110 km 2 , with a spacing of 60 km between the windows.The flexure-derived Moho shown in Figure 4b is in excellent agreement with the gravity-derived Moho.The low values of the residual Moho (Figure 4c; -2.5 to +2.5 km) supports this.For positive values, the flexural Moho is shallower than the gravity Moho, and vice versa.We have used standard values for the elastic parameters (Table 1): e.g., for the Young modulus (100 GPa), the Poisson ratio (0.25), and the mantle and crustal densities (3,350 kg/m 3 and 2,800 kg/m 3 , respectively).The variation in the elastic thickness shows low values along the CLR (1-6 km) compared to the Mascarene Pateau (ca.4-20 km) (Figure 4d).Along the Mascarene Plateau we obtained a variation of the Te values from south (Réunion, 20 km) to north (Seychelles, 4 km).

Discussion
The elastic thickness values obtained in our analysis show considerable variations along the CLR and the Mascarene Plateau (Figure 4d).The Te values along the CLR are low (1-6 km) whereas along the Mascarene Plateau, the Te values are high (ca. 4-20 km).There is also spatial variation in the Te values along the Mascarene Plateau, where they decrease from 20 km in the south, to 8 km in the north.The seafloor spreading ages are shown in Figure 2a.The Te values along the Mascarene Plateau vary with age.The Te value is ca.12 km at the Seychelles, where the age is ca.65 Ma.As we progress southwards along the Mascarene Plateau, the age decreases, while the Te values increase; e.g., the Te values are ca.4 km to 9 km at Saya de Malha Bank (63 Ma), 14 km to 18 km at Nazareth Bank (34 Ma), and 20 km at Réunion (0 Ma) These data are somewhat similar to the Te values computed in the spectral domain (using free-air admittance analysis) by Tiwari et al. [2007].A comparative study of the elastic thickness by these two methods is shown in Table 2.However, one clear advantage of our analysis is that the spatial variation of Te is available over the whole region, unlike with the data of Tiwari et al. [2007], where only one Te value is computed for some discrete blocks.
The low Te values along the CLR imply its proximity to a spreading ridge, whereas the high Te values at Réunion strongly suggest an intraplate volcanism [Watts 2001].The low Te values along the CLR at the western continental margin of India suggest that the lithosphere was weak and therefore younger at the time of emplacement of the volcanic material.In other words, the CLR was created by the Réunion hotspot when it was close to the spreading ridge (CIR) axis.These low values are quite similar to the results obtained by Ashalata et al. [1991] using transfer function analysis of gravity and bathymetry data along seven onedimensional profiles.Chaubey et al. [2008] also reported similar low results (Te ca.2-3 km) from admittance analysis along the Laccadive ridge.This hypothesis also finds support in the strong spreading symmetry along the Carlsberg Ridge.On the contrary, the spatial variation of the Te values with age along the Mascarene Plateau implies intraplate volcanism.At around 65 Ma, as a consequence of the Deccan flood basalt event, the Indian plate started migrating northwards, which caused the Central Indian ridge to pass over the Réunion hotspot in the Indian Ocean.With the northward progress of the Indian plate, the Réunion hotspot influenced the adjacent areas, initiating the arrangement of the CLR and the reorganization of the nearby spreading centers.Thus, to summarize, the estimation of Te in the spatial domain has produced results comparable to the results in the spectral domain.This provides strong support in favor of the hotspot theory in the genesis of the lithosphere of the southwest Indian Ocean.Moreover, a database of the Moho depth for the region has been created by gravity as well as flexural inversion.

Conclusions
By the combined analysis of gravity and flexure inversion, we have obtained the Moho depth for the southwest Indian Ocean.The low values of residual Moho indicate that the Moho depths obtained by the two approaches are in excellent agreement with each other.By using the flexure inversion process, we have obtained a high resolution spatially varying Te map for the whole region.The estimated Te values are in close agreement with the values obtained in the spectral domain.The Te distribution shows low values along the CLR (1-6 km), implying its proximity to a spreading center at the time of its creation.On the other hand, the Te values along the Mascarene Plateau show variation with seafloor age, implying intraplate volcanism.The Te along the CIR is relatively high, but there are thermal effects and lateral crustal density changes that must be taken into account, and which influence the outcome of the Te values.However, this study would go beyond the scope of the present report and would require some more constraining data.

Figure 1 .
Figure 1.Bathymetry map of the study area, with the major tectonic features shown.The black rectangle indicates the area over which the calculations were performed.The locations of Réunion, Chagos Bank and Maldives (short name of Chagos-Maldives-Laccadives Ridge), and Saya de Malha Bank, Nazareth Bank and the Seychelles are shown.
Figure 2. The maps of the study region.(a) Bathymetry map.(b) Free-air gravity map.(c) Bouguer gravity map.(d) Sediment thickness map.

Figure 4 .
Figure 4. Moho depth and elastic thickness.(a) Moho depth from the constrained gravity inversion.(b) Moho depth from the flexural inversion.(c) Residual Moho (gravity-inversion-derived Moho -flexural-inversion-derived Moho).(d) Spatial variation of the effective elastic thickness.The white squares in (d) show the blocks used byTiwari et al. [2007] to compute the Te by the free-air admittance method.

Table 2 .
Comparisons between the estimated Te values in the spectral and spatial domains.