# Water, oceanic fracture zones and the lubrication of subducting plate boundaries-insights from seismicity

David Schlaphorst, J. -Michael Kendall, Jenny S. Collier, James P. Verdon, Jon Blundy, Brian Baptie, Joan L. Latchman, Frederic Massin, Marie-Paule Bouin

<figure><img src="image_1.png" /></figure>

To cite this version:

David Schlaphorst, J. -Michael Kendall, Jenny S. Collier, James P. Verdon, Jon Blundy, et al.. Water, oceanic fracture zones and the lubrication of subducting plate boundaries-insights from seismicity. Geophysical Journal International, 2016, 204, pp.1405-1420. 10.1093/gji/ggv509 . insu-03581674

HAL Id: insu-03581674

https://insu.hal.science/insu-03581674

Submitted on 20 Feb 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

# Water, oceanic fracture zones and the lubrication of subducting plate boundaries—insights from seismicity

David Schlaphorst,<fnref n="1" /> J-Michael Kendall,<fnref n="1" /> Jenny S. Collier,<fnref n="2" /> James P. Verdon,<fnref n="1" /> Jon Blundy,<fnref n="1" /> Brian Baptie,<fnref n="3" /> Joan L. Latchman,<fnref n="4" /> Frederic Massin<fnref n="5" /> and Marie-Paule Bouin<fnref n="5" />

Accepted 2015 November 23. Received 2015 November 22; in original form 2015 March 10

## SUMMARY

We investigate the relationship between subduction processes and related seismicity for the Lesser Antilles Arc using the Gutenberg–Richter law. This power law describes the earthquake-magnitude distribution, with the gradient of the cumulative magnitude distribution being commonly known as the *b*-value. The Lesser Antilles Arc was chosen because of its along-strike variability in sediment subduction and the transition from subduction to strike-slip movement towards its northern and southern ends. The data are derived from the seismicity catalogues from the Seismic Research Centre of The University of the West Indies and the Observatoires Volcanologiques et Sismologiques of the Institut de Physique du Globe de Paris and consist of subcrustal events primarily from the slab interface. The *b*-value is found using a Kolmogorov–Smirnov test for a maximum-likelihood straight line-fitting routine. We investigate spatial variations in *b*-values using a grid-search with circular cells as well as an along-arc projection. Tests with different algorithms and the two independent earthquake catalogues provide confidence in the robustness of our results. We observe a strong spatial variability of the *b*-value that cannot be explained by the uncertainties. Rather than obtaining a simple north–south *b*-value distribution suggestive of the dominant control on earthquake triggering being water released from the sedimentary cover on the incoming American Plates, or a *b*-value distribution that correlates with on the obliquity of subduction, we obtain a series of discrete, high *b*-value ‘bull’s-eyes’ along strike. These bull’s-eyes, which indicate stress release through a higher fraction of small earthquakes, coincide with the locations of known incoming oceanic fracture zones on the American Plates. We interpret the results in terms of water being delivered to the Lesser Antilles subduction zone in the vicinity of fracture zones providing lubrication and thus changing the character of the related seismicity. Our results suggest serpentinization around mid-ocean ridge transform faults, which go on to become fracture zones on the incoming plate, plays a significant role in the delivery of water into the mantle at subduction zones.

**Key words:** Spatial analysis; Earthquake dynamics; Seismicity and tectonics; Subduction zone processes; North America.

## 1 INTRODUCTION

Slab lubrication plays an important role in the mechanical behaviour of subduction zones. Potential lubricants include water contained within the pore spaces of ocean floor sediments or structurally bound in hydrous minerals in hydrated crust and uppermost mantle por-

tions of the slab (Reynard 2013). Once these materials enter the subduction zone, their water is released by a combination of porosity reduction and mineral dehydration reactions in response to changing pressure and temperature along the subduction geotherm. Water release is likely to be continuous along the subduction path, albeit with spikes related to the breakdown of specific mineral phases (Schmidt

<footnote n="1">¹School of Earth Sciences, University of Bristol, Bristol, United Kingdom. E-mail: david.schlaphorst@bristol.ac.uk</footnote>

<footnote n="2">²Department of Earth Science & Engineering, Imperial College London, London, United Kingdom</footnote>

<footnote n="3">³British Geological Survey, Edinburgh, United Kingdom</footnote>

<footnote n="4">⁴Seismic Research Centre, The University of the West Indies, St. Augustine, Trinidad & Tobago</footnote>

<footnote n="5">⁵Institut de Physique du Globe de Paris, Paris, France</footnote>
Poli 2014). The estimated depth range for appreciable water release is between 50 and 150 km (Anderson et al. 1978; Kawakatsu & Watada 2007). Water released from the slab can ascend into the overlying peridotite mantle wedge. Close to the trench, where the wedge is relatively cool, water ingress produces serpentine (Bostock et al. 2002). Conversely, water that enters the hotter part of wedge at some distance from the trench will lower the peridotite solidus temperature and trigger copious melting to produce primitive hydrous arc basalts (Pichavant et al. 2002; Skora & Blundy 2010; Grove et al. 2012). This mantle melting in response to water flux in turn leads to the generation of volcanic arcs (Tatsumi 1989; Schmidt & Poli 1998; Hattori & Guillot 2003; Rüpke et al. 2004; Hacker 2008; Grove et al. 2012) typically located 65–130 km above the slab (England et al. 2004).

Serpentine (Mg₃Si₂O₅(OH)₄) is a hydrous mineral that can host considerable structurally bound amounts of water in its structure (up to 13 wt per cent) and has an upper thermal stability of ~720 °C (Ulmer & Trommsdorf 1995). Consequently, serpentine breakdown within the subducting slab has been widely implicated in water supply to the mantle wedge beneath volcanic arcs (Scambelluri et al. 2004). There are two potential sources of the serpentine on the incoming plate: the first is that formed at mid-ocean ridges by hydration of mantle peridotites exposed within core complexes (Schmidt & Poli 1998; Smith et al. 2006; Kawakatsu & Watada 2007). Serpentinized mantle peridotite is particularly abundant on seafloor formed close to transform faults at slow-spreading ridges such as the Mid Atlantic Ridge (Alt et al. 2013). The second potential source of the serpentine forms along faults that are generated as the plate bends as it approaches the trench (Ranero et al. 2003). In the first case the occurrence of serpentine would be expected to correlate with large-offset fracture zones on the downgoing plate, in the second case it would be expected to be evenly distributed along the length of the arc. In both cases subduction is able to carry serpentine to considerable depths below the arc, in principle until the slab temperature exceeds serpentine thermal stability.

The presence of fluids will also affect seismicity associated with subduction zones by raising the pore pressure leading to an increase in the number of small earthquakes (Wiemer & Benoit 1996). In most tectonic settings the earthquake-magnitude distribution follows the Gutenberg–Richter power law, $log\{N(M)\} = a-bM$, where N is the number of events larger than magnitude M and a and b are constants (Gutenberg & Richter 1954). The gradient of the slope, called the *b*-value, gives a measure of the distribution of earthquake magnitudes. A higher *b*-value implies a higher abundance of smaller events.

Several studies have sought to map *b*-values using various methods. Studies range from regional studies of a number of subduction zones (Wiemer & Wyss 1997, 2000; Katsumata 2006; Amorèse et al. 2010; Vorobieva et al. 2013; Nishikawa & Ide 2014; Tormann et al. 2015) to more local studies of specific volcanic areas (Wyss et al. 1997) and during the aftershocks of big events (De Gori et al. 2014). Depending on the comprehensiveness of the seismicity catalogue, variation in *b*-value can be investigated laterally, vertically and/or temporally. Microseismic data from hydraulic fracturing (Williams & Le Calvet 2013) and from laboratory experiments (Scholz 1968) has also been used to better understand factors that affect magnitude distributions.

In general, the *b*-value lies in a range of 0.5–2.0, with the average value for global catalogues in the range of 0.79–1.25 (Frohlich & Davis 1993). It has been previously suggested that the absolute *b*-value varies with fault type: on average, thrust faulting producing a value as low as 0.7, normal faulting a value around 1.1 and strike-slip

movement a value around 0.9 (Schorlemmer et al. 2005). However, caution has been raised about comparing absolute *b*-values (e.g. Frohlich & Davis 1993; Wiemer & Benoit 1996), as they can change when using different calculation methods or different seismicity catalogues.

Over a broad range of scales, *b*-values have an inverse relationship to the differential stress ($b \sim 1/\Delta\sigma$). This has been tested using laboratory experiments (Scholz 1968) and fractal simulations (Huang & Turcotte 1988), as well as seismicity induced by mining (Urbancic et al. 1992), volcanic activity (Wyss et al. 1997) or hydraulic fracturing (Wyss 1973). Using global seismicity catalogues, Schorlemmer et al. (2005) concluded that this relationship holds true on both microscopic and macroscopic scales. Therefore, knowledge of the *b*-value can give direct information about the dynamic process in an area of subduction. The dehydration of the downgoing slab is likely the main trigger for high *b*-values, as it increases the pore pressure and lowers the effective stress (Wiemer & Benoit 1996; Wiemer & Wyss 1997), which is linked to the differential stress. In effect, the presence of water facilitates slip at lower differential stress, leading to an increased fraction of small earthquakes. Conversely, the absence of water inhibits slip except at higher differential stresses, leading to relatively few small events and more large magnitude events and so a smaller *b*-value.

In addition to differential stress, the *b*-value can also be influenced by the heterogeneity of the material (Mogi 1962a,b). For example, based on an experimental study, Sammonds et al. (1992) find a negative correlation between the *b*-value and the combined effects of the stress and crack length.

In this work, we map *b*-values along the Lesser Antilles Arc (LAA, see Fig. 1) in an effort to correlate any observed variation with structural features. We identify a number of possible influences, the first group relating to the broad-scale tectonics and the second water-delivery. Within the former group, the first factor is the curvature of the arc, which increases the amount of strike-slip movement towards the north and south of the arc. Drawing on Schorlemmer et al. (2005), if this feature influences the seismicity it would lead to gradual change to higher *b*-values towards the north and south. The second factor is variation in slab dip, which steepens from an average value of 50° in the central arc to nearly vertical close to Venezuela in the south (Bouysse & Westercamp 1990). If this feature influences the seismicity we would expect to see a gradual change in *b*-value towards the south. The third factor we identify within the broad-tectonic group is the possible presence of a slab gap at depth in the centre of the arc (van Benthem et al. 2013). Again, this feature might be expected to directly influence the *b*-values at this part of the system. Within the water-delivery group of possible influences, the first factor is the nature and abundance of subducted sediment. This has a bimodal distribution along the length of the arc (Plank & Langmuir 1998) from sediment-rich subduction in the south to sediment-poor subduction in the north. If sediment thickness is a key factor in the supply of water to the subduction zone then one would expect a gradual change along the arc with higher *b*-values in the south. Conversely, if the supply of water to the arc is controlled by the structurally bound H₂O in serpentine, the second factor in this group, then the variation along the arc would be irregular, correlating with the presence or absence of ocean-floor features associated with serpentinization, such as fracture zones. Thus, the two prime sources of water that enter the subduction zone have distinctive spatial patterns, potentially giving us the opportunity to discriminate between them by the observed seismicity pattern.
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Map of the Lesser Antilles Arc. (a) Bathymetry map (Gebco 2008) of the study area. The strong bathymetric contrast between the southern and northern parts of the forearc is due to the sediments delivered to the Atlantic seafloor from the South American continent being trapped by prominent basement ridges such as the Tiburon Ridge on the incoming plate and the development of a thick accretionary prism near Barbados (Ba) in the south (Bouysse & Westercamp 1990). Also shown is the location of the trench (toothed line) marking the subduction of the North and South American plates beneath the Caribbean plate together with the 50-, 100-, 150 and 200-km slab-depth contours (Gudmundsson & Sambridge 1998). These slab contours are independent of the locations of the local earthquake catalogues used in this study. The islands are (from north to south): S, Saba; StE, St Eustatius; StK, St Kitts; Mo, Montserrat; Gu, Guadeloupe; D, Dominica; Mq, Martinique; StL, St Lucia; StV, St Vincent; Gr, Grenada. The triangles mark the location of the three-component broadband stations of the SRC (yellow) and the IPGP (blue) networks used in this study. (b) Location of the study area in a broader tectonic setting. The position of map (a) is indicated by the red outline. (c) Vertical gravity gradient of the latest satellite altimetry-derived grid (Sandwell et al. 2014) for the Central Atlantic from the Mid-Atlantic Ridge (black line) in the east to the Lesser Antilles in the west. To map the most significant fracture zones (black lines, named) we followed the methods of Matthews et al. (2011) for the satellite data combined with an inspection of shipboard magnetic data. An uninterpreted image of a larger area is given in Fig. S1. Projections of these fracture zones into the subduction zone are shown by dashed black lines. We use the geometries of the fracture zones still preserved on the African Plate to predict the trends of those on the subducting American plates. In performing the projections we assume that fracture zones are vertical and therefore not deflected as the plate bends.</figcaption></figure>

We test these five hypotheses, concluding that variations in *b*-values along the LAA are predominantly controlled by fracture zones and the associated water released from serpentine emplaced at the Mid-Atlantic Ridge. We conclude by discussing the implications this has for arc volcanism.

## 2 TECTONIC SETTING

The LAA extends from the eastern end of the Greater Antilles in the north to the South American continent and results from subduction

of the North and South American Plates beneath the Caribbean Plate (Wadge & Shepherd 1984; Bouysse & Westercamp 1990; DeMets et al. 1990). The boundary between the two American plates forms a ~200 km wide zone of mild north-south compressional structures such as the Baracuda Ridge and the Tiburon Rise (Patriat et al. 2011). Recent seismic tomography has suggested the presence of a slab gap at this boundary with depth (van Benthem et al. 2013). The subduction rate is relatively slow, at 18–20 mm yr⁻¹ (DeMets et al. 2000), in an almost arc-perpendicular direction (~95°). The slow subduction rate could be the reason the seismicity is relatively low
for an active margin (Bouysse & Westercamp 1990). The dip of the subduction zone is believed to be roughly constant at around 50° along the entire arc (Bengouhou-Valerius *et al.* 2008), but it may be slightly steeper north of Martinique (50–60°) compared to the south (45–50°), but becomes almost vertical close to Venezuela in the extreme south (Bouysse & Westercamp 1990). Due to the curvature of the arc the strike-slip component of subduction increases to the north and the south.

The volcanic arc itself is constructed on thickened oceanic crust of the Caribbean Plate (Mauffret & Leroy 1997). Beneath the arc, estimates of crustal thickness vary between 24 and 35 km (Boynton *et al.* 1979; Christeson *et al.* 2008; Kopp *et al.* 2011). The arc splits into two branches north of Guadeloupe. The inner arc, made up of St Kitts, St Eustatius and Saba, is still active today. Remnants of the inactive outer arc include the islands of Antigua and St Martin (Bouysse *et al.* 1990). The westwards shift of the northern arc occurred at approximately 9 Ma (Briden *et al.* 1979). To the west the LAA is bordered by the Grenada Basin and the Aves Ridge, a remnant arc of Mesozoic age (Kearey 1974; Neill *et al.* 2011).

Fed by the Orinoco River, the incoming plate in the south is rich in continental clastic sediment. Towards the north where the sediment supply becomes blocked by several submarine banks, including the Tiburon Ridge, the incoming plate is dominated by pelagic marine sediments. Traditionally this spatial variation in the amount and nature of sediment on the incoming plate has been assumed to play a significant role in LAA subduction processes, including the systematic north to south variations in arc magma geochemistry (White & Dupre 1986; Carpentier *et al.* 2008). However, other authors have suggested that significant volumes of the incoming sediments in the south are scraped off to form the Barbados accretionary prism (e.g. Faugères *et al.* 1993) and therefore do not enter the subduction zone.

There are five major fracture zones meeting the arc between 12°N and 18°N, Fifteen-Twenty, Marathon, Mercurius, Vema and Doldrums. Seafloor mapping at the modern-day Mid-Atlantic Ridge axis shows that first-order ridge segments between the transform faults that later become these fracture zones can be roughly divided into thirds, where only the middle part does not show core complexes (Smith *et al.* 2006). Within the core-complex zones the exposure of ultramafic rocks at the seabed can be as high as 65–90 per cent (Alt *et al.* 2013 and references therein) and their serpentinization can result in the acquisition of up 13 per cent H₂O (Boschi *et al.* 2013). Therefore the bands of oceanic lithosphere around the fracture zones represent significantly hydrated zones.

## 3 DATA

We use two seismicity catalogues to estimate *b*-values (Fig. 2). The first is provided by the *Seismic Research Centre* (SRC) of The University of the West Indies (The UWI) from the years 2008 to 2013 (hereafter referred to as SRC2008). It uses the Duration Magnitude (Md) for all events with values rounded into 0.1 bins.

The second earthquake catalogue is from the two *Observatoires Volcanologiques et Sismologiques (de la Guadeloupe: OVSG; and de la Martinique: OVSM)* of the Institut de Physique du Globe de Paris (IPGP) from the years 1996 to 2012 (hereafter referred to as OVSM1996). In order to make a direct temporal comparison with SRC2008, we also analyse a subset of this data between 2008 and 2012 (OVSM2008), results for which are given in the supplementary figures. This catalogue does not use magnitude bins for the events, but reports magnitudes for each event to the second decimal place.

The catalogue has a lower detection threshold than the SRC2008, but only covers the region north of approximately 13.5°N. In both catalogues we exclude crustal events in our analysis (events located above a depth of 40 km), thus focusing on slab events (Figs 3 and 4).

To optimize coverage of the entire study area, we will use both catalogues. This needs to be handled with care as previous workers have shown that different absolute *b*-values can be catalogue-dependent, and might represent changes due only to different magnitude values from different stations or magnitude estimations (Frohlich & Davis 1993). In our case, the seismicity maps, histograms and average *b*-values suggest a clear difference between the two catalogues used in this study. However, while we recognize that each catalogue may produce different absolute *b*-values, it is the relative variations across the arc that we are most interested in. The relative *b*-value change within the same catalogue is unaffected by the features of the other catalogue.

## 4 METHODS

We test a number of different statistical approaches to *b*-value calculations and their associated uncertainties. These include both the maximum-likelihood method (Aki 1965; Utsu 1965; Bender 1983) and an explicit line-fitting of the data in log space together with either a Kolmogorov–Smirnov (KS) test or a Lynden-Bell algorithm (Lynden-Bell 1971). This creates four different method combinations that will be discussed in more detail in Section 4.1.

In Section 4.2 we detail the criteria used for dividing the catalogue into subsets based on the location and time of the events. For clarity spatial groups of events will be called ‘cells’, temporal groups will be referred to as ‘time bins’ and groups based on event magnitudes as ‘magnitude bins’. Cells can be either circular, using a grid-search over latitude and longitude, or sections along a pre-defined arc-line (seen in Fig. 3a). Furthermore, the ‘size’ of a cell can be defined by its spatial dimensions (radius of cells or length along the arc-line) or by the number of events per cell.

In Section 4.3, we discuss some of the challenges and issues associated with calculating and interpreting *b*-values.

### 4.1 Calculating *b*-values

In this study we test different methods to calculate *b*-values. First we use a maximum-likelihood method to fit a straight line to the data in M-log(N) space. To establish whether the quality of this best fit straight line is acceptable we then apply a KS-test. This test compares the cumulate probability density functions (PDF) of the data and the constructed straight line fit. The test uses a significance level to determine the maximum accepted value for the vertical difference between these PDFs. Typical values for this level range from 1 to 20 per cent, with higher significance levels making the test more rigorous. For this reason we use a high level of 20 per cent. In our analysis we calculate uncertainty boundaries as the point at which the KS-test between the data and straight line fits with gradually altered gradients are significantly different. Here, a lower significance level increases the uncertainty. We use a level of 5 per cent. Choosing both a high significance level for the main calculation and a low significance level for its uncertainties generates a stricter evaluation of significant *b*-value variation.

The Magnitude of Completeness ($M_C$; also called ‘minimum magnitude’ or $M_{min}$) is chosen automatically to be the smallest magnitude at which the similarity between the data curve and the
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Overview of the full seismicity catalogue data used in this work regardless of the event depth. (a) Seismicity maps of the Lesser Antilles Arc for the OVS catalogue from 1996 to 2012 (left-hand column) and the SRC catalogue containing data from 2008 to 2012 (right column). The black lines (A–A' and B–B') indicate cross-sections shown in Fig. 4. Note that the SRC2008 data are collated in 0.1 magnitude intervals. The arc has a region of low seismicity around the island of St Vincent (~12.5°N–13.5°N), visible in the SRC2008 catalogue. Note that OVS1996 is incomplete south of 13.5°N. (b) Histograms of the magnitude distributions for the different catalogues. Note the large number of events without magnitude information in the SRC catalogue, which have been excluded from our analysis. (c) Average *b*-values for the two catalogues. The green and blue circles and lines represent the cumulative catalogue. The blue circles show all the events used to calculate the *b*-value (solid red line) and its uncertainties (dashed red lines). Therefore, the minimum magnitude is shown by the colour change of the circles from green to blue.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Overview of the events used to calculate the *b*-values. The presentation is similar to Fig. 2. In (a), the maps show the arc-line, built as a smooth curve through all of the arc islands. The line is divided into 100-km-intervals that are used as the *y*-axis of the along-arc calculations shown in Figs 8 and 9.</figcaption></figure>

straight line fit is sufficient, based on a manually preset significance level, so that there is no need to restrict the seismicity catalogue to events above minimum magnitude. The maximum-likelihood method was developed with the assumption that there is no magnitude uncertainty and no magnitude upper bound. However, it has been tailored using different constraints. If data are grouped into magnitude bins, the results are biased, but this bias is negligible

for an interval as small as 0.1 (Bender 1983), which is the case for SRC2008.

We also try alternative methods, using the KS-test, to compare *b*-value variation along the arc (see supplementary figures); rather than using the maximum-likelihood method it is possible to calculate the *b*-value from an explicit fitting of the data in log space (see Figs S7 and S12).
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Cross-sections perpendicular to the arc showing event locations in OVS1996 (left-hand column) and SRC2008 (right-hand column). The cross-sections follow the lines A–A' (top) and B–B' (bottom) shown in Fig. 2. Further cross-sections can be found in Bengouhou-Valerius *et al.* (2008). The influence of station distribution (Fig. 1) of the two catalogues is clearly seen, with the clustering of stations in the central part of the arc for OVS array causing the paucity of events in the southern B–B' profile compared to the SRC catalogue. Our analysis of *b*-values includes all events deeper than 40 km (dashed line), the majority of which clearly follow the line of the subducting slab, but with a subset of events within the overlying wedge. Since we focus on the variation in latitude and longitude, our *b*-value results do not discriminate between these various depth regions.</figcaption></figure>

Additionally to the KS-test, we use a method introduced by Williams & Le Calvez (2013): they replace the commonly used Magnitude of Completeness by a conservative detection limit using the Lynden-Bell algorithm (LB). All events with magnitudes greater than the minimum are assumed to be included in the catalogue. This approach generates a new minimum number of events required per cell. Like the KS-test approach, we use the LB-algorithm with either the maximum-likelihood method or the explicit line-fitting in log space (see Figs S8 and S9 and S12).

Using the KS-test as well as the LB-algorithm together with a maximum-likelihood approach and an explicit line-fitting of the data we are able to test for the stability of observed *b*-value variation over a range of methods.

## 4.2 Setup of spatial cells

We apply different approaches to group the events, described in the following paragraphs as Case 1 to Case 4. We applied these various approaches to ensure that any *b*-value pattern observed is not an artefact of the calculation process.

Case 1 and Case 2: we use a grid-search to map the variation in *b*-value. For this each 0.1 degree in latitude and longitude forms the centre-point of an individual cell.

Case 1: for every point on the grid we use circular cells with a fixed radius and select all events within each cell (cf. Wiemer & Wyss 1997). The *b*-value is mapped at the centre of the circle and cells with less events than a given threshold are rejected. The cell size was chosen to be 100 km for OVS1996 and OVS2008 and 200 km for SRC2008 due to this catalogue's sparser event density. Smaller cell sizes lead to a higher number of cells that have fewer than the minimum number of events required for the analysis. Larger cell sizes tend to smooth out small-scale patterns; we note that previous studies have shown that significant variations

can already be observable at scales of 10 km (Wiemer & Benoit 1996). The results were therefore also compared to different cell sizes to check for stable patterns. To get a higher resolution in areas with a high density of events, each event was used as a cell centre in addition to the pre-defined gridpoints. We use a minimum number of 50 events (cf. Wiemer & Wyss 1997). Schorlemmer et al. (2005) use at least 100 events per cell, but are more confident about results with at least 200 events. Bender (1983) uses 25 events, which leads to a standard deviation of at least $\sigma b \sim 0.25b$; this value decreases to $\sigma b \sim 0.1b$ for simulated data with $N = 100$.

Case 2: the cells are circular but we fix the number of events per cell to 100. The distance of the furthest event therefore varies according to the density of events in that area (cf. Wiemer & Benoit 1996; Wyss et al. 1997). The maximum allowed search distance was set to 100 km for OVS1996 and OVS2008 and 300 km for SRC2008. Note that as the grid-search method is set up using *x–y* coordinates, it does not explicitly consider depth (*z*) when calculating distances. However in our case since most of the events follow the downgoing slab the method serves as a basic tool to map the variation in *b*-values with depth.

Case 3: we use a 2-D analysis of the *b*-value with the *x*-axis being the arc contour ('along-arc'). The arc contour is shown in Fig. 3. For this purpose the individual events have been projected from their original location onto the arc using the direction of motion of the subducting plate. The width of the cell is kept constant. For this technique the distance of the event to the arc also becomes an adjustable parameter. We include all events with a maximum distance of 200 km to the arc-line, but different distances are chosen to compare the results.

Case 4: we use a similar approach to Case 3, but in this analysis the number of events per cell is kept constant.

A graphical summary of all methods applied can be seen in Fig. 5 and Fig. S3.
**(a) Legend:**

Colour: Catalogue

OVS1996

OVS2008

SRC2008

Texture: Cell Geometry

Grid-Search

<figure><img src="image_8.png" /></figure>

Along-Arc

Shape: Calculation method

Boundary Style: Cell Size

KS-test, max. likelihood

<figure><img src="image_9.png" /></figure>

Constant cell size

<figure><img src="image_10.png" /></figure>

KS-test, expl. line fitting

<figure><img src="image_11.png" /></figure>

Constant no. events

LB-alg., max. likelihood

Number: Case

LB-alg., expl. line fitting

e.g., Case 1
(grid-search with
constant cell size)

<figure><img src="image_12.png" /></figure>

**(b) Figures**

<figure><img src="image_13.png" /><figcaption>**Figure 5.** Graphical summary of the techniques used to estimate *b*-values and a key to how they are applied in each figure. The legend is shown in (a). The shape (square or circle) is used to indicate the calculation method (see text for detail). A combination of colour, texture and boundary style is used to indicate, respectively, the seismicity catalogue, the cell geometry and the cell size. The approach used in each figure, is indicated in (b); for example, Fig. 6(a) uses the SRC2008 catalogue, the cells are segments along the arc and the cell size is constant (referred to as Case 3 in the text). Note that a KS-test using the maximum likelihood method is applied to calculate the individual *b*-values. See Fig. S3 for a similar graphical summary, which includes results for the entire range of calculation methods.</figcaption></figure>

## 4.3 Challenges working with *b*-values

**Incomplete earthquake data sets:** Artefacts can be introduced by incomplete earthquake data sets, so that different minimum magnitudes have an effect on the results (Frohlich & Davis 1993; Vorobieva et al. 2013). Following large main shocks the seismicity catalogues can be incomplete due to missing aftershocks (Enescu et al. 2009). Rather than using an overall $M_{min}$ we use the full dataset,

calculating individual $M_{min}$ for each *b*-value estimation with the significance level of the KS-test. Our tests show that a truncation at $M_{min} = 2.5$ does not change the result.

**Temporal variation:** Previous studies have shown that the *b*-value can vary over time (see Smith 1981 for a detailed discussion). These temporal changes in the magnitude range could bias the results (Frohlich & Davis 1993; Wiemer & Wyss 2000). Palacios et al. (2006) observe multiple slope gradients from temporal displacements of the probability function, which relates the event to its magnitude and which is truncated at the magnitude of $M_{min}$. In our study, due to the relatively short time interval of SRC2008 and the lack of any unusual periods of seismicity, we do not expect large temporal variation. For OVS1996 we check the temporal stability by dividing the catalogue into time bins. To match the time interval of SRC2008 we additionally use the OVS catalogue for the time interval between 2008 and 2012 ('OVS2008'), results for which are given in the supplementary material. For areas in California, Wiemer & Wyss (1997) conclude that the temporal changes of the *b*-value are less significant than spatial ones even though they may exist. Furthermore, the work of Urbancic et al. (1992) favors the use of mapping *b*-value changes in a spatial rather than a temporal environment. We can confirm these observations with our results obtained using different time bins of OVS1996. The changes are smaller than the uncertainties of our results.

**Earthquake pattern:** Frohlich & Davis (1993) address the question of whether the *b*-value of fore- and aftershock events differs from the *b*-value of main shocks. Tormann et al. (2015) observe a short (3 months after the event) temporal change in *b*-value around the 2011 Tohoku earthquake. However, impacts of large earthquakes are found to significantly subside towards a background level of spatio-temporal fluctuations. A correlation on a temporal basis has been found in some areas (e.g. Smith 1981; Latchman et al. 2008) but equally there are cases where none could be found (Shi & Bolt 1982). We did not find a correlation between our results and the location of earthquakes with a magnitude greater than 5.

**Spatial variation:** Spatial grouping of the data using fixed radii can lead to problems. If a cell includes different clusters that have different average *b*-values, then different values become superimposed on each other. This effect causes multiple slope gradients in different magnitude intervals in the magnitude-frequency domain (Wiemer & Wyss 2000; Vorobieva et al. 2013; Williams & Le Calvez 2013). In these cases, the power law decay will not normally follow a straight

<figure><img src="image_14.png" /><figcaption>**Figure 6.** Relation of *b*-values to the number of events for SRC2008 using the KS-test with a maximum-likelihood method (see supplementary figures for other calculation approaches and catalogues). (a) Along-arc *b*-value versus number of events. (b) Grid-search *b*-value versus number of events. The uncertainties in *b*-values have been omitted for clarity.</figcaption></figure>
line in the M-log(N) plot; it will instead form a convex curve or have two different gradients at different magnitude intervals. Therefore, we examine individual measurements to observe potential multiple gradients and compare the behaviour of individual results where the drop-off does not depict a straight line. Furthermore, we used different cell sizes and different numbers of earthquakes per cell to estimate the stability of the pattern observed. In our study, the individual measurements do not show multiple slope gradients, except for regions around 17°N. In all of these cases the line is fitted to the gradient around the lower magnitude interval. The second gradient always produces a higher b-value. As a result of this, the uncertainties become very high so that the b-value lies within uncertainties of the lower b-values in the surrounding cells.

**Network coverage:** A change in network coverage can affect the recorded magnitudes and the overall detection of events (Frohlich & Davis 1993). This problem is often analysed on a temporal basis with an attempt to include historical data. In this work we concentrate on recent data, which we check for temporal stability. Still, a similar shortcoming exists on a regional basis; areas with low station coverage can give a biased seismicity map due to unlocated smaller events. Therefore, we have less confidence in the results at the boundaries of our search area. The majority of our results, however, is unaffected by this, since we have at least one seismic broadband station on every major island along the arc.

**Different catalogues:** The use of different catalogues has an influence on the results (Frohlich & Davis 1993). Generally, regions without good station coverage have a high magnitude detection threshold (Wiemer & Wyss 2000). Since we do not compare the absolute b-values this problem does not affect our study.

**Different methods:** Results of b-value studies can be ambiguous if none, or only a single, statistical test method is used (Amorèse et al. 2010). Using different methods described before enables us to check for the stability of our results.

**Minimum number of events:** The use of different fitting techniques can lead to very different b-value estimates, especially if the number of events for the calculation is small (Bender 1983; Frohlich & Davis 1993). For more than *N* = 100 events and magnitudes spanning more than two units, the maximum-likelihood and the weighted least squares methods give nearly the same result (Wiemer & Benoit 1996), whereas the cumulative least squares method tends to give estimates of the b-value that are too high with increasing *N* (Bender 1983). In general, previous studies have shown that the maximum-likelihood approach yields a more constant mean and a more symmetric distribution of all b-values around this mean even if the number of large events changes (Shi & Bolt 1982; Bender 1983; Wiemer & Wyss 1997). In our work, we set a minimum number of events that is large enough to assure for stable results.

**Uncertainties on b-values:** The uncertainties have to be taken into account when b-value patterns are observed; small changes in b-value are often explained by their uncertainties rather than having an actual physical meaning (Frohlich & Davis 1993; Amorèse et al. 2010). With this is mind, we set a significance level of 5 per cent to ensure only large variations to be significant. Using this technique the observed b-value changes are unlikely to be caused by measuring uncertainties alone.

In summary, there are many factors that can have an effect on the reliability of b-values. With our method, however, we are able to exclude these concerns from our study, thus gaining confidence in the significance of our results.

# 5 RESULTS

Significant changes in b-value from both catalogues can be seen in the grid-search and along-arc approaches. Confidence in the robustness of the results stems from the fact that the number of events per cell does not correlate with the b-value (Fig. 6). This becomes even more apparent where variations in b-value are observed when the number of events per cell is fixed. Altering the cell size or width, the prescribed number of events or the maximum distance of event locations to the arc-line does not change the overall pattern. Even though the two catalogues have different mean b-values and different magnitude distributions (see Figs 2 and 3), a similar pattern of variation in b-values can be observed along the arc. Regions of high b-values are found regardless of the method used.

The results of the grid-search calculation (Fig. 7) and the along-arc method (Figs 8 and 9) show that the b-value varies between values of 0.6 and 2.0 for SRC2008 and 0.6 to 1.6 for OVS1996. This range is too large to be explained by uncertainties, which nowhere exceed 0.5. Furthermore, independent of the method or the catalogue used (see supplementary figures), a consistent pattern is observed, with a series of discrete, high b-value ‘bull’s-eyes’ along strike. The along-strike wavelength of the bull’s-eye pattern is around 150–200 km, and there is no overall north-south gradient in absolute b-value. The average event depth around these bull’s eyes is 100 km.

In total the study has revealed four b-value bull’s eye highs. The strongest area is found around 15°N (around 600 km in along-arc distance), approximately between the islands of Martinique and Dominica. Three smaller regions with higher b-values can be observed at 17°N (800 km), around the island of Antigua; at 13.8°N (350 km), just south of St Lucia and at 12.3°N (200 km), just north of Grenada.

Even with a conservative approach the main bull’s eye shows a significant and stable b-value elevation with both catalogues, varying cell sizes or number of events and different calculation methods (Figs 7–9 and supplementary figures). The centre of this anomaly is slightly further to the north for OVS1996 than SRC2008, such that it lies between Dominica and Guadeloupe. Using a fixed number of events, which reduces the cell size in areas of higher event density and therefore increases the resolution, this bull’s-eye structure gets distorted having a more serrated appearance (compare Figs 7d to c for example). Still, this structure keeps a bigger spatial extent than the other three.

The northernmost bull’s eye is especially prominent in the along-arc observation of the OVS1996 catalogue. Due to the large number of events in that catalogue it is possible to either decrease the fixed cell size or to increase the fixed number of events per cell. This allows for a better detectability of smaller structures and for much lower uncertainties (see also Figs S4 and S5). Since this catalogue has a much higher number of events in that area, this raises confidence in the significance of this structure.

The two southern bull’s-eyes of elevated b-values are most prominent in the grid-search calculation using a fixed cell size for SRC2008 (Fig. 7a), although the southernmost one is also particularly prominent using a fixed number of events as well (Fig. 7b). In the along-arc calculation with fixed cell widths only one bull’s eye can be seen (Fig. 8b at around 13°N or 300 km). Since its location is near the middle between the locations of the two bull’s eyes seen with the grid-search it is possible that the difference is the effect of lower resolution due to a smaller number of events in that area.

Unfortunately, the seismicity gap between the islands of Grenada and St Lucia hinders a more detailed mapping of that area. Furthermore, it does not allow for an along-arc calculation with a fixed
<figure><img src="image_15.png" /></figure>

**Figure 7.** Results for the grid-search observation of SRC2008 (a, b) and OVS1996 (c, d) using the KS-test with the maximum-likelihood approach, which gives the most stable result even in areas of sparser event density. (a, c) Calculations using a fixed cell size. (b, d) Using a fixed number of events per cell. Example earthquake-magnitude plots for cells with high and low *b*-values are shown in Fig. S2. Uncertainties in *b*-values for these maps are shown in the Figs S4 and S5. The grey lines show the positions of major fracture zones shown in Fig. 1. Note the change in scale between the catalogues. Also, note the prominent bull's-eye pattern of elevated *b*-values, most notably in the vicinity of Martinique, but to a lesser extent near Montserrat and Grenada. These are best seen in the SRC catalogue (a, b). In the OVS catalogue the pattern is less pronounced since the variation in *b*-values there is smaller.

number of events, since the lower limit of events per cell is exceeded. Also, the two southern structures are invisible to the OVS1996 catalogue due to the insufficient event coverage south of 13.5°N.

## 6 DISCUSSION

The variation in seismicity across the LAA allows us to test between the various tectonic and water-delivery hypotheses outlined in the introduction. Based on our results, we rule out *b*-values being associated with the change in subduction to a more strike-slip dominated regime in the south and north of the curved arc. Additionally, the slab gap does not seem to have a major effect on the behaviour of the *b*-value variation. Furthermore, the north-south contrast in terms of earthquake frequency and slab dip (Wadge & Shepherd 1984; Bouysse & Westercamp 1990) does not manifest itself in

the observed *b*-value pattern. We do not see a general increase in *b*-values towards the southern end of the arc, which suggests little correlation with the amount of subducted sediments. Instead, we argue that there is a correlation between *b*-value and the location of fracture zones being subducted at the LAA.

In several locations the higher *b*-value bull's-eyes appear to coincide with places where fracture zones reach the arc. This pattern is particularly well observed with SRC2008; where the large bull's-eye at 15°N corresponds to the incoming Marathon and Mercurius Fracture Zones, the bull's-eye at 17°N Fifteen-Twenty FZ, the bull's-eye at 13.8°N Vema FZ and the bull's-eye at 12.3°N Doldrums FZ (Fig. 7a). The strength and robustness of the 15°N bull's-eye may be due to the combined effect of water bought by these two closely-spaced fracture zones. The relationship is demonstrated in Fig. 10, where *b*-value versus distance from the closest fracture zone is
<figure><img src="image_16.png" /><figcaption>**Figure 8.** Results for the along-arc observation using fixed cell sizes. The observation was done using OVS1996 (a) and SRC2008 (b). Column 1 (subimage on the left-hand side): Positions of volcanoes (red triangles) and islands (black rectangles) along the arc. The horizontal areas in grey represent the areas around the subducted fracture zones including the associated core complexes. Column 2: Calculated *b*-values for one cell width with its uncertainties. The average *b*-values in areas where of subducted fracture zones and the remaining areas are calculated using all results in the grey shaded areas (‘frac’) and the white gaps (‘nofrac’). Note the change of the *b*-value upper limit between the two catalogues. Column 3: Results without uncertainties for a range of different cell widths. Note that the range of values is generally less than the calculated uncertainties for a single cell size. This implies that cell size is not a critical parameter in this study. Column 4: Number of events per cell. The dashed black line shows the minimum required number for the *b*-value calculation (50). Note that there is no set maximum number of events per cell; we only cut the figure at 500 events per cell for clarity.</figcaption></figure>

plotted for the entire dataset. Even given the uncertainties and errors in both the *b*-value calculation and fracture zone positioning within the subduction zone, the relationship is clear. An explanation can be found in water delivered to the subduction zone via

fracture zones. The spatial correlation between fracture zones and higher *b*-values is consistent with a recent proposal of Manea *et al.* (2014), using a combination of arc magma geochemistry and geodynamic modelling to show that fracture zones play a fundamental
<figure><img src="image_17.png" /><figcaption>**Figure 9.** Results for the along-arc observation using fixed numbers of events per cell. The observation was done using OVS1996 (a) and SRC2008 (b). Column 1 (subimage on the left): positions of volcanoes (red triangles) and islands (black rectangles) along the arc. The horizontal areas in grey represent the areas around the subducted fracture zones within location uncertainties. Column 2: calculated *b*-values for one number of events per cell with its uncertainties. The average *b*-values in areas where of subducted fracture zones and the remaining areas are calculated using all results in the grey shaded areas ('frac') and the white gaps ('nofrac'). Column 3: required size of the cell to reach the defined number of events. The dashed black line shows the maximum allowed cell width for the *b*-value calculation (100 km).</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>**Figure 10.** Relationship of *b*-values to the distance from the nearest fracture zone for SRC2008 using the KS-test with a maximum-likelihood method (see Figs S13 and S14 for other calculation approaches and catalogues). The density of the results is shown with grey shading, the average *b*-values (bold line) and their uncertainties (fine lines) were calculated in 10 km steps. The smaller the distance from a fracture zone the greater the *b*-value.</figcaption></figure>

role in supplying water for arc magmatism. The relative variation in chemical tracers of serpentine, such as boron (e.g. Manea *et al.* 2014) or of halogens (e.g. Scambelluri *et al.* 2004; Kendrick *et al.* 2012), in volcanic rocks silicate melt inclusions could be used in future studies to establish a geochemical link between water released serpentine and the geochemistry of arc magmas.

The *b*-value variation correlates with the fracture zones in the LAA. This supports the assumption that a significant amount of water is transported along the fracture zones and influences the subduction mechanics. By the same token, an increased amount of water can trigger enhanced degrees of melting in the overlying mantle wedge. The generally high *b*-value around the island of Dominica at 15°N is accompanied by the highest density of volcanoes in the LAA (Wadge & Shepherd 1984).

The large high *b*-value region is located close to the marked kink in the volcanic arc, where the Limestone Caribees and Volcanic Caribees diverge. The region appears to underlie the island of Dominica, one of the largest and most active volcanic islands in the Lesser Antilles with at least six major eruptive centres active since the upper Pleistocene. Dominica has the highest rate of magma production in the Lesser Antilles within the past 100 000 yr with around 250 km³ of eruptive products (terrestrial and submarine). Some of the most voluminous pyroclastic deposits in the Caribbean originated from Dominican volcanoes and at least eight volcanic centres are considered to be potentially active today (Smith *et al.* 2013).

The amount of magma generated beneath any particular volcanic island is related ultimately to the amount of melting in the mantle wedge, modulated to some extent by the efficiency of the differentiation process, that is the ratio between magma intruded into the crust and magma extruded at the surface. Assuming that this ratio is broadly constant along the arc, it would appear that the elevated *b*-values beneath the central part of the arc are related to increased magma production in the mantle wedge, that itself is in large part a consequence of the amount of fluid added. Thus, our approach

suggests that there may be a link between seismic *b*-values, H₂O supply and magma productivity. This could be investigated at other arcs where variations in volcanic island (or edifice) dimensions vary along-arc.

Comparing our results with studies from the Pacific region is not straightforward: oceanic lithosphere in the Atlantic was formed at slow spreading rates compared to the intermediate-to-fast rates found in the Pacific and subduction beneath the Caribbean arc is slow compared with much higher subduction velocities in the Pacific. Alt *et al.* (2013) estimate a total volume of serpentinite production of 0.39–0.86 km³ y⁻¹ for all slow spreading ridges compared to only 0.08 km³ yr⁻¹ for intermediate and fast spreading ridges even though the total oceanic crustal production there is four times as high. Oceanic lithosphere formed at slow spreading ridges is therefore wetter and has a stronger along-strike variation in water content than oceanic lithosphere formed at intermediate-to-fast spreading ridges. In detail it is likely that the water distribution in the downgoing plate varies considerably from arc to arc, depending on how the water is incorporated into the slab (e.g. fracture zones versus plate bending), the thickness of sediment and the thermal structure. Nonetheless, our methodology is readily transferable to other subduction zones where the earthquake catalogue is similarly extensive.

## 7 CONCLUSIONS

Water appears to play an important role in the regulation of subduction dynamics both in terms of induced seismicity and in terms of arc volcanism. Using a set of two individual local seismicity catalogues we were able to map the *b*-value variation along the entire Lesser Antilles Arc for the first time. A relationship between *b*-values and the subducting fracture zones can be observed. This relationship is independent of the seismic catalogue or the calculation method. We attribute the increase in *b*-values at places where fracture zones approach the arc to the presence of greater amounts of liberated water that serve to lubricate the subduction zone and enable failure at lower differential stresses. Additional water supply may also be related to increased magma production in the wedge that, in turn, is related to volcanic output and volcano size. Further investigations in the Lesser Antilles should seek to investigate the geochemical link between water released serpentine and the geochemistry of arc magmas. Furthermore, establishing the extent of serpentinization associated with fracture zones in downgoing plates is an important next step in quantifying the water supply to arc magmatism. Such an enterprise will require continued interaction between seismology, geochemistry and petrology.

## ACKNOWLEDGEMENTS

The authors would like to thank the Seismic Research Centre of the University of The West Indies and the Institut de Physique du Globe de Paris for providing the seismic catalogues used in this study, as well as the British Geological Survey and Natural Environment Research Council for the funding of this project. We would also like to thank the anonymous reviewers for their insightful and helpful comments. JB acknowledges ERC Advanced Grant CRITMAG and Wolfson Research Merit award. JC thanks Javier Escartin and Mathilde Cannat for discussions.
## REFERENCES

Aki, K., 1965. Maximum likelihood estimate of b in the formula log N = a - bM and its confidence limits, *Bull. Earthq. Res. Inst.*, **43**, 237-239.

Alt, J.C. et al., 2013. The role of serpentinites in cycling of carbon and sulfur: seafloor serpentinization and subduction metamorphism, *Lithos*, **178**, 40-54.

Amorèse, D., Grasso, J.-R. & Rydelek, P.A., 2010. On varying b-values with depth: results from computer-intensive tests for Southern California, *Geophys. J. Int.*, **180**, 347-360.

Anderson, R., Delong, S. & Schwarz, W., 1978. Thermal-model for subduction with dehydration in downgoing slab, *J. Geol.*, **86**, 731-739.

Bender, B., 1983. Maximum-likelihood estimation of b values for magnitude grouped data, *Bull. seism. Soc. Am.*, **73**, 831-851.

Bengouhou-Valerius, M., Bazin, S., Bertil, D., Beuducel, F. & Bosson, A., 2008. CDSA: a new seismological data center for the French Lesser Antilles, *Seismol. Res. Lett.*, **79**, 90-102.

Boschi, C. et al., 2013. Serpentinization of mantle peridotites along an uplifted lithospheric section, Mid Atlantic Ridge at 11 degrees N, *Lithos*, **178**, 3-23.

Bostock, M.G., Hyndman, R.D., Rondenay, S. & Peacock, S.M., 2002. An inverted continental Moho and serpentinization of the forearc mantle, *Nature*, **417**, 536-538.

Bouysse, P. & Westercamp, D., 1990. Subduction of Atlantic aseismic ridges and late Cenozoic evolution of the Lesser Antilles island-arc, *Tectonophysics*, **175**, 349-380.

Bouysse, P., Westercamp, D. & Andreieff, P., 1990. The Lesser Antilles island arc, *Proc. Ocean Drill. Prog. Sci. Results*, **110**, 29-44.

Boynton, C., Westbrook, G., Bott, M. & Long, R., 1979. A seismic refraction investigation of crustal structure beneath the Lesser Antilles arc, *Geophys. J. R. astr. S.*, **58**, 371-393.

Briden, J.C., Rex, D.C., Faller, A.M. & Tomblin, J.F., 1979. K-Ar geochronology and palaeomagnetism of volcanic rocks in the Lesser Antilles island arc, *Phil. Trans. R. Soc. A*, **291**, 485-528.

Carpentier, M., Chauvel, C. & Mattielli, N., 2008. Pb-Nd isotopic constraints on sedimentary input into the Lesser Antilles arc system, *Earth planet. Sci. Lett.*, **272**, 199-211.

Christeson, G., Mann, P., Escalona, A. & Aitken, T., 2008. Crustal structure of the Caribbean-northeastern South America arc-continent collision zone, *J. geophys. Res.: Solid Earth*, **113**, doi:10.1029/2007JB005373.

De Gori, P., Akinci, A., Lucente, F.P. & Kilic, T., 2014. Spatial and temporal variations of aftershock activity of the 23 October 2011 M-w 7.1 Van, Turkey, earthquake, *Bull. seism. Soc. Am.*, **104**, 913-930.

DeMets, C., Gordon, R., Argus, D. & Stein, S., 1990. Current plate motions, *Geophys. J. Int.*, **101**, 425-478.

DeMets, C., Jansma, P., Mattioli, G., Dixon, T., Farina, F. & Bilham, R. 2000. GPS geodetic constraints on Caribbean-North America plate motion, *Geophys. Res. Lett.*, **27**, 437-440.

Enescu, B., Mori, J., Miyazawa, M. & Kano, Y., 2009. Omori-Utsu law c-values associated with recent moderate earthquakes in Japan, *Bull. seism. Soc. Am.*, **99**(2A), 884-891.

England, P., Engdahl, R. & Thatcher, W., 2004. Systematic variation in the depths of slabs beneath arc volcanoes, *Geophys. J. Int.*, **156**, 377-408.

Faugères, J., Gonthier, E., Griboulard, R. & Masse, L., 1993. Effects of the free surface on shear wave trains, *Mar. Geol.*, **110**, 110-142.

Frohlich, C. & Davis, S., 1993. Teleseismic b-values - or, much ado about 1.0, *J. geophys. Res.: Solid Earth*, **98**, 631-644.

Gebco, 2008. The GEBCO 08 Grid, version 20081212.

Grove, T., Till, C., Krawczynski, M. & Jeanloz, R., 2012. The role of H₂O in subduction zone magmatism, *Annu. Rev. Earth planet. Sci.*, **40**, 413-439.

Gudmundsson, O. & Sambridge, M., 1998. A regionalized upper mantle (RUM) seismic model, *J. geophys. Res.: Solid Earth*, **103**, 7121-7136.

Gutenberg, B. & Richter, C.F., 1954. *Seismicity of the Earth, and Associated Phenomena*, 2nd ed., Princeton Univ. Press.

Hacker, B., 2008. H₂O subduction beyond arcs, *Geochem. Geophys. Geosyst.*, **9**, doi:10.1029/2007GC001707.

Hattori, K. & Guillot, S., 2003. Volcanic fronts form as a consequence of serpentinite dehydration in the forearc mantle wedge, *Geology*, **31**, 525-528.

Huang, J. & Turcotte, D., 1988. Fractal distributions of stress and strength and variations of b-value, *Earth planet. Sci. Lett.*, **91**, 223-230.

Katsumata, K., 2006. Imaging the high b-value anomalies within the subducting Pacific plate in the Hokkaido corner, *Earth Planets Space*, **58**, e49-e52.

Kawakatsu, H. & Watada, S., 2007. Seismic evidence for deep-water transportation in the mantle, *Science*, **316**, 1468-1471.

Kearey, P., 1974. Gravity and Seismic-Reflection Investigations into Crustal Structure of Aves-Ridge, Eastern Caribbean, *Geophys. J. R. astr. Sci.*, **38**, 435-448.

Kendrick, M.A., Woodhead, J.D. & Kamenetsky, V.S., 2012. Tracking halogens through the subduction cycle, *Geology*, **40**, 1075-1078.

Kopp, H. et al., 2011. Deep structure of the central Lesser Antilles Island Arc: relevance for the formation of continental crust, *Earth planet. Sci. Lett.*, **304**, 121-134.

Latchman, J.L., Morgan, F.D.O. & Aspinall, W.P., 2008. Temporal changes in the cumulative piecewise gradient of a variant of the Gutenberg-Richter relationship, and the imminence of extreme events, *Earth-Sci. Rev.*, **87**, 94-112.

Lynden-Bell, D., 1971. Method of allowing for known observational selection in small samples applied to 3CR quasars, *Mon. Not. Roy. Astron. Soc.*, **155**, 95-118.

Manea, V.C., Leeman, W.P., Gerva, T., Manea, M. & Zhu, G.Z., 2014. Subduction of fracture zones controls mantle melting and geochemical signature above slabs, *Nat. Commun.*, **5**, doi:10.1038/ncomms6095.

Matthews, K., Muller, R., Wessel, P. & Whittaker, J., 2011. The tectonic fabric of the ocean basins, *J. geophys. Res.: Solid Earth*, **166**, doi:10.1029/2011JB008413.

Mauffret, A. & Leroy, S., 1997. Seismic stratigraphy and structure of the Caribbean igneous province, *Tectonophysics*, **283**, 61-104.

Mogi, K., 1962a. Study of elastic shocks caused by the fracture of heterogeneous materials and its relations to earthquake phenomena, *Bull. Earthq. Res. Inst. Univ. Tokyo*, **40**, 125-173.

Mogi, K., 1962b. Magnitude-frequency relation for elastic shocks accompanying fractures of various materials and some related problems in earthquakes, *Bull. Earthq. Res. Inst. Univ. Tokyo*, **40**, 831-853.

Neill, I., Kerr, A.C., Hastie, A.R., Stanek, K.-P. & Millar, I.L., 2011. Origin of the Aves Ridge and Dutch-Venezuelan Antilles: interaction of the Cretaceous "Great Arc" and Caribbean-Colombian Oceanic Plateau?, *J. Geol. Soc.*, **168**, 333-348.

Nishikawa, T. & Ide, S., 2014. Earthquake size distribution in subduction zones linked to slab buoyancy, *Nat. Geosci.*, **7**, 904-908.

Palacios, P., Molina, I. & Segovia, M., 2006. *The Gutenberg-Richter Law: Assumptions, Limitations and Interpretations, Special Publication of IAVCEI*, Vol. 1, pp. 115-127, Geological Society, London.

Patriat, M., Pichot, T., Westbrook, G., Umber, M., Deville, E., Benard, F., Roest, W.R. & Loubrieu, B., 2011. Evidence for Quaternary convergence across the North America-South America plate boundary zone, east of the Lesser Antilles, *Geology*, **39**, 979-982.

Pichavant, M., Mysen, B.O. & Macdonald, R., 2002. Source and H₂O content of high-MgO magmas in island arc settings: an experimental study of a primitive calc-alkaline basalt from St. Vincent, Lesser Antilles arc, *Geochim. Cosmochim. Acta.*, **66**, 2193-2209.

Plank, T. & Langmuir, C.H., 1998. The chemical composition of subducting sediment and its consequences for the crust and mantle, *Chem. Geol.*, **145**, 325-394.

Ranero, C., Morgan, J., McIntosh, K. & Reichert, C., 2003. Bending-related faulting and mantle serpentinization at the Middle America trench, *Nature*, **425**, 367-373.

Reynard, B., 2013. Serpentine in active subduction zones, *Lithos*, **178**, 171-185.

Rüpke, L., Morgan, J., Hort, M. & Connolly, J., 2004. Serpentine and the subduction zone water cycle, *Earth planet. Sci. Lett.*, **223**, 17-34.
Sammonds, P.R., Meredith, P.G. & Main, I.G., 1992. Role of pore fluids in the generation of seismic precursors to shear fracture, *Nature*, **359**, 228–230.

Sandwell, D., Mller, R., Smith, W., Garcia, E. & Francis, R., 2014. New global marine gravity model from CryoSat-2 and Jason-1 reveals buried tectonic structure, *Science*, **346**, 65–67.

Scambelluri, M., Fiebig, J., Malaspina, N., Muntener, O. & Pettke, T., 2004. Serpentinite subduction: implications for fluid processes and trace-element recycling, *Int. Geol. Rev.*, **46**, 595–613.

Schmidt, M.W. & Poli, S., 1998. Experimentally based water budgets for dehydrating slabs and consequences for arc magma generation, *Earth planet. Sci. Lett.*, **163**, 361–379.

Schmidt, M.W. & Poli, S., 2014. Generation of mobile components during subduction of oceanic crust, in *Treatise on Geochemistry*, Vol. 3, 2nd ed., eds Holland, H.D. & Turekian, K.K., Elsevier.

Scholz, C., 1968. The frequency-magnitude relation of microfracturing in rock and its relation to earthquakes, *Bull. seism. Soc. Am.*, **58**, 399–415.

Schorlemmer, D., Wiemer, S. & Wyss, M., 2005. Variations in earthquake-size distribution across different stress regimes, *Nature*, **437**, 539–542.

Shi, Y. & Bolt, B., 1982. The standard error of the magnitude-frequency b-value, *Bull. seism. Soc. Am.*, **72**, 1677–1687.

Skora, S. & Blundy, J., 2010. High-pressure hydrous phase relations of radiolarian clay and implications for the involvement of subducted sediment in arc magmatism, *J. Petrol.*, **51**, 2211–2243.

Smith, W., 1981. The b-value as an earthquake precursor, *Nature*, **289**, 136–139.

Smith, D., Escartin, J., Schouten, H. & Cann, J., 2006. Fault rotation and core complex formation: significant processes in seafloor formation at slow-spreading mid-ocean ridges (Mid-Atlantic Ridge, 13 degrees-15 degrees N), *Geochem. Geophys. Geosyst.*, **9**, doi:10.1029/2007GC001699.

Smith, A.L., Roobol, M.J., Mattioli, G.S., Fryxell, J.E., Daly, G.E. & Fernandez, L.A., 2013. The Volcanic Geology of the Mid-Arc Island of Dominica, Vol. 496, Geological Society of America.

Tatsumi, Y., 1989. Migration of fluid phases and genesis of basalt magmas in subduction zones, *J. geophys. Res.: Solid Earth Planet*, **94**, 4697–4707.

Tormann, T., Enescu, B., Woessner, J. & Wiemer, S., 2015. Randomness of megathrust earthquakes implied by rapid stress recovery after the Japan earthquake, *Nat. Geosci.*, **8**, 152–158.

Ulmer, P. & Trommsdorff, V., 1995. Serpentine stability to mantle depths and subduction-related magmatism, *Science*, **268**, 858–861.

Urbancic, T., Trifu, C.I., Long, J. & Young, R., 1992. Space-time correlations of b values with stress release, *Pure appl. Geophys.*, **139**, 449–462.

Utsu, T., 1965. A method for determining the value of b in a formula log n = a – bM showing the magnitude-frequency relation for earthquakes (with English summary), *Geophys. Bull. Hokkaido Univ.*, **13**, 99–103.

van Benthem, S., Govers, R., Spakman, W. & Wortel, R., 2013. Tectonic evolution and mantle structure of the Caribbean, *J. geophys. Res.: Solid Earth*, **118**, 3019–3036.

Vorobieva, I., Narteau, C., Shebalin, P., Beauducel, F., Nercessian, A., Clouard, V. & Bouin, M.-P., 2013. Multiscale mapping of completeness magnitude of earthquake catalogs, *Bull. seism. Soc. Am.*, **103**, 2188–2202.

Wadge, G. & Shepherd, J., 1984. Segmentation of the Lesser Antilles subduction zone, *Earth planet. Sc. Lett.*, **71**, 297–304.

White, W.M. & Dupre, B., 1986. Sediment subduction and magma genesis in the Lesser Antilles—isotropic and trace-element constraints, *J. geophys. Res.: Solid Earth Planets*, **91**, 5927–5941.

Wiemer, S. & Benoit, J., 1996. Mapping the b-value anomaly at 100 km depth in the Alaska and New Zealand subduction zones, *Geophys. Res. Lett.*, **23**, 1557–1560.

Wiemer, S. & Wyss, M., 1997. Mapping the frequency-magnitude distribution in asperities: an improved technique to calculate recurrence times?, *J. geophys. Res.: Solid Earth*, **102**, 15 115–15 128.

Wiemer, S. & Wyss, M., 2000. Minimum magnitude of completeness in earthquake catalogs: examples from Alaska, the western United States, and Japan, *Bull. seism. Soc. Am.*, **90**, 859–869.

Williams, M. & Le Calvez, J., 2013. Reconstructing frequency-magnitude statistics from detection limited microseismic data, *Geophys. Prospect.*, **61**, 20–38.

Wyss, M., 1973. Towards a physical understanding of earthquake frequency distribution, *Geophys. J. R. astr. Soc.*, **31**, 341–359.

Wyss, M., Shimazaki, K. & Wiemer, S., 1997. Mapping active magma chambers by b values beneath the off-Ito volcano, Japan, *J. geophys. Res.: Solid Earth*, **102**, 20 413–20 422.

## SUPPORTING INFORMATION

Additional Supporting Information may be found in the online version of this paper:

**Figure S1.** Vertical gravity gradient in the North-Atlantic with (B) showing the interpretation of Mid-Atlantic Ridge, transform faults and fracture zones in white. See Fig. 1 in the main paper for details. The fracture zones have been identified using both the gravity data and shipboard magnetic data. Only first-order ridge discontinuities have been mapped, which are assumed to be symmetrical on either side of the ridge. Previous workers have suggested that the NAmsAm-Af triple junction migrated northwards during the time period shown, today residing close to the Fifteen-Twenty FZ, with small-scale compression across the SAmsAm plate boundary (Pichot et al. 2012). Such details are not important at the scale of our study.

**Figure S2.** Individual plots of earthquake-magnitude distribution using the SRC2008 catalogue. The top panel shows a cell with a high b-value, the bottom panel shows a cell with a low b-value. The blue data points indicate those used in the b-value calculation, with the red solid line being the best-fit straight line and the dashed lines the uncertainties to 5 per cent significance.

**Figure S3.** Graphical summary of the techniques used to estimate b-values and a key to how they are applied in each figure. The legend is shown in A. The shape (square or circle) is used to indicate the calculation method (see text for detail). A combination of color, texture and boundary style is used to indicate, respectively, the seismicity catalogue, the cell geometry and the cell size. The approach used in each figure, is indicated in B; for example, Figure 4(A) uses the SRC2008 catalogue, the cells are circular on a grid search and the cell size is constant (referred to as Case 1 in the text).

**Figure S4.** Results for the grid-search observation for SRC2008 using the KS-test with the maximum-likelihood approach. (A) Calculations using a fixed cell size. (B) Using a fixed number of events per cell. Left column: Calculated b-values. Right column: Uncertainty of b. Note that since we do not have symmetrical uncertainties, we use the lower uncertainty boundary for results above the median b-value, whereas for b-values below the median we use the upper uncertainty boundary.

**Figure S5.** Results for the grid-search observation for OVS2008. The calculations were done using a KS-test with the maximum-likelihood approach. See Suppl. Fig. 4 for further explanations.

**Figure S6.** Results for the grid-search observation for OVS1996. The calculations were done using a KS-test with the maximum-likelihood approach. See Suppl. Fig. 4 for further explanations. The disrupted appearance of the uncertainty map for the fixed cell size calculation (top right) is not a result of any unstable calculation method but rather of the low uncertainties of the b-values.

**Figure S7.** Results for the grid-search observation using OVS1996. The calculations were done using a KS-test with explicit line fitting in log space. See Suppl. Fig. 4 for further explanations. The instability of the explicit line fitting technique shows in the disrupted appearance of the b-value maps (left).

**Figure S8.** Results for the grid-search observation using OVS1996. The calculations were done using a LB-algorithm
with maximum-likelihood method. See Suppl. Fig. 4 for further
explanations.

**Figure S9.** Results for the grid-search observation using OVS1996. The calculations were done using a LB-algorithm with explicit line fitting in log-space. See Suppl. Fig. 4 for further explanations.

**Figure S10.** Results for the along-arc observation using fixed cell sizes. This figure extends Fig. 8 of the main paper by including OVS2008 (B). OVS2008 is a subset of OVS1996 (A) decimated to match the temporal range of SRC2008 (C). Tests on this dataset were made to ensure the *b*-value results were not biased by the earthquake catalogue used. Figure details as main paper Fig. 8, with the horizontal areas in grey representing the areas around the subducted fracture zones including the core complexes. Column 1 (sub image on the left): Positions of volcanoes (red triangles) and islands (black rectangles) along the arc. Column 2: Calculated *b*-values for one cell width with its uncertainties. The average *b*-values in areas which subducted fracture zones are calculated using all results in the grey shaded areas ('frac') with the remaining areas the white gaps ('nofrac'). Column 3: Results without uncertainties for a range of different cell widths. Column 4: Number of events per cell. The dashed black line shows the minimum required number for the *b*-value calculation (50).

**Figure S11.** Results for the along-arc observation using fixed num-
bers of events per cell. This figure extends Fig. 9 of the main paper
by including OVS2008 (B). All other details as the previous Figure.
Here, column 3 shows the required size of the cell to reach the de-
fined number of events. The dashed black line shows the maximum
allowed width for the *b*-value calculation (100 km). Both catalogues
show a significantly higher *b*-value around the island of Dominica.
For OVS1996 further significantly higher *b*-values can be identified
around the island of Montserrat and the South of Martinique.

**Figure S12.** Results for the along-arc observation for OVS1996 using different methods. The upper row (A–D) uses a fixed cell size of 100 km; the lower row (E–H) uses a fixed number of 300 events. The left column shows the positions of volcanoes (red trian-

gles) and islands (black rectangles) along the arc. The methods are
Kolmogorov-Smirnov (KS) test with a maximum-likelihood method
(A, E); KS-test with explicit line fitting (elf) in log-space (B, F);
Lynden-Bell (LB) algorithm with ml-method (C, G); LB-algorithm
with elf (D, H). Even though the average *b*-value and the strength of
the bulls eyes regions of elevated *b*-values change with the calcula-
tion methods they can be observed in all cases and are significant
in most of them.

**Figure S13.** Results for the calculations using a fixed cell size with OVS1996. Relation of *b*-values and their uncertainties to the cell size for the along-arc calculation (left column) and the distance to the nearest fracture zone of the grid-search (right column). In the right column the density of the results is shown with grey shading, the lines mark the average *b*-values and their uncertainties which were calculated in 10 km steps. The calculation methods are: KS-test with ml-method (A, E); KS-test with elf (B, F); LB-algorithm with ml-method (C, G); LB-algorithm with elf (D, H). The consistency of the results obtained by the different methods gives confidence to the results.

**Figure S14.** Results for the calculations using a fixed number of events with OVS1996. Relation of *b*-values and their uncertainties to the cell size for the along-arc calculation (left column) and the distance to the nearest fracture zone of the grid-search (right column). In the right column the density of the results is shown with grey shading, the average *b*-values and their uncertainties were calculated in 10 km steps. The calculation methods are: KS-test with ml-method (A, E); KS-test with elf (B, F); LB-algorithm with ml-method (C, G); LB-algorithm with elf (D, H).

(http://gji.oxfordjournals.org/lookup/suppl/doi:10.1093/gji/ggv509/-/DC1).

Please note: Oxford University Press is not responsible for the con-
tent or functionality of any supporting materials supplied by the
authors. Any queries (other than missing material) should be di-
rected to the corresponding author for the paper.