F San

• A Gaujard

S Angélibert

S Mastrorillo

R Céréghino
email: cereghin@cict.fr

F Santoul

A Gaujard

Gravel pits support waterbird diversity in an urban landscape

Keywords: Artificial wetlands, Probability of presence, Rare species, Species richness, Self-Organizing Maps, W aterbirds

W e assessed the benefit of 11 gravel pits for the settlement of waterbird communities in an urbanized area lacking natural wetlands. Gravel pits captured 57% of the regional species pool of aquatic birds. W e identified 39 species, among which five were regionally rare. We used the Self-Organizing Map algorithm to calculate the probabilities of presence of species, and to bring out habitat conditions that predict assemblage patterns. The age of the pits did not correlate with assemblage composition and species richness. There was a positive influence of macrophyte cover on waterbird species richness. Larger pits did not support more species, but species richness increased with connectivity. As alternative

Introduction

While human activity has resulted in the destruction of natural wetlands [START_REF] Hull | The pond life project: a model for conservation and sustainability[END_REF], artificial pools, such as farm ponds, rice fields, etc., became important alternative habitats for the pond biota (Declerck et al., 2006;Céréghino et al., 2008). It also becomes increasingly accepted that man-made ecosystems are likely to support biodiversity while they provide resources that have economie values, calling for more attention on the importance of these ecosystems to both wildlife and people [START_REF] Odling-Smee | Dollars and sense[END_REF]. Throughout the world, gravel pits contribute to local economies while they constitute new wetlands for species of conservation interest, notably dragonflies, amphibians and birds (Frochot & Godreau, 1995). In France, gravel pits cover an area of about 90,000 ha, and about 5,000 ha are still created each year to satisfy the demands of the construction trade (Barnaud & Le Bloch, 1998). These new wetlands are colonized by waterbirds, creating a situation which raises new management concerns. We performed a test to verify the hypothesis that gravel pits contribute to waterbird diversity in urban areas. We recorded the species occurring in 11 gravel pits in the suburbs of the city of Toulouse (SW France), and the regional species pool of aquatic birds (gamma diversity, SW France 57,000 km 2 ) was used to assess the benefit of gravel pits for the diversity of waterbird communities. In arder to maximize the information extracted from "simple" presence-absence data, we used the Self-Organizing Map (SOM, neural network) to calculate the probabilities of the presence of each species in the various clusters of pits. Subsequently, environmental variables were introduced into the SOM trained with biological variables to interpret the variability of waterbird communities with respect to habitat features.

Methods

Several gravel pits were excavated in the River Garonne floodplain around the city of Toulouse, SW France (Fig. 1 ), a highly populated area ( over 1 million inhabitants) from which natural wetlands were largely eliminated by drainage during the nineteenth and twentieth centuries. Eleven gravel pits were characterized with four environmental variables: age (years), surface area (ha), % macrophyte caver ( estimated from point 30 samples using an Eckman grab, each year at the end of the spring period), and an index of connectivity C, calculated as C = 2::7= 1 Si x Di where S (surface) and D (distance) are divided into classes defined in [START_REF] Oertli | Diversité biologique et typologie écologique des étangs et petits lacs de Suisse[END_REF]. The use of simple variables was intended to keep models broadly applicable for management applications. All gravel pits had a similar mean depth of 3-4 m.

Bird censuses were carried out weekly from October 1996 to October 1998 using binoculars (8 x 30) and telescopes (20 x 60). Waterbirds were selected according to [START_REF] Gillier | Les comptages d'oiseaux d'eau hivernant en France: actualisation des connaissances, effectifs moyens, critères numériques d'importance internationale et nationale[END_REF]. The adequacy of sampling was assessed by plotting the cumulative frequency of species against sampling effort (samplerarefaction curve with 500 randomizations) [START_REF] Colwell | Iaterpolating, extmpolaûng. aud wmpared ゥ ョ 」 ゥ 、 」 ョ セ 「 。 ウ ・ 、 @ species accumulation c:urves[END_REF]). The regional species pool (after [START_REF] Maurel | Liste de collemboles déterminants[END_REF]) consists of 68 waterbird species occurring in SW France (see Table 1).

A full description of the SOM modellng procedure was given in Céréghino et al. (2008). The structure of the SOM consisted of two layers of neurons connected by weights (connection intensities): the input layer constituted by 39 neurons (one by species) connected to the 11 gravel pits, and the output layer constituted by 15 neurons (visualized as hexagonal cells) organized on an array with five rows and three columns. The SOM plots the similarities of the data in a 2D grid, by grouping similar data items together through an iterative learning process. At the end of the training, the connection intensity between input and output layers calculated during the learning process can be considered as the probability of occurrence of each species in the area concemed (see [START_REF] Céréghino | Using self-organizing maps to investigate spatial patterns of non-nal:i.ve species[END_REF]. The occurrence probabilities of each species were visualized on the SOM map in grey scale, and allowed us to analyse the effect of each species on the patterning input dataset (sites). The SOM was clustered using Ward's algorithm [START_REF] Leroy | Morphogenetic trends in morphological, optical, biochemical features of phyllodes in Acacia mangium Willd (Mirnosaceae)[END_REF]. In arder to bring out relationships between biological and environmental variables, we introduced the environmental variables into the SOM previously trained with bird occurrences

+ +

The 39 species recorded in our 11 gravel pits (Toulouse suburbs) appear in bold, + = presence. Information on rarity and commonness at the regional scale is given (after [START_REF] Maurel | Liste de collemboles déterminants[END_REF] (see [START_REF] Park | Applications of artificial neural networks for patterning and predicting aquatic insect species richness in running waters[END_REF]. Ali mean values of environmental variables assigned on the SOM map were visualized in grey scale.

Results

Fifty-seven percent of the regional species pool was captured by our 11 gravel pits (Table 1 ). The number of species per pit ranged from 12 to 38. Among the 39 species recorded, five species were regional! y rare, and the remaining ones were common. Accumulation of new species reached its asymptote (Fig. 2); we could th us consider that our sampling was satisfactory. Mter training, the SOM with species occurrences, Ward's algorithm helped to derive three clusters of sites (gravel pits) according to waterbird assemblages (Fig. 3a). The ordinate of the SOM represented the number of species, from low (top areas of the SOM) to high (bottom) (Fig. 3b ). Eleven species (Podiceps auritus, Platalea leucorodia, Aythya marila, Netta rufina, Mergus albellus, Glareola pratincola, Tringa totanus, Tringa nebularia, Actitis hypoleucos, Limosa lapponica, and Calidris alpina) only occurred in cluster C (Fig. 4). Five species (Podiceps cristatus, Phalacrocorax carbo, Ardea cinerea, Egretta garzetta, and Anas platyrhynchos) occurred in ali the gravel pits, and the remaining 23 species occurred in two clusters of sites. When environmental variables were introduced into the SOM (Fig. 5), the ordinate on the SOM showed a gradient of connectivity and macrophyte cover [from low (top area) to high (bottom)]. When clusters of gravel pits were compared to their distribution on the geographical map of the study area, the species hosted by sorne pits within the same sub-areas tended to be similar (gravel pits 7 and 9 in cluster B; 5, 6 and 10, and 2, 3 and 11 in cluster A, 1 and 4 in cluster C, see Fig. 1) and those characteristics tended to differ when sites belonged to more distinct areas. However, gravel pits from distinct areas also had very similar assemblages (e.g., 8 and 4; 3 and 5-6; 11 and 10), thus suggesting that local factors (e.g., macrophyte caver) interacted with connectivity to shape bird assemblages.

"' 40

Q) • • • '(3 • • & 35 • • "' • 0 30

Discussion

While the potential pool of colonists in our study region was made of 68 waterbird species, 11 gravel pits allowed the presence of 39 species in the Toulouse city suburbs. Our observations therefore highlight how a small set of artificial wetlands may sustain an important fraction (57%) of the larger regional species pool in landscapes where natural wetlands are lacking. Other observations in France and UK showed increases in population densities of several waterbirds during the past decades (Aythya fuligula, Anas strepera, Phalacrocora:x carbo, Podiceps cristatus, Stema hirundo), as a result of the creation of gravel pits and reservoirs [START_REF] Frochot | L'expansion réc:ente des oiseaux d'e1111[END_REF][START_REF] Fuller | Birds and habitat change in Britain. Part I: a review of losses and gains in the twentieth century[END_REF]. Species richness and the presence of rare species are frequently cited criteria for site selection by conservationists [START_REF] Myers | Biodiversity hotspots for conservation priorities[END_REF]. If man-made habitats only attract the common (widespread) species, one may argue that they do not make a significant contribution to biodiversity. Conversely, rare species are of special interest to environmental managers [START_REF] Rey-Benayas | Commonness and rarity: theory and application of a new model to Mediterranean montane grasslands[END_REF], and it was recently demonstrated that areas which carry rare species may also concentrate an important fraction of the regional biodiversity (Cucherousset et al., 2008). For instance, sorne species were exclusive to gravel pit 4, which also hosted the highest species richness. Such artificial wetlands might therefore benefit from

P. cristatus A cinerea P. carbo E. garzella A p/atyrhynchos セ Z N @ 0 A clypeata A purpurea P. leucorodil' C. alpins L. tapponica L. cachinnans G. pratincola A lerina セ Z N @ 0 A querquedula A pene/ope A strepera F. atra Q P N セ Y @ セ o N m T @ 0.69 A hypoleucos N. rufina T. rulicollis T. totanus T. nebularia P. auritus A mari/a 1 0.888 m r -8 m 1 0.567 セ P N T T T @ セ @ セ M P N T @ • セ P R X T @ 0 0 0 N . nycticorax B. ibis セ @ FA セ Z N @ セ セ @ セ Z Z Z Q @ セ Z @ 0 0 0.727 0.281 0.4 10962 セ P N W S R @ 0.501 A creees A BCU/8 H. h imantopus G. grus G. chloropus L. ridibundus FI FI セ Z N @ 1 セ Z Z Z Z @ セ @ 0.562 0 0 0 10.9 セ P N W @ 10.983 セ P N U S S @ 0.0833
セ Z Z N ュ @ Em [ 1 セ ] @ セ P N T Y U @ 0
Dg. 4 Gradient malysis of the pmbability of occ:uroence of each watetbild speciea cm the 1IlÙDed SOM (see Fig. 3a), with v.i1malization in shading scale (tklrk = bigh probability of occummce, light = low probability of occurrence). l!ach map bigher management priorities. We thus support the idea that urban landscapes containing man-made wetlands make a significant contribution to freshwater biodiversity. Larger gravel pits did not suppon more waterbird specles. This absence of species-erea relationship suggests tbat larger gravel pits were not more easily colonized by immigrants, and/ar were not likely to show a higher diversity of ecological niches facilitating the coexistence of a larger number of species. Smaller but well-connected gravel pits ( セ T @ ha) bad tbe greatest SW!Ceptibility to host mme taxa, and more rare species. They potentially bad higher conservation value for waterbirds than Iarger gravel pits, although smal1 surface may become a limiting factor if the carrying ca.pacity becomes insufficient for waterbirds. Thus, gravel pits were certainly attractive to waterbirds when they acted as stepping stones that ensured connectivity between ia superimposed on the map n:pacnting the distribution of sites poesented in. Fig. 3a. See Table 1 for fuD apecies 118DleS and authorities Iarger natural and/or artificial wetlands separated in space [START_REF] Boumaud | L'espace étang daDs ses rapporta avec l'avifawle en période de nidüicatiOD[END_REF]. Species-poor assemblages were subsets of richer assemblages, suggesting nested patterns of waterbird assemblages. We assume that such patterns would be 」 ッ ャ ッ ョ ゥ コ 。 エ ゥ ッ ョ セ カ ・ ョ @ because most waterbirds did not live at the gravel pits throughout the year. The geographicallocation of the study area near the Pyrenees mountainous barrier maloes the region important as a stop-over for migrant birds (Royer, 1994). Rare species were specifically present during the stop-over period, and therefore they pmbably preferred those connected sites which allowed them to find quieter sites in case of disturbance. There was a positive inftuence of the extent of tbe macrophyte cover on the number of waterbird species present at a site. Aquatic plants over large areas are attractive to many birds, such as charophytes, and the many invert.ebmtes living on them provide food for du.cks and other waterbirds [START_REF] Knapton | Changes in distribution and abundance of submerged macrophytes in the inner bay at Long Poing, Lake Erie: implications for foraging waterfowl[END_REF][START_REF] Blindow | Lcmg-teml waterfowl ductuali01111 in relation to altemative states in two sballow laloes[END_REF]. In tum, aquatic birds may have a negative impact on macrophyoe abondance, and may be important determinants of the aquatic system dynamics (reviewed in [START_REF] Lodge | Impact of herbivory on plant standing crop: comparisons among biomes, between vascular and among freshwater herbivory taxa[END_REF]. The age of the gravel pit was not a structuring factor for waterbird assemblages, at least within the range of 3-22 years considered in this study (no gravel pits older than 22 years in the area). Nevertbeless, Frochot & Godreau (1995) emphasized that gravel pits over 30 years are less attractive for watel'birds, due to their more homogeneous and dense habitats, although sorne species such as herons may benefit from the development of trees.

In conclusion, gravel pits should be seen as a network of habitats integratecl within the broader network of natural and artificial wetlands. We must enlarge the spatial scale needed to manage/survey these particular habitats, and we should primarily pay attention to welland networks rather than attempting to target some species and/or bodies of water for parti.cular mana.gement actions. 'lbe lack of fundamelllal knowledge needed to implement managenumt plans usually limits the conservation potential of gravel pits. Restoration is generally car.ried out to transform gravel pits into recreational areas or ftshing lakes, and oDly rarely into habitats favorable to waterbirds. Owing to the continuing Joss of natural wetlands, there is a need to enhance the contribution of artüicial wetlands such as gravel pits for future conservation of waterbirds, and probably other taxa such as amphibians, insects, or plants. Therefore, further understanding of the distribution of biologi.cal diveuity in non-natural systems may facilitate the adoption of positive solutions for wildlife, with limitecl costs for human activities.

Fig. 1

 1 Fig. 1 Location of the 11 gravel pits (1-11) around the city of Toulouse, SW France

Fig. 2

 2 Fig. 2 Sample-based rarefaction curves representing the number of species accumulated by sampling the 11 gravel pits

  Fig. 3 a Distribution of gravel pits on the self-organizing map (SOM) according to the presence or absence of 39 waterbird species. Solid lines show the cluster boundaries (i.e., for clusters A, B, C), delineated according to Ward's algorithm. Gravel pits that are neighbors within clusters are expected to have similar waterbird assemblages. Codes (1-11) correspond to gravel pits. b Mean number of species (±SE) per cluster

  Fig. 5 Visualization of enviroDIDCD1al variabks on the SOM trained with watelbird specics. The mean value of each varlabl.e was calculated in each output ueumn of tbe tniDed SOM. Dari J:epJaelda a high value. wbile light ia 1ow

Table 1

 1 List of the waterbird species occurring at the regional scale (SW France)

	Gravel pit

Aythya nyroca

(Gûld, 1769)

Rare

Bubulcus ibis

(Linnaeus, 1758)

Co mm on

+ + + + + + + + + Bucephala clangula

(Linnaeus, 1758) Co mm on Calidris alpina (Linnaeus, 1758) Co mm on + Calidris minuta (Leisler, 1812) Co mm on Charadrius dubius Scopoli, 1786 Co mm on + + + + Chlidonias hybridus (Pallas, 1811) Co mm on Chlidonias niger (Linnaeus, 1758) Co mm on Egretta alba (Linnaeus, 1758) Co mm on Egretta ganetta (Linnaeus, 1758) Co mm on + + + + + + + + + + + Eudromias morinellus (Linnaeus, 1758) Rare Fulica atra Linnaeus, 1758 Co mm on + + + + + + + + + + Tadorna tadorna (Linnaeus, 1758) Common + + + + Tringa erythropus (Pallas, 1764) Common Tringa glareola Linnaeus, 1758 Common Tringa nebularia (Gunnerus, 1767) Common + + Tringa ochropus Linnaeus, 1758 Common Tringa totanus (Linnaeus, 1758) Common + + Vanellus vanellus (Linnaeus, 1758) Common


P A Cowie

C H Scholz

G P Roberts

J P Faure

P Steer

Viscous roots of active seismogenic faults revealed by geologic slip rate variations

During the earthquake cycle viscous flow at depth contributes to elastic strain accumulation along seismogenic faults 1 . Evaluating the importance of this contribution to fault loading is hampered by uncertainty about whether viscous deformation mainly occurs in shear zones or by distributed flow. Furthermore, viscous strain rate has a power-law dependence on applied stress 2 but few estimates exist for the power-law exponent applicable to the long term in situ behaviour of active faults. Here we show that measurements of topography and whole-Holocene offsets along seismically active normal faults in the Italian Apennines can be used to derive a relationship between stress and strain rate (averaged over 15±3 kyrs). This relationship, which follows a welldefined power-law with an exponent in the range 3.0-3.3 (1σ), is used to infer the rheological structure of the crust and constrain the width of active extension across the Apennines. Our result supports the idea that the irregular, stick-slip movement of upper crustal faults, and hence earthquake recurrence, are controlled by down-dip viscous flow in shear zones over multiple earthquake cycles.

Earthquakes in the crust occur down to depths of approximately 15km in most regions because below this depth temperature-and time-dependent creep (aseismic) deformation processes become progressively more important. It is therefore generally accepted that the upper crustal, seismogenic, portion of a fault is rooted down dip into a ductile (mylonitic) shear zone and that the transition from frictional stick-slip to viscous flow is temperature and strain rate dependent. At sufficiently high temperatures, distributed ductile deformation may also occur in the lower crust and upper mantle. Both localised flow in shear zones and distributed flow lead to elastic strain accumulation in the upper crust and thus loading of faults to failure but currently there is disagreement as to which dominates 3, [START_REF] Freed | Earthquake triggering by static, dynamic and postseismic stress transfer[END_REF] . Experimental work, field data and theory indicate the flow law for the lithosphere at tectonic strain rates should be that of dislocation creep in which strain rate, ė, is proportional to stress raised to an exponent n, where n is typically in the range 2 to 4 5, [START_REF] Newman | The dynamics of extensional sedimentary basins: constraints from subsidence inversion[END_REF] :

      - = RT Q A e n exp σ  Equation 1
Here σ is driving stress, A is a material property, Q is activation energy, R is the molar gas constant and T is absolute temperature. Geodetic observations of post-seismic relaxation reveal temporal and spatial variations in effective viscosity that are most easily explained by power law creep with n ≈ 3 [START_REF] Freed | Evidence of power law flow in the Mojave Desert mantle[END_REF] . However geodetic data generally do not permit discrimination between contributions of bulk flow of the upper mantle, of the lower crust, or plastic creep within a shear zone [START_REF] Bürgmann | Rheology of the lower crust and upper mantle: Evidence from rock mechanics, geodesy and field observations[END_REF] . Moreover, it is not clear that the rheological properties indicated by postseismic transients are applicable to longer term behaviour of the coupled frictionalviscous fault system 8 .

Here we show that extensional strain rates derived from slip on seismogenic normal faults in the actively uplifting and extending central and southern Italian Apennines can be used to address this issue. The strain rates are measured at the surface using published structural data 9,11 (Fig. 1 (Fig. 2a). Gravity admittance data indicate that the long wavelength topography is supported by buoyancy variations in the uppermost mantle 12 . Regional surface uplift rates [START_REF] D'anastasio | Short-term vertical velocity field in the Apennines (Italy) revealed by geodetic levelling data[END_REF] increase in magnitude inwards from the Adriatic and Tyrrhenian coasts, mimicking in shape the long wavelength topography [START_REF] D'anastasio | Short-term vertical velocity field in the Apennines (Italy) revealed by geodetic levelling data[END_REF] .

The extensional strain rates, averaged over the whole Holocene (15±3 kyrs), correlate with average topographic elevation along the length of the central and southern Apennines 10 (Fig. 1c). This observation is confirmed by the map view distribution of active faults relative to topographic contours (e.g., Fig. 2a). Geodetic data also show that the highest contemporary strain rates coincide with the highest elevation area in the central Apennines [START_REF] D'agostino | Evidence for localized active extension in the central Apennines (Italy) from global positioning system observations[END_REF] . A power law regression between the strain rate, ė, and elevations, h, (in transects 90 km across strike by 30 km wide along strike) reveals a well-defined relationship with power law exponents in the range 3.0-3.3 (1σ) (see Methods and Supplementary material). Data from two independent sets of 30km transects show that the result is not location dependent (Fig. 1d). Varying transect width (from 5 km to 60 km) shows that over all scales the exponent lies within the range 2.7-3.4 and 2.3-4.0 at 95% and 99% confidence intervals respectively. These variations in strain rate cannot be attributed to thermal structure as heat flow increases gradually from < 40 mWm -2 along the Adriatic coast to > 60 mWm -2 along the Tyrrhenian coast, independent of elevation and distance along strike 15 .

To interpret our data (Fig. 1d) in terms of Equation 1 we need to demonstrate that h and σ are directly proportional. Previous workers (e.g., ref 16) made the connection by approximating the lithosphere as a homogeneous thin viscous sheet. However, where thicker than average crust (40 -50km) overlies thinned mantle lithosphere, as it does in the central and southern Apennines, the vertical velocity field is unlikely to be continuous at the scale of the entire lithosphere [START_REF] Flesch | The relationship between surface kinematics and deformation of the whole lithosphere[END_REF] . Furthermore, the topography varies by 100's of meters over wavelengths < 100 km in which case approximations made in the thin sheet model break down [START_REF] Naliboff | The effects of lithospheric thickness and density structure on Earth's stress field[END_REF] . To avoid making these approximations we use observational constraints to relate h to σ by noting that (1) the upper crust is at or close to the threshold for brittle failure, i.e., "at yield" 19 and (2) earthquake focal mechanisms and fault kinematic data along active faults [START_REF] Walker | Comparison of earthquake strains over 10 2 to 10 4 year timescales: Insights into variability in the seismic cycle in central Apennines, Italy[END_REF][START_REF] Roberts | Spatial and temporal variations in growth rates along active normal fault systems: an example from The Lazio-Abruzzo Apennines, central Italy[END_REF] indicate that the maximum compressive stress, σ 1 , is vertical and the least compressive stress σ 3 is parallel to the principal extensional strain orientation (NE-SW in Fig. 2). In an elastic-brittle upper crust at yield, σ 3 , is directly proportional to σ 1 , compatible with incipient frictional failure on optimally oriented planes 20 (Fig. 3). Thus the differential stress is also proportional to σ 1 , e.g., (σ 1 -σ 3 ) ≈ 2σ 1 /3 if Byerlee friction constants are assumed. Below the base of the seismogenic zone (~14-17 km depth in this region 19,21 ), where viscous flow dominates, differential stress is less but we assume there is no stress discontinuity across this transition over long time scales (Fig. 3). Additional topographic loads that result from surface uplift relative to sea level increase σ 1 , and hence (σ 1 -σ 3 ), driving deformation (by a depth and temperature dependent combination of frictional slip and viscous flow; Fig. 3) such that differential stress in the upper crust is relaxed to re-establish the "at yield" condition (Fig. 3c). As both frictional and viscous components of the fault system undergo the same overall strain, they operate in parallel (A, B in Fig. 3b). As long as surface uplift since extension began is proportional to elevation, h, which as explained in the Methods is a reasonable assumption 10, [START_REF] D'agostino | Interactions between mantle upwelling, drainage evolution and active normal faulting: an example from the central Apennines (Italy)[END_REF] , the increase in differential stress is independent of depth and simply proportional to ρgh (ρ is crustal density; g is acceleration due to gravity). This reasoning, which requires that buoyancy forces rather than plate boundary forces are the dominant control on upper crustal deformation, explains the spatial variation in regional strain rates, provides the link between the relationship shown in Fig. 1(d) and Equation 1, and thereby allows us to constrain the exponent n ≈ 3. Other studies that compare topography with strain rates derived from fault slip and geodetic data [START_REF] England | Late Quaternary to decadal velocity fields in Asia[END_REF] assume n ≈ 3 to explain their data but here, for the first time, we use such data to constrain its value.

Our surface strain rate measurements are derived from slip along faults so they do not represent deformation of a continuum. For example, at the scale of individual fault blocks (20km x 20km; Fig. 2a) the spatial gradients in strain rate and mean elevation (Fig. 2b) covary, consistent with the regional relationship (Fig. 1), but the correlation between ė and h is poor (R 2 < 0.5) because of the heterogeneous pattern of brittle faulting. However, is viscous deformation at depth also likely to be heterogeneous, as proposed in Fig. 3(a)? Deforming non-linear viscous materials (i.e. n > 1), in general, show a tendency to localise strain and the development of a brittle fault up-dip provides a geometrical discontinuity that influences where in the viscous regime localisation preferentially develops [START_REF] Huismans | Roles of lithospheric strain softening and heterogeneity in determining the geometry of rifts and continental margins[END_REF] . Brittle-frictional faults extending down to depths of at least 10km in this area are revealed by earthquake aftershocks [START_REF] Chiarabba | and 28 co-authors. The 2009 L'Aquila (central Italy) M w 6.3 earthquake: Main shock and aftershocks[END_REF] . During major earthquakes, cataclasis, hydrous alteration and shear heating together contribute to grain size reduction and material weakening 8 , processes associated with localisation at the frictional-viscous transition [START_REF] Fusseis | Networking of shear zones at the brittle-to-viscous transition (Cap de Creus, NE Spain)[END_REF] (Fig. 3) and enhanced shear zone development within the viscous regime down-dip of seismogenic faults [START_REF] Ellis | Imposed strain localization in the lower crust on seismic timescales[END_REF] . We also know that fine grained, poly-mineralic, mylonitic shear zones remain weaker [START_REF] Fusseis | Networking of shear zones at the brittle-to-viscous transition (Cap de Creus, NE Spain)[END_REF] and deform at higher strain rates than surrounding coarser grained less deformed rock (red line, Fig. 3a) and contribute to lowering the effective viscosity (Fig. 1d) [START_REF] Mehl | Plagioclase recrystallization and preferred orientation in layered mylonites: Evaluation of flow laws for the lower crust[END_REF] .

In the mid-to lower crust quartzo-feldspathic mylonites form a fabric of mineral segregated layers parallel to shear so that their strength is controlled by the weakest phase:

quartz. Using a flow law for wet quartz calibrated for mylonitic rocks [START_REF] Hirth | An evaluation of quartzite flow laws based on comparisons between experimentally and naturally deformed rocks[END_REF] to fit measured strain rates across brittle fault zones (~5 km wide) [START_REF] Walker | Comparison of earthquake strains over 10 2 to 10 4 year timescales: Insights into variability in the seismic cycle in central Apennines, Italy[END_REF] , we estimate a lower bound on the temperature, T, of the viscously deforming material to be 710±120K (~440 ±120°C). These temperatures are reached just below the base of the seismogenic zone (25±7.5 km), as constrained by surface heat flow data [START_REF] Vedova | Deep temperatures and surface heat flow distribution[END_REF] and upper crustal seismicity [START_REF] Boncio | The instrumental seismicity of the Abruzzo Region in Central Italy (1981-2003): seismotectonic implications[END_REF][START_REF] Chiarabba | and 28 co-authors. The 2009 L'Aquila (central Italy) M w 6.3 earthquake: Main shock and aftershocks[END_REF] (see Methods). At the 99% confidence level the exponent we derive is consistent with this flow law (Fig. 1d) as well as with recent theoretical predictions of shear zone rheology [START_REF] Platt | Grainsize evolution in ductile shear zones: Implications for strain localisation and the strength of the lithosphere[END_REF] and textural evidence 8,23,25,26,27 from poly-mineralic mylonitic rocks. We suggest therefore that the rate of viscous flow in shear zones dominates over distributed flow within the lower crust and/or upper mantle and regulates the slip rate we measure over long timescales on the up-dip seismogenic part of the fault even if on shorter time scales rupture depends on elastic loading and mechanical instabilities 3,4 .

The relationship between ė and h (Fig. 1d) also permits the controls on regional seismic hazard to be re-evaluated. Uplift of the Italian Peninsula has been on-going since at least the Quaternary, fast enough in places to raise Plio-Pleistocene marine sediments by over 400 m since the Early Pleistocene 10,28,29 . Contemporary surface uplift, documented by geodetic levelling lines across both the footwalls and hanging-walls of active faults, indicates regional uplift of 0-0.5 mm/yr close to the coasts, increasing to 1.0-1.5 mm/yr inland [START_REF] D'anastasio | Short-term vertical velocity field in the Apennines (Italy) revealed by geodetic levelling data[END_REF] . Furthermore, slip rates along the largest normal faults in the highest elevation areas of central Italy (Fig. 2a) increased significantly ~0.8-1Ma [START_REF] Roberts | Fault Slip-Rate Variations During Crustal-Scale Strain Localisation, Central Italy[END_REF] . Some hanging-wall basins became internally drained 12 in part because of high rates on basin-bounding faults, consistent with progressive strain localisation 14, [START_REF] Roberts | Fault Slip-Rate Variations During Crustal-Scale Strain Localisation, Central Italy[END_REF] . Although recent analysis of geodetic data and instrumental seismicity suggests that active deformation may be localised in a zone only ~50km wide [START_REF] D'agostino | Evidence for localized active extension in the central Apennines (Italy) from global positioning system observations[END_REF] (Fig. 4), historical earthquake shaking records (since 1349 A.D.) and the distribution of Holocene scarps imply a broader active zone 80-90km wide (Fig. 4). These apparently contradictory observations become consistent when viewed in the light of the derived relationship (Fig. 1d) as it predicts that strain rate decreases rapidly across strike from the highest elevation areas towards the coast (Fig. 4) and smaller strains are more readily resolved when observed over a longer period of observation (10 2 -10 4 yr versus < 10 1 yr in this case).

In summary, rates of seismicity on frictional faults in this example are regulated by rates of localised viscous flow at depth constrained for the first time by observations averaged over the timescale of multiple earthquake cycles and likely representing the in situ and longterm (10 4 years) mechanical properties of deforming mylonitic shear zones.

Methods:

Strain rate calculation: Reference 9 presents the calculations used to convert field measurements of the direction and amount of Holocene fault slip (since 15±3 ka) into strainrates within grid cells, the dimensions of which can be specified and thus varied. Fig. 2(a) shows principal extensional strain rates (blue bars) calculated using this approach for 20 km x 20 km grid cells for the central Apennines. The strain rate calculations follow established methods. For details see Supplementary material. To obtain the strain rate data shown in Fig. 1 the principal extensional strain rate is first calculated using a 5 km × 5 km grid oriented 8 NE-SW, approximately parallel to the principal strain orientation (043°-223°), extending a distance 90 km across strike. From these data we calculate the strain rate in adjacent rectangular regions 5km x 90km (see Fig. S1 in Supplementary material). We combine together adjacent 5km x 90km regions to obtain independent estimates of the principal strain rate in transects of a given width (10km, 20km, 30km etc.; Fig. S2) over the whole study area. Transects 30km wide suppress short length scale variations due to fault displacement gradients and constrain well the regional variations in strain rate along strike along the Italian Apennines from NW to SE (see Fig. S2). To characterise shorter wavelength spatial variations in strain rate and elevation across strike across the central Apennines we calculate Δė / Δx and Δh / Δx (Fig. 2b), by taking the difference in strain-rate (Δė) and mean elevation (Δh) between adjacent 20 km x 20 km cells (Fig. 2a) along the extension direction (043°-223°), from NE to SW, and dividing by Δx = 20 km.

Topographic data and stress:

To calculate elevation (Fig. 1) topographic profiles located along the centre of 5 km x 90 km regions (see strain rate calculation) were constructed from SRTM 90 m DEM data using GeoMapApp. Each of the topographic data profiles are orientated NE-SW and are separated along-strike by 5 km intervals. Spot heights along the topographic profiles are sampled approximately every 850 m. The 5 km width transects were combined to calculate the mean elevation within wider transects (10km, 20km, 30km etc.) to derive regional elevation variations along strike. The 95% confidence intervals of the mean elevation are calculated using the assumption of a normal distribution in the topographic spot heights [START_REF] Walker | Relationship between topography and strain rate in the actively extending Italian Apennines[END_REF] . Remnants of a flat palaeolandscape, formed close to sea-level during the Pliocene and now identified at high elevations, plus preservation of uplifted marine deposits and terraces indicate that the present day topography has mainly formed since extension began [START_REF] Walker | Relationship between topography and strain rate in the actively extending Italian Apennines[END_REF] .

Subsequent erosion/deposition is minimal and many of the high elevation hanging-wall basins remain unincised [START_REF] D'agostino | Interactions between mantle upwelling, drainage evolution and active normal faulting: an example from the central Apennines (Italy)[END_REF] . See Supplementary material for more details. Stress (MPa) and effective viscosity (η = σ/2ė) (top axis in Fig. 1d) are derived using σ = ρgh where ρ = 2800 kgm -3 and h = elevation in meters.

Correlation between strain rate and elevation: Grey lines in Fig. 1d are best fit regressions through data obtained by sampling mean elevation, h, and strain rate, ė, in 30 km wide transects (90km across strike). Data from two transect positions (offset by 15 km) are presented to demonstrate there is no selection bias. Regressions lines (Fig. 1d) are given by: ė = 10 -17.8 h 3.2 with R 2 = 0.8 (1σ) and ė = 10 -17.9 h 3.3 with R 2 = 0.9 (1σ). The estimate of the power law exponent depends on transect width; 30 km wide transects best constrain its value (see Fig. S2). The six different transect positions at this scale yield a value of the exponent in the range 3.0 to 3.3 and a pre-factor of 10 -17.6+0.3-0.9 . We use ė = 10 -17.6 h 3.2 in Figure 4 to predict the variation in strain rate from topography.

Temperature-depth calculation:

We use values for Q, log A, R and n (Eqn. 1) from a flow law for wet quartz [START_REF] Hirth | An evaluation of quartzite flow laws based on comparisons between experimentally and naturally deformed rocks[END_REF] and solve for the temperature T that predicts strain rates of similar magnitude to those measured across 5 km wide fault zones [START_REF] Walker | Comparison of earthquake strains over 10 2 to 10 4 year timescales: Insights into variability in the seismic cycle in central Apennines, Italy[END_REF] . Q = 135 ± 15 kJ/mol, log A = -11.2 ± 0.6 MPa -n /s, R = 8.314472 m 2 kg s -2 K -1 mol -1 , n = 4. We obtain T = 710±120 K, i.e., ~440 ±120°C. Average surface heat flow of 50 mWm -2 in this area [START_REF] Vedova | Deep temperatures and surface heat flow distribution[END_REF] is used to derive temperature versus depth through the upper crust assuming a surface temperature of 10°C, crustal heat production = 1 x 10 -6 Wm -3 and thermal conductivity = 2.5 J s -1 m -1 K -1 .

Following the standard approach, crustal heat production decreases exponentially with depth (characteristic depth = 10 km). The inferred depth range of viscous flow implied by the temperature range is 25±7.5 km, i.e., below the depth extent of upper crustal seismicity (14-

, 2 ;

 2 see Methods) along active normal faults, characterised by bedrock scarps that exhibit striated fault planes and offset dated Holocene sediments and geomorphic surfaces 11 . These faults have developed in the last 2-3 My since thrusting in this region diminished as westward subduction of the Adriatic plate beneath the Italian peninsula slowed and slab tearing/detachment initiated ~6 Ma 10,12 . Present day topographic elevation increases inland reaching elevations up to 2900 m locally along the footwall crests of major extensional faults. Short wavelength (10-20 km) topographic variations due to faulting are superimposed on long wavelength (100-150 km) topography aligned NW-SE along the axis of the Italian Peninsula

Figure 1 .

 1 Figure 1. Correlation between Holocene-averaged regional extensional strain rates and

Figure 2 .

 2 Figure 2. Spatial variation in strain rate with elevation across the central Apennines

Figure 3 .

 3 Figure 3. Rheology and loading of a coupled frictional-viscous fault system. (a)

Figure 4 . 2 x 10 - 8 )

 42108 Figure 4. Contemporary strain accumulation across Abruzzo. Width of high strain rate zone (along line shown in Fig. 2a) implied by geodesy 14 (red arrow) (current resolution: ė ≥ 2 x 10 -8 ) versus fault scarp mapping (blue arrow) (estimated resolution: ė ≥ 5 x 10 -9 ). Strain rate variation (grey line) predicted from long wavelength topography 14 (thick black line) using ė = 10 -17.6 h 3.2 (see Methods). Thin black line: SRTM topography. Stars indicate locations of

Figure 2 Y

 2 Figure 2

12

17 km 19,21 ). As mantle lithosphere is thinned in this region, extrapolating this temperature depth profile may underestimate lower crust/upper mantle temperatures.

Acknowledgements

. N. D'Agostino supplied the long wavelength topography data used in Figs. 2 and 4. This work was supported by NERC grants: NER/S/A/2006/14042, NE/E01545X/1 and NE/I024127/1. Financial support was also provided by the Statoil Earth System Modelling project (P.S.) and the Statoil-University of Bergen Akademia agreement (P.C.). We thank R. Huismans for useful discussions and P. Molnar, M. Handy and an anonymous reviewer for their comments.

Supplementary Information:

Section 1: Summary of main calculation steps used in strain rate calculations plus evaluation of sources and magnitudes of uncertainties.

Section 2: Estimate of power law exponent relating strain rate to elevation using different transect widths and positions including an analysis of confidence intervals on these estimates.

Author contribution statement:

PAC led the interpretation of the scaling exponent in terms of mid-crustal shear zones.

CHS contributed to understanding the behaviour of coupled frictional-viscous fault systems.

GPR provided the structural data and analysed the strain rate vs elevation relation. JFW performed the strain rate calculations and quantified data uncertainties.

PS contributed to understanding stress and strain rate variations in a layered lithosphere.