**Origin of the lozenge-shaped topographic highs. **One possibility is that these features represent intra-transform volcanism40. However, the observed low backscatter intensity of these structures8 does not support recent volcanism. A more likely scenario is that these features represent crustal blocks exhumed via transpressional stresses underlain by altered mantle, i.e., features typically referred to as positive flower structures. Similar morphological features, are frequently observed and fault plane locations are sometimes available to substantiate a transpressional origin. This is supported by the highly tectonised mafic and ultra-mafic compositions of rocks from features with similar morphology and gravity signatures at the nearby St. Paul TF that are also thought to have a transpressional origin41.

**Seismic catalogue and focal mechanisms. **The analysis of the broadband ocean bottom seismic deployment revealed 972 events, 812 of which are located within the area surrounded by the OBS network31, which belongs to the Chain TF and the adjacent ridge spreading centres (Fig. 1b). The location of the events is performed by NonLinLoc software42 using a 1-D velocity model of the Chain TF region based on CRUST1.043 and earthquake arrival times31. In order to focus on the TF seismicity, we discard the events that unequivocally occurred at the ridge segments and at the inside corners after considering their location and focal mechanisms (north-south striking, normal faulting events). This selection leaves 626 events along the Chain transform valley, within a local magnitude (ML) range between 1.1 and 5.4. We further remove events below the completeness threshold MC=2.3 (see below), resulting to a dataset of 370 events along the transform valley, which is used for our analysis.

The epicentral coordinates of these events are sufficiently constrained with median lateral uncertainties of 2.5 km. Vertical uncertainties are larger (median ~18 km) and several events are above 20 km. For this reason, the hypocentre depths as well as focal mechanisms of well-recorded events are re-evaluated using the Grond software44, which carries out a Bayesian-bootstrap time-domain deviatoric moment tensor inversion. We calculated 119 focal mechanism solutions. For the stress inversion we select the 47 focal mechanisms classified as having good-fits (Fig. 1). We include 42 additional events (89 events total) with sufficiently well-determined depth (mean vertical error of 6km).

We supplement our seismicity data recorded by the OBS with historical earthquakes from the Global Centroid Moment Tensor catalogue45. Due to the relatively large epicentral uncertainties, we only consider the strong events (MW≥5.6) after 1993 relocated by Shi et al.26.

**OBS seismicity with respect to bathymetry and rMBA. **The bathymetry and residual mantle Bouguer anomaly (rMBA) data acquisition and analysis are described in Harmon et al.8. Here, we investigate whether there is preference for seismic events in the local OBS data occurring at sites with particular rMBA/ bathymetry values. In doing so, we compare the rMBA/bathymetry values at the epicentral location of events with the background rMBA/bathymetry data. Background data are defined as all the rMBA/bathymetry values calculated within a 15km wide zone, centred on the transform valley axis and bounded in the along-strike direction by our seismicity data (roughly between the ridge spreading centres). In such way, we consider only the values corresponding to the active deformation zone along the transform valley, where the vast majority of seismicity in the local OBS data occurs, discarding the essentially aseismic areas elsewhere. If the earthquake epicentres show no preference for particular parameter values, then the red and blue histograms in Fig. 3d and Fig. 3e should be similar. We see, however, that there is a clear preference of events occurring at negative rMBA values and an event deficit for positive rMBA (Fig. 3d). There is also a preference of events occurring at very shallow ocean depths (~<3700m, Fig. 3e), i.e., directly beneath the ELFS.

**Magnitude distribution. **We investigate the seismicity magnitude distribution dependence on bathymetry and rMBA, after assigning a bathymetry and rMBA value to each earthquake by nearest neighbour interpolation. Then the events are sorted by bathymetry and rMBA values in an ascending order and the b-values from the corresponding event magnitudes are calculated, together with standard error. Calculations are performed for overlapping 75-event windows, sliding by 1 event (different windows of 50 to 100 events were also tested, as shown in Extended Data Animation 1).

The exponentiality of magnitude distribution is investigated by the Anderson-Darling Test46,47 and further established by goodness of fit test48,49, applying the Maximum Likelihood Estimation (MLE) of b-value50. Both techniques suggest that the magnitude distribution can be sufficiently modelled by an exponential distribution for ML≥2.3. A total of 370 events with ML≥2.3 comprise the complete data for the study area (circles, Fig. 1b), which roughly corresponds to an average rate of ~1 event/day. However, the MLE returns considerably unstable b-values, even well above Mc, fluctuating between 0.76 and 0.89 for 2.3≤ML≤3.0 (Extended Data Table 1, Extended Data Fig. 7). This is incompatible with the unique b-value predicted by the GR law and may lead to artifacts and misinterpretation of the derived outcome. For this reason, we apply alternative methods which are more suitable for small datasets and are less sensitive to MC selection51,52. The comparison of the results from the aforementioned methods, assuming different MC are shown in Extended Data Fig. 7. We choose to estimate b-values by applying the repeated medians53 (RM), a non-parametric regression technique, which provides more robust results for small data sets, since it is highly resistant to observational uncertainties and outliers51. According to the RM approach, we consider n magnitude intervals between two points, i and j, having Ni and Nj events, respectively, and calculate n-1 slopes (b-values) between these points as:

With Mj≠Mi. For each point, i, we calculate the median slopes, thus n median values. Then, the bRM is estimated as the median of the n median values:

bRM=-median(median(bij))

The standard errors, SE, of the corresponding bRM are derived by bootstrap resampling process as described in Amorèse et al.51. The significance of the b-value difference for two datasets, e.g. A and B, can be estimated by a bootstrap t-test, with the t statistic defined as:

A t=1.96 corresponds to a significance at 0.05 level. This approach returns b=0.83±0.09, with this value being insensitive to the applied magnitude cut-off (b=0.82-0.84 for 2.0≤M≤3.0) as shown in Extended Data Fig. 7 and Extended Data Table 1. We therefore apply the RM technique for all b-value calculations in this study, since it provides the most consistent results, regardless the MC choice.

The highest b-values are highest beneath the ELFS (1.1-1.2), they are lowest (<0.8) in the central parts of the fault (~14.7°E – 13.8°E) and intermediate (0.9 – 1.1) in the western area (Extended Data Animation 1). Shallower water depths generally correspond to higher b-values in comparison to the deeper bathymetry areas, where b-values gradually decline with increasing depth at ~0.15 units/km (Fig. 3a). There are higher b-values (b=1.04±0.13) in regions of low rMBA (<-3mgal) than regions of high rMBA (<-3mgal) where b=0.72±0.11 (Fig. 3b, c) with the difference between these two b-values being statistically significant at 0.05 level, indicated by the calculated t statistic equal to 1.96. Events also occur preferentially in areas with negative rMBA (Fig. 3d). These areas are mainly located in the regions of shallowest bathymetry (<3700 m), i.e., beneath the ELFS (Fig. 1; Fig. 3e).

**Seismic Vs Predicted Moment. **We compare the seismic moment release, M0, in each asperity (A1, A2) and barrier (B1, B2, B3) using the formula54 , where M0 is given in N∙m. We express the results in terms of percentage (a) of seismically released moment (SM0s) over the potential moment release assuming full seismic coupling (SM0p), i.e., a= SM0s/SM0p. In doing so, we estimate SM0p by summing the seismic moment of all MW≥5.6 events occurred since 1993 in each patch (barrier or asperity). We then estimate M0p as5:

M0p=GLwsDT

where, G is the shear modulus [40-60 GPa], L is the patch length, w, is the patch width from sea floor to the 600 °C isotherm, s is the slip rate [28-33mm/yr] and DT is the time period duration since 1993, equal to 28 years. Our results range from a=0.25 to a=0.46 with an average value of a=0.32 (for G=50GPa and s=32mm/yr). This range agrees well with previous findings (a=0.33) for the entire Chain Transform7. However, we find the proportion of seismically released moment varies along the length of the transform from 0.02-0.06 in barrier regions to 0.67-0.78 in asperity regions (Fig. 1c; Fig. Extended Data Fig. 2).

**Stress Inversion. **The applied stress inversion algorithm55,56 is used to invert focal mechanisms for a set of earthquakes into tectonic stress, under the assumptions that a uniform regional stress field applies and that the earthquakes occur on faults with varying orientations, but they do not interact with each other. Given that there are only 2 ML>5.0 events, with the stronger earthquake being ML=5.4, the effects of co-seismic stress changes can be neglected57,58 in our study. Under these assumptions, we invert focal mechanisms to estimate the orientation of the principal stress axes (and consequently, the principal focal mechanism) and the shape ratio, R, which determines the relative amplitude of the principal stresses. A typical problem in stress inversion from focal mechanisms concerns the selection of the actual fault plane from the two nodal planes, a choice which may severely bias the results57 mainly the shape ratio determination. To avoid this effect, an algorithm which inverts jointly for stress and fault orientations is considered59. The fault orientation determination is based on the evaluation of fault instability index, I:

where τ and σ are the normalized shear and normal tractions, respectively and m, the friction coefficient. The fault plane is selected as the one having the highest I value, after applying an iterative process. The iteration process performs a linearized inversion55, each time introducing a different noise realization, which randomly rotate the given focal mechanisms. We applied 1000 such realizations in each case tested. The final stress tensor is derived as the mean of the results from all realizations. The uncertainty of the principal stress axes orientations is demonstrated by the scatter of the derived values (Extended Data Fig. 4a), whereas the shape ratio uncertainty is determined as the 95% confidence interval. The instability of principal stress axes shows the impact of the uncertainties from random perturbation of the solution, however it is not clear whether this instability arises exclusively from the spatial stress variations or uncertainties in the focal mechanism determination.

We use the available focal mechanisms along the Chain transform valley to determine the stress field. The vast majority of the reported GCMT solutions45 suggest strike-slip faulting along Chain TF. However, our 1-year seismicity analysis of the OBS data reveals that a proportion of smaller magnitude events (3.0<ML<5.0) demonstrate a considerable vertical component, plausibly related with the positive flower structures (Fig. 1a). We derived 47 best constrained focal mechanisms that are considered for the stress inversion. The principal focal mechanism has strike=275°, dip=62° and rake=165° (Extended Data Fig. 4b). We also derived a shape ratio, R=0.86, with 95% confidence interval [0.69-0.89]. Dividing the region into smaller areas for a spatial analysis would significantly lower the significance of the results due to the decreasing sample size, therefore, we are not able to quantify the spatial variation of the stress along Chain.

The high R value also suggests that the intermediate, s2 and minimum, s3 principal stresses, both deviate considerably from the vertical direction (plunge is 60° for s2 and 30° for s2) and they have similar amplitudes (Extended Data Fig. 4). Therefore, they are virtually indistinguishable and might switch with each other. With s1 being clearly sub-horizontal (2° plunge), this physically means that under such stress field the simultaneous activation of both strike-slip and reverse faults is plausible60. This suggests a non-negligible reverse movement, related to the flower structures, especially within the ELFS, indicating the existence of active transpressional features along Chain.

**S-wave model from local Rayleigh and ambient noise tomography. **We compare our result to the S-wave velocity structure derived from local Rayleigh wave and ambient noise tomography29. Group velocity measurements were inverted from 16-40 s period to determine 2-D velocity maps for each period of interest. Then the group velocities at each pixel across all of the maps were inverted for S-wave velocity structure. Additional details can be found in the original work. The lateral smoothing or correlation length in the group velocity map is 100 km, and checkerboard tests from the study suggest that anomalies at this length scale are well resolved by the method, particularly in the region around the Chain Transform Fault as it is within the centre of the array, with many crossing raypaths with a wide range of azimuths. The formal resolution for the inversion for S-wave velocity with depth is 0.1 for 1 km thick layers in the upper 10 km beneath the seafloor and 20 km at greater depths. This suggests that the S-wave velocity inversion can uniquely resolve the average velocity over the upper 10 km of the model beneath the seafloor, the average velocity from 10-30 km beneath the seafloor and the average value from 30-50 km beneath the seafloor. The lateral and vertical resolution suggests that the variation in lithospheric thickness/velocity we observe across the transform fault at ~100 km scale (WFS1-WF3, A2 vs ELFS) and over 50 km depth is robust. The broader thermal structure of the lithosphere in the region is well within the resolving power of the model presented here. The low P-wave velocities just beneath the Moho of the ELFS from the active source model do not appear in the S-wave velocity model derived from the surface waves, although this is expected given the different resolution of the waveforms. The scale of the crustal and shallow mantle variations observed in our P-wave refraction model is~25-30 km laterally along strike of the transform, 10-20 km across the transform and occurring over ~ 5-10 km in depth. This is below the resolution of the tomographic model used here as described above. In other words, the velocity anomalies in the refraction model would be averaged out in the S-wave velocity model with the relatively faster lithosphere on either side of the fault zone.

The 600°C isotherm derived from the S-wave velocity derived model is closer to 30 km depth across the entire transform, deeper and different in shape than that of simple thermal models (Fig 2). The scale of this variation, 100 km laterally, is well resolved as described above. Isotherms with shapes that deviate from those of simple thermal models can be caused by upwelling in the centre of the transform, as predicted by geodynamic modelling with a visco-elastic-plastic rheology6. While such a rheology cannot be precluded here, we also require an additional factor, given that our isotherms are also deeper than predicted by those models.

**S-wave model translated into thermal model.*** *We examine two methods for translating seismic S-wave velocity to temperature based on the assumption of a peridotitic mantle. In the main text we use the empirical fit to predictions for a pyrolite mantle30:

where Vs is S-wave velocity in km/s, P is pressure in GPa and T is temperature in Kelvin. We also investigate the effects of different composition assumptions on the conversion from S-wave velocity to temperature. We use the Abers and Hacker61 MATLAB toolbox to calculate the predicted S-wave velocity for three end member compositions. We use a depleted upper mantle (harzburgite), an undepleted mantle (lherzolite) and, for reference, an idealized unmelted mantle composition (pyrolite). The specific mineral modes we use are as follows: harzburgite (fo 72.48% volume, fa 7.52 % , en 18.24%, fs 1.76%), lherzolite (fo 45.41%, fa 5.07%, en 22.48%, fs 2.50% di 20.03, sp 4.50%) and pyrolite (alm 2.10%, gr 1.10%, py 10.80%, fo 54.50%, fa 6.70%, en 14.70%, fs 1.50%, di 7.50%, hed 1.2%). We do not account for phase changes of the aluminious phases in the calculations although the spinel – garnet transition is encompassed somewhat between the lherzolite and pyrolite models. We apply a frequency dependent attenuation correction to the S-wave velocity output from the codes using both the approach and the 1-D attenuation structure62. The velocities were calculated for a range of temperatures at a given pressure corresponding to the depth in the model and the temperature corresponding to each velocity in the model was determined by interpolation.

The overall effect of choosing a different starting composition shallows the isotherms predicted from the S-wave velocity model relative to the Stixrude and Bertelloni29 parameterization, particularly the 600°C isotherms, although the effect is not substantial. With the Abers and Hacker61 parameterization the 700°C isotherm lies at nearly the same depth as the 600°C isotherm from the Stixrude and Bertelloni30 parameterization. The difference between the compositions from the Abers and Hacker61 calculations is small, typically the isotherms are within a few km of each other. If the compositions and calculations from Hacker are used all the well-located events are shallower than the 700°C isotherm.

**Seismic tomography of active source airgun profiling. **We used data collected at the easternmost section of the Chain Fracture Zone as part of the ILAB-SPARC experiment63 (Extended Data Fig. 5, Extended Data Fig. 6). Seismic data were acquired in October of 2018 aboard the French research vessel Pourquoi Pas?. A 140 km long seismic profile was acquired at the eastern Ridge-Transform Intersection running from the pressure ridge into the adjacent fracture zone (Extended Data Fig. 5). Five OBS spaced every 12 km recorded airguns shots fired with the Pourquoi Pas?; two sub-arrays containing eight G-guns each provided a total volume of 82 litres and were towed at a depth of 10 m.

Travel times of first arrival P-waves and secondary wide-angle-refection interpreted to be crust/mantle boundary reflection (PmP) arrivals are hand-picked. In general, picking uncertainties were 20-30 ms for short-offset P-waves (Pg) and reach 60 ms for far-offset P-waves (Pn) and secondary PmP arrivals.

Seismic refraction travel time data were used to derive 2-D velocity models using a seismic tomography approach64. The method employs a hybrid ray‐tracing scheme combining the graph method with further refinements utilizing ray bending with the conjugate gradients method. Smoothing and damping constraints regularize the iterative inversion. Details of the procedure are given elsewhere23. Picking errors and starting velocity models may control inversion results. We therefore chose a nonlinear Monte Carlo-type error analysis to derive model uncertainties. The approach consists of randomly perturbing the velocity values of an initial average 1-D model to create a set of 100 2-D reference models, providing a well constrained seismic velocity model (Extended Data Fig. 6).

In the centre of the flower structure the <7.0 km/s contour, corresponding to the Moho, is uplifted. On the eastern end of the flower structure, in the mantle beneath the Moho velocities may be as low as 7.0 – 7.5 km/s, i.e., 4 – 10 % slower in comparison to velocities >7.8 km/s on the western side of the ELFS.

**References**

40. Grevemeyer, I., Rüpke, L.H., Morgan, J.P., Iyer, K, and Devey, C.W. Extensional tectonics and two-stage crustal accretion at oceanic transform faults. *Nature***591**, 402–407, (2021).

41. Maia, M., Sichel, S., *et al.* (2016) Extreme mantle uplift and exhumation along a transpressive transform fault. *Nature Geosci***9, **619–623.

42. Lomax, A., Virieux, J., Volant, P., & Berge-Thierry, C. Probabilistic earthquake location in 3D and layered models. *Advances in Seismic Event Location*, **18**, 101– 134 (2000).

43. Laske, G., Masters, G., Ma, Z. & Pasyanos M. Update on CRUST1.0—a 1-degree global model of Earth’s crust. *Geophys. Res. Abstr*. **15**, 2658 (2013).

44. Heimann, S., Isken, M., Kühn, D., Sudhaus, H., Steinberg, A., Vasyura-Bathke, H., et al. *Grond - A probabilistic earthquake source inversion framework*. V. 1.0. Potsdam: GFZ Data Services (2018).

45. Ekström, G., Nettles, M. & Dziewoński, A. The global CMT project 2004–2010: centroid-moment tensors for 13,017 earthquakes. *Phys. Earth. Planet. Int.***200–201**, 1–9 (2012).

46. Marsaglia, G. & Marsaglia J. Evaluating the Anderson-Darling distribution. *J. Stat. Soft.***9**, 1-5 (2004).

47. Leptokaropoulos, K. Magnitude distribution complexity and variation at The Geysers geothermal field. *Geophys. J. Int. ***222**, 893-906 (2020).

48. Wiemer, S. & Wyss M. Minimum magnitude of completeness in earthquake catalogs: Examples from Alaska, the Western United States, and Japan. *Bull. Seismol. Soc. Am.***90 (4)**, 859–869 (2000).

49. Leptokaropoulos K., Karakostas, V., Papadimitriou, E., Adamaki, A., Tan, O. & İnan S. A Homogeneous Earthquake Catalogue for Western Turkey and Magnitude of Completeness Determination. *Bull. Seismol. Soc. Am.***103****(5)**, 2739-2751 (2013).

50. Aki, K. Maximum likelihood estimate of b in the formula logN = a - bM and its confidence limits. *Bull. Earthq. Res. Inst. Tokyo Univ*. **43**, 237–239 (1965).

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

52. van der Elst, N. J. B-positive: A robust estimator of aftershock magnitude distribution in transiently incomplete catalogs. *J. Geophys. Res.: Solid Earth***126**, e2020JB021027 (2021).

53. Siegel, A. F. Robust regression using repeated medians. *Biometrika***69**, 242–244, (1982).

54. Hanks, T. C. & Kanamori H. A moment magnitude scale. *J. Geophys. Res.***84****(B5)**, 2348–50 (1979).

55. Michael, A. J. Determination of stress from slip data: faults and folds. *J. Geophys. Res.***89**, 11,517-11,526 (1984).

56. Vavryčuk, V. Iterative joint inversion for stress for stress and fault orientations from focal mechanisms. *Geophys. J. Int.***199**, 69-77 (2014).

57. Martínez-Garzón, P., Ben-Zion, Y., Abolfathian, N., Kwiatek, G. & Bohnhoff, M. A refined methodology for stress inversions of earthquake focal mechanisms. *J. Geophys. Res. Solid Earth***121**, 8666–8687 (2016).

58. Fojtíková, L. & Vavryčuk, V. Tectonic stress regime in the 2003-2004 and 2012-2015 earthquake swarms in the Ubaye Valley, French Alps. *Pure Appl. Geophys.***175**, 1997-2008 (2018).

59. Vavryčuk, V. Moment tensor decompositions revisited. *Pure Appl. Geophys.***19**, 231-252 (2015).

60. Hallo, M., Opršal, I., Asano, K. & Gallovič, F. Seismotectonics of the 2018 northern Osaka M6.1 earthquake and its aftershocks: joint movements on strike‑slip and reverse faults in inland Japan. *Earth Planets Space***71:34** (2019).

61. Abers, G. A., & Hacker, B. R. A MATLAB toolbox and Excel workbook for calculating the densities, seismic wave speeds, and major element composition of minerals and rocks at pressure and temperature. *Geochem. Geophys. Geosyst.*,**17**, 616–624 (2016).

62. Saikia, U., Rychert, C., Harmon, N. & Kendall, J.M. Seismic attenuation at the equatorial Mid-Atlantic Ridge constrained by local Rayleigh wave analysis from the PI-LAB experiment. *Geochem. Geophys. Geosyst.*, **22**, e2021GC010085. (2021).

63. Marjanović, M., Singh, S. C., Gregory, E. P. M., Grevemeyer, I., Growe, K., Wang, Z., et al. Seismic crustal structure and morphotectonic features associated with the Chain Fracture Zone and their role in the evolution of the equatorial Atlantic region. *J. Geophys. Res. Solid Earth*, **125**, e2020JB020275 (2020).

64. Korenaga, J., Holbrook, W. S., Kent, G. M., Kelemen, P. B., Detrick, R. S., Larsen, H. C., et al. Crustal structure of the southeast Greenland margin from joint refraction and reflection seismic tomography. *J. Geophys. Res.: Solid Earth*, **105**, 21591-21614 (2000).