Can the boundary profiles at 26N be used to extract buoyancy-forced AMOC signals?

The AMOC circulation is driven both by direct wind stresses and by the buoyancy-driven formation of North Atlantic Deep Water over the Labrador and Nordic Seas. In many models low frequency density variability down the western boundary of the Atlantic basin is linked to changes in the buoyancy forcing over the Atlantic Sub-Polar Gyre (SPG) region, and this is found to explain part of the geostrophic AMOC variability at 26N. In this study, using different experiments with 10 an OGCM, we develop statistical methods to identify characteristic vertical density profiles at 26N at the western and eastern boundaries which relate to the buoyancy-forced AMOC. We show that density anomalies due to anomalous buoyancy forcing over the SPG propagate equatorward along the western Atlantic boundary, through 26N, and then eastward along the equator, and poleward up the eastern Atlantic boundary. The timing of the density anomalies appearing at the eastern and western boundaries at 26N reveals a propagation speed leading to ~2-3 years lags between boundaries with maxima along deeper levels 15 (2600-3000m). Time record required to capture those vertical density profiles in the model is ~26 years. Results suggest that depth structure, and the lagged covariances between the boundaries at 26N, may both provide useful information for detecting density anomalies of high latitude origin in more complex models, and potentially in the observational RAPID array. However, time filtering will be required together with the continuation of the RAPID program in order to extend the time period.


Coupled experiment
We also analyse 120 years of monthly-mean data from a control run of the high-resolution coupled ocean-atmosphere model 125 HadGEM3-GC2 (hereafter GC2, Williams et al 2015). (~60Km in mid-latitudes) and 85 levels. This model was used in the Met Office seasonal and decadal prediction systems (GloSea5 and DePreSys3 respectively). The model has been used to study the North Atlantic variability and its predictability 130 (Menary et al 2015;Williams et al 2015;Ortega et al 2017;Robson et al., 2016). Results are shown in Section 6.

Model evaluation
We use the RAPID array (McCarthy et al., 2015;Smeed et al 2017) to evaluate the boundary densities in the model. We use the merged profiles at the Western Boundary (26.52N, 76.74W) and Eastern Boundary (26.99N, 16.23W) for the period April 2004 to February 2017. 135 NEMO1 and GC2 are both able to capture important aspects of the observed boundary density profiles such as the mean vertical density gradients (N 2 , Fig. S1a). On the WB the profiles are similar between 1500m to 4500m but the model stratification is stronger between 300-700m. The EB profiles are similar at all levels below 500m (Fig. S1a). However, the NEMO1 model underestimates the density variance at all levels, especially at the WB, while GC2 has a more realistic variance 140 on the WB at depth (Fig. S1b).
The AMOC at 26N in NEMO1 has a time mean (12Sv) and maximum (18 Sv) at a depth of 1000m in the CTRL experiment.
The mean AMOC is higher in the BUOY experiment by ~2 Sv. The AMOC measured at 1000m has a prominent trend in the BUOY experiment (+3.2 Sv in 52 years), but the trend is not significant for the CTRL (-0.2 Sv in 52 years). The AMOC 145 seasonal cycle (not shown) in the CTRL presents a maximum in boreal winter with a secondary peak in boreal summer, which is also reproduced in the BUOY experiment. The annual cycle defined as the difference between the maximum (in boreal winter) and the minimum (in boreal spring) is 3.9 and 4.9 Sv for the CTRL and BUOY experiments. After removing the linear trend and the seasonal cycle the standard deviation of both experiments is similar: 1.89 and 1.43 Sv for CTRL and BUOY respectively (see also Fig. S2 for the AMOC distributions). The AMOC at 26N in the RAPID observations presents a mean 150 and maximum of 17 and 31 Sv respectively from 2004-2014, with monthly standard deviation of 4.35 Sv. Trends have been reported for RAPID data of -0.6Sv/year (Smeed et al 2014;2018) which could be part of a longer variation cycle (Smeed et al 2018, Jackson et al 2016. Results from the modes of variability at the western boundary density profile are shown in section 5. https://doi.org/10.5194/os-2020-8 Preprint. Discussion started: 21 February 2020 c Author(s) 2020. CC BY 4.0 License.

Statistical analysis 155
Model experiments are first sampled at the western and eastern Atlantic boundaries at 26N to simulate the monthly-mean density profiles from the RAPID array. Empirical Orthogonal Function (EOF) analysis of these density profiles is used to obtain vertical modes of density variability and related timeseries of Principal Components (PCs) which together represent the largest fractions of the total variance (Bretherton et al., 1992).
Before calculating the EOFs, the data are processed to remove the seasonal cycle and linear trends. Unlike in PA14, density anomalies are weighted by the thickness of each layer to ensure that all points are appropriately represented for the total density variability. The EOF analysis is computed for the individual boundaries; Western Boundary (WB) and Eastern Boundary (EB), and also for the combined anomalies at both boundaries. Finally, we have explored the combined EOFs by time lagging the eastern boundary variability in order to understand related signals from both boundaries. 165 Regression analysis of the PC time series associated with these EOFs on other fields (e.g. 3D density) allows us to detect spatial patterns and depth structures of the propagating modes associated with the EOFs at 26N. We show regression and the correlation coefficients where they are statistically significant at the 95% confidence level, according to a Student's t-test for the effective number of degrees of freedom (Metz, 1991). 170

Spectral analysis
In order to remove the high frequencies in the timeseries in section 3, we have used a one-year running mean filter. This filtered timeseries is obtained by taking the average of a data subset (13 months) which is centered in a monthly time step (von Storch and Zwiers, 1999).

175
Spectral analysis is used to decompose timeseries to show signals that lie within different frequency bands. The analysis is performed in order to identify the frequencies involved in the propagation of density anomalies at different depths. Power spectra of the time series are obtained using the multi-taper method, which provides more degrees of freedom and therefore more significance (Thomson 1982). The power spectra are tested against the hypothesis that the signals are generated by a first-order autoregressive process AR(1) with the same time-scale as the original, yielding a red noise spectrum, and the 95% 180 confidence limit for the rejection of the red noise hypothesis is applied. Additionally, when we have a large internal variability, we use the decomposition in order to filter some of the timeseries using a Lanczos (1956) filter. This is done in particular for the control GC2 run simulation.
We have used the Radon Transform (RT) function (Dean 1983) in order to estimate the phase speed of propagation of density 185 anomalies. The angle of the maximum RT standard deviation determines the propagation phase speed. We calculate the RT https://doi.org/10.5194/os-2020-8 Preprint. i.e. from the BUOY experiment, as noted by PA14. We focus here on developing the vertical density fingerprint of this signal and using it to identify the buoyancy forced AMOC signal as it appears in the CTRL run. Figure 1a shows the monthly AMOC (defined as total AMOC minus Ekman component) variability, defined as the integral of 195 the meridional transport at 26N down to 1028m, for both the CTRL and BUOY experiments (as in Fig. 1a in PA14). There is a prominent decadal signal in both CTRL and BUOY with peaks in 1975, 1985and 1995 although CTRL also shows additional monthly and inter-annual variability. The monthly-mean timeseries correlate at 0.43, but the correlation rises to 0.62 when using a 1 year running mean filter (Fig. 1b). Wind forced inter-annual variability explains most of the remaining differences; when the 1-year smoothed AMOC anomalies from BUOY and WIND are summed the correlation with CTRL rises to 0.86 200 (SUM in Fig 1b).
The majority of the boundary density variability is also recreated in BUOY and WIND. The correlations between boundary density anomalies in the CTRL and SUM are shown as a function of depth in Fig. 1c-d. For the WB, most of the variability is linearly reproduced by SUM from 1800m to 4000m ( Fig. 1c and Fig. S3). For the EB SUM explains most of the variability 205 seen in the CTRL experiment at all depths ( Fig. 1d and Fig. S3). A 1yr low pass filtering does not influence the correlations for the WB, although for the EB filtering reduces the correlation at some depths. We now relate this density variability with the AMOC signals in Fig 1a, b. Figure 2 shows the Principal Component time series of the first EOF computed using monthly density profiles on the western 210 and eastern boundary at 26N for the CTRL (blue) and BUOY (red) experiments. As density profiles near the surface contain significant noise, we calculate the PCs for both full depth profiles (f0m, Fig. 2a-b) and only from 800m downwards (f800m, a 1-yr filter to the PC-CTRL timeseries increases the correlations with the CTRL AMOC variability to 0.47 (Fig. 2c). However, for the BUOY experiment the temporal filter has little impact on the correlation between PC and AMOC timeseries (0.89).

EOFs of boundary density profiles
Figures 2b and d show the corresponding eastern boundary density EOF timeseries for full depth and below 800m variability. 220 The full depth PC is rather different to that on the western boundary, and also to the AMOC variability itself. However, when only the deep variability is retained, the inter-annual variability in BUOY is more similar to the western boundary buoyancy driven variability. In CTRL there is still considerable high frequency wind generated variability, however, when a 1-yr filter is introduced the buoyancy forced AMOC related signal also becomes visible on the eastern boundary, and the correlation with the AMOC rises to 0.42 (Fig. 2d). All correlations are summarised in Table 1. 225 In order to understand the gain of truncating the density profiles in the EOFs and theirs limits, Fig. 3a-b summarises correlations between leading PCs at the boundaries and AMOC timeseries in the CTRL experiment by increasing the truncation level for the density profiles. For the WB, maximum correlation is found simultaneously at all truncations especially for deeper levels ( Fig. 3a). In contrast for the EB simultaneous correlations are always low even with a 1 year filter but correlations increase 230 greatly when a lag to the AMOC is applied (Fig. 3b). Correlations still require a 1 year filter to remove noise but now peak at 0.7 with a lag around 2-3 years, with the deeper signals also showing the longer lags. Similar lag increases of EB correlation is seen for the BUOY experiment (not shown). The nature of this lag in the EB-AMOC correlation is related to the link between boundaries and it will be discussed along in the next section.
235 Figure 3c-d shows the vertical structure of density anomalies associated with the leading EOF modes for both boundaries using the 800m depth truncation and monthly-mean data for the three experiments CTRL, BUOY and WIND. The leading EOFs are very similar between CTRL and BUOY on the WB, with maxima between 1200m and 4000m. Note that in PA14 Fig. 3b, their PC with only 500m depth truncation was substantially different and was mainly wind-driven. The WIND experiment has much less variability at greater depths and the PC is uncorrelated with the PC-WB from CTRL (r=0.09, Fig. 3c). Therefore on the 240 WB, Fig.s 2c and 3c show that the EOF analysis successfully extracts the buoyancy-forced signal related to the AMOC in the CTRL.
On the EB the leading EOFs, even below 800m, show more correspondence between CTRL and WIND ( Fig. 3d, with PC correlations r=0.87). This explains why further filtering and the application of a lag to the PC timeseries in Fig 2d, 3b is needed 245 to extract the weaker buoyancy-forced AMOC-related signal on the EB. The relationship between these buoyancy driven density variations at both boundaries is now explored further. Figure 4a shows lead-lagged correlations between the WB leading PCs (computed from 800m depth) and the EB at different depth truncations in the CTRL experiment. Dashed (solid) lines indicate PCs without (with) the application of 1year running 250 mean filtering. For truncations deeper than 1600m the highest correlation is found when EB is lagged, up to 30 months ( Fig.   4a), revealing longer links between boundaries for deeper levels. However this lag is much less clear when shallower depth are retained and substantial WB-EB density correlations can still be seen with lag0, unlike in Fig 3b. Even in the BUOY experiment only when the EB EOFs are truncated to below 1600m is a strong lag clearly seen between the boundaries, again reaching up to 30 months for the deepest signals (Fig. 4b). 255 In order to reduce upper level noise at the EB and to bring out the deep density signal connecting the boundaries more clearly we compute the combined EOF while truncating the EB to below 1600m, which shows a maximum in the WB-EB correlations for CTRL (Fig. 4a). The new combined EOF shows similar WB structure for all lags (Fig. 4c) with a deeper signal around 2000m on the EB (Fig. 4d). Figure 4e shows the timeseries, which is very similar to the PC1-WB in Fig. 2c. The new combined 260 EOF explains more variance in CTRL (43% compared to 40%) and a slighly longer lag between WB and EB (18 months in CTRL and 25 months in BUOY, not shown).

Relationship between boundaries
We conclude that the deep densities on the EB contain a very clear signal of the buoyancy forced AMOC variability but that this signal plays no detectable role in the direct (lag 0) control of the AMOC. The EB lag can also be seen in relation to PC-265 WB density although the signal is less clearly lagged, probably reflecting the noise still present in the upper layer densities on the WB.
The WB clearly contains the core density information on the buoyancy-forced AMOC changes at low frequencies. Using PC1-WB, we will now identify the propagating density signal connecting the boundaries (section 4) and search for similar signals 270 in RAPID data (section 5), and in the higher resolution GC2 model (section 6).

Propagation of the buoyancy-driven signals
Motivated by the lagged signal at the EB, we analyse the i) spatial coherence of the anomalies at deeper levels (~3000m), and ii) the propagation fingerprints of the connecting signals. BUOY, we are identifying the same signal in both experiments. The density signal in the Labrador Sea appears from lag -30 ( Fig. 5a, d) and intensifies and propagates down the WB to the Equator (Fig. 5b, e), and then across the equator and poleward 280 at the EB (Fig. 5c,f).

Spatial regression patterns 275
Density regressions at shallower levels (900-1300m) against PC1-WB in the CTRL experiment (suggested by the maximum in the EOF profiles in Fig. 3) also finds anomalies beginning near the Labrador Sea, leading PC1-WB by 30 months, and a pattern of equatorial Kelvin and Rossby-waves as in Johnson and Marshall (2002) from lag -30 to lag ~0 (Fig. S5). However, 285 the absence of such tropical signals in BUOY at long lags (lag -30) shows that these could be wind driven Ekman pumping signals (see supplementary material in Fig. S6). Therefore, we concentrate on the deeper signal which is clearly related to buoyancy. Figure 6a shows the wave-track defined following the signals in Fig. 5 and using the topography at the 3269 m model level. 290

Wave path and phase speed
The wave-track starts in the Labrador Sea (60N) and proceeds southwards to the equator. We plot it along the equator and then North along the eastern boundary to 55N. The path avoids entering into the Gulf of Mexico but at these depths this is not expected.
The Hovmöller diagram along this path (Fig. 6b) shows the propagation related to peaks in the deep water formation at high 295 latitudes for the CTRL experiment (BUOY is very similar, not shown). Density anomalies propagate continuously along the track from the Labrador Sea around to the British Isles. The propagation shows density maxima in 1975, 1985 and 1995, also seen as peak-AMOC years in Fig.s 1-2. Additionally, the propagation speed from the Radon transform is similar for both experiments (0.41 and 0.31m/s for the BUOY and CTRL respectively). This phase speed is consistent with the lags found between boundaries (i.e. density anomalies in CTRL will take ~25 months to travel between WB and EB following the defined 300 track). Background currents at 2700-3000m level for the CTRL is shown in supplementary material (Fig. S4). Importantly these propagating buoyancy-related signals are clearly seen in the CTRL experiment where wind and buoyancy forcing are both applied, again suggesting that the analysis is extracting the same buoyancy forced processes. Therefore, the diagnostic methods developed should allow identification of similar signals in other models and in the observations. In the next section we look for similar buoyancy-forced signals in the RAPID observations.

Modes of WB density variability in RAPID data
The RAPID timeseries dataset is considerably shorter than the datasets analysed for buoyancy signals in the NEMO model. Nevertheless Fig. 8 uses the same EOF analysis on the WB density profile for the RAPID array and for the NEMO experiments for the common period 2004-2010. The first two leading EOFs and their PCs indicate inter-annual variability in RAPID data (black lines in Fig. 8a, c). EOF1 has a maximum density anomaly at ~1000m (with a maximum of 0.04 kg/m3 in Fig. 8b), 315 while the second mode describes variations at deeper levels (~3000m, Fig. 8d).
The density EOFs of the WB in both the CTRL (blue lines) and WIND (green lines) experiments for the common period 2004-2010 look quite similar to the EOFs from the RAPID array, in particular with EOF1 peaking at around 1000m and EOF2 peaking much deeper at ~3000m (Fig. 8b, d). However in BUOY (red lines) the absence of inter-annual wind variability 320 means that the deep 3000m peak shows up as EOF1 which is more comparable to EOF2 from RAPID, while BUOY EOF2 shows very little clear density signal.
The PC1 time series in CTRL and WIND also look very similar to PC1 from RAPID, reaching a peak in 2007 and declining to 2010 (Fig. 8a). For these common 6 years of simulation 2004-2010, the wind-forced inter-annual density variability on the 325 WB is remarkably well captured by the model. Figure 9 shows the statistics of the correlation between the leading WB mode for CTRL and both, WIND (Fig. 9a) and BUOY ( Fig. 9b) experiments using different numbers of years for each sub-period (x-axis). Figure 9a (first box) shows the test of sampling 6 years periods (as we have in 6 years in the common period in Fig. 8a). For all 6 year periods selected between 1958 330 and 2010 we also found that EOF1 and PC1 in CTRL and WIND agreed well (Fig. 9a, first box), even in periods when the buoyancy signal was known to be changing rapidly. This dominance of the wind forcing over short time periods is not surprising and was noted previously (PA14).
We find that typically 16 years of data are needed to find a significant (above 0.35) correlation between CTRL and BUOY 335 (Fig. 9b), however the leading mode in CTRL may still be a mix of wind and buoyancy forcing (as seen from Fig. 9a). The best extraction of buoyancy forcing signals occurs when we have periods longer than 35 years, when the wind-forced signal nearly disappears (Fig. 9a). Figure 9c shows the correlation between the PC1-WB and the AMOC timeseries for the CTRL and different periods. This would represent how much variance of AMOC can be explained by the PC1-WB. For periods with more than 25 years, PC1-WB in CTRL experiment is able to explain more than 25% of the AMOC variance (r>0.5), therefore, 340 it is still able to extract an AMOC-related signal at 26N (Fig. 9c). Although EOF2 from CTRL and RAPID both represent deeper density variability (Fig. 8d) they are still dominated by wind forcing over this short time period, and we note that EOF2 from WIND shows the same deep density peak. The PC2s from CTRL and WIND show very similar time series and even the RAPID PC2 shows a considerable level of agreement with CTRL 345 and WIND (Fig. 8c). However PC1 from BUOY, which has a similar EOF but represents only the buoyancy forced component, has only lower-frequency changes with no relationship to the other timeseries (Fig. 8a). The BUOY PC2 timeseries (Fig. 8c) also shows no comparable variability. Therefore, the short record of the RAPID array would not allow us to follow buoyancy-forced signal from the NEMO1 model. 350 Hence, it appears likely that the variability seen here in the RAPID record, both at shallow and deeper depths, is mainly related to wind-forcing (both PC1 and PC2 in CTRL experiment and RAPID show agreement in Fig. 8). Note that the same EOF analysis using the longer record now available for RAPID data (2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017)(2018), but not for these model results, still gives similar density profiles to Fig. 8 (not shown).

355
The NEMO1 model results suggest that the lower frequency buoyancy forced signals from higher latitudes may start to dominate over the wind forced signals after ~25 years of RAPID data have been collected, when their leading density variability should show up at deeper depths ~3000m. In the next section we look at the ability of the analysis to extract buoyancy forced signals and their propagation in higher resolution HadGEM3-GC2 coupled model run, which is the current UK operational coupled model. 360 6 Boundary density in a high resolution coupled model Figure 10a shows a wave track for the GC2 model (section 2.2) bathymetry using the 3138m level boundary, and figure 10b shows the Hovmöller diagram of the density anomalies for 120 years along this track , after applying a 2 year low pass filter as GC2 is much noisier at higher frequencies (see methods). We notice that the high frequency density anomalies (<2years 365 period) in a similar Hovmöller are not presenting propagating signals but very noisy signals (not shown), therefore we suggest that the high frequency signal is dominated by local Ekman pumping.

Density propagation in GC2
Anomalies propagate down the WB from ~40N (track point 200) to the equator and across to the EB. These signals can be traced back across the Gulf Stream to subpolar latitudes (points 100-200), but only appear as lower frequency decadal 370 variations in the subpolar gyre and into the Labrador Sea (points 0-100). The Radon Transform phase speed is ~2m/s, which is faster than the phase speed calculated in NEMO1 at the same depth (Fig. 6), and closer to the theoretical and observed Kelvin wave propagation speed (Polo et al 2008). Therefore, as the density anomalies propagate down the deep western boundary, we would expect to find this deep density variability signal using the EOF analysis. https://doi.org/10.5194/os-2020-8 Preprint. Discussion started: 21 February 2020 c Author(s) 2020. CC BY 4.0 License. Figure 11a shows the vertical density profiles associated with the first two EOFs from 800m downwards at the 26N WB in the GC2 experiment. The profile location is the closest grid-point to the western wall at the Bahamas, which reaches the bottom at 3200m, and the first 2 EOFs explain more than 70% of the total variance. The EOF1 shows an equivalent barotropic vertical structure peaking near the bottom ~3000m (blue line) while EOF2 changes sign between 900m and 3000m (red line). 380

EOFs at the WB in GC2
The PC1 timeseries associated with EOF1 is plotted in Fig. 11b (blue line). Unlike in NEMO1 this PC1 shows high frequency variability which is nevertheless still correlated with the AMOC-Ekman (i.e. the peak AMOC stream function at 900m after the variability due to Ekman has been removed) without filtering at r=0.45, rising to r=0.49 with high-pass (<2 years) filtering.
However, PC1 becomes less correlated with the AMOC after 2 year low-pass filtering, r=0.25 (Table 2). In contrast, for PC2 385 the unfiltered correlation with the AMOC-Ekman is low (r=0.13) but this increases with 2 years low pass filtering (r=0.32, Fig. 11c, Table 2). As density anomalies at deep levels are able to be excited by wind alone (already seen in both NEMO1 and the RAPID observations in Fig. 8), the two EOFs in GC2 both capture some wind and buoyancy forcing, with variability signals in both PC1 and PC2.

390
Although the PC1 and PC2 are orthogonal by construction and thus correlation between timeseries is zero, after time filtering the PCs timeseries, the modes are correlated (r= 0.41) and the lead-lag correlations with AMOC present a cycle between vertical profiles modes, with a 16 months between the peaks (Fig. 11d). This indicates the limitations of extracting low frequency AMOC-related signals in complex environments using linear methods.

395
High frequency PC1-WB (<2years) represents high frequency density signal at 26N, which could be wind-forced. It is correlated with the AMOC and is independent of the number of years used to identify it (Fig. 11e). Low frequency PC1-WB (>2years) represent a small part of the low frequency density signal at 26N (r<0.2) and it is independent of the number of years (not shown). PC2-WB correlates better with a lower-frequency AMOC signal (>2years) and is also independent of the number of years used to identify it (Fig. 11f). 400 If we filter the density anomalies prior to performing the EOF, then the leading mode corresponds to PC2-WB seen here. This confirms that we cannot isolate the low frequencies by identifying a deep density signature as works well in NEMO1, therefore time-filtering is needed to identify the inter-annual buoyancy signal in GC2 as a leading mode. Nevertheless, the buoyancy signal is still traceable emerging from the Labrador Sea, and is well captured in PCs, representing relevant information on 405

Discussion
In this work we have used model output and statistical methods to identify vertical density profiles along the boundaries that are consistent with the buoyancy-forced variability. We have shown that the most relevant profile at 26N is found at the WB 440 using EOF analysis after truncating the density profile from 800m (PC1-WB). This truncation is very effective in emphasising low frequency (decadal time-scales) signals and in NEMO1 negates the need for temporal filtering, which can also add spurious signals or lead to excess smoothing. Caveats that warrant further discussion include the differences between the EOFs in NEMO1 and GC2, and the role of the EB.

445
We note that vertical profile of the WB EOFs in GC2 and NEMO1 are different, especially in the top 1500m, and the frequencies of the dominant variability in boundary signals in GC2 are higher than the decadal signature seen in NEMO1. We noticed that shallow (<1500m) density signal related to PC1-WB in NEMO have different timing in CTRL and BUOY experiments at the tropics (Fig. S5-6), suggesting that at 26N wind forcing is modulating the buoyancy-forced density signal in CTRL. This is also an argument to suggest that in GC2 the shallow signal is more probably wind-forced signal. 450 Although PCs GC2 contains the low-frequency AMOC-related variance it is perhaps not surprising that the details of the WB EOFs are different between GC2 and NEMO1 given the range of AMOC variability in models (Biastoch et al., 2008;Cabanes et al., 2008;PA14;Ortega et al., 2017). The 1° horizontal resolution, and even the ¼° model, may still be too coarse to correctly capture propagating boundary signals (Johnson and Marshall, 2002;Getzlaff et al., 2005;Hodson and Sutton, 2012) from the 455 Labrador Sea. Therefore, we may not expect the exact details of the boundary density EOFs, or the phase speeds nor the phase lags identified from the boundary and Labrador Sea signals to be very realistic. Nevertheless, the methods reveal, in two very different modelling environments, boundary signals consistently related to the geostrophic AMOC at 26N simultaneously.
It is worth notice that time filtering is needed to see a clear signal due to noise from high frequency wind-forcing in GC2. 460 Wind forcing can be projected onto density anomalies at deep levels (as it is seen in RAPID data and GC2), it is possible that in the observations similar time-filtering would be needed in order to extract the buoyancy-forced signals. Understanding the differences between these boundary density EOFs between models, would be useful to interpret observations as the record becomes longer, and it should be a focus of further work.

465
In Section 3.3, we tested the value of using WB and EB together in a single EOF. The fact that the combining the WB and EB at different lags in a singular structure does not improve the explained variance of the EOF WB alone reveals that the WB contains most of the variability at decadal time-scales. This is in agreement with previous results by PA14 and also with the propagating signal from the Labrador Sea that can reach the subtropics along the WB (as in shown by Jackson et al., 2016), but it is not that obvious for the EB. However, we must be clear that we are not saying that EB observations do not make up 470 an important component of the RAPID array observations. Indeed, it has already been shown in observations (Kanzow et al, 2010;Duchez et al 2011), and in models (PA14), that the EB is important for understanding the wind-forced variability in the observed AMOC at 26N from sub-annual to inter-annual time-scales. Indeed, in our own study the use of the EB was important to isolate the density propagation at 3000m depth which is an AMOC-related signature. Therefore, the EB observations are still important in understanding the role of decadal-buoyancy forced variability. 475 We have found density variability in RAPID which the models suggest could be buoyancy-forced at high latitudes (PC1-WB in NEMO1). However, the temporal variability suggests that wind forcing is still dominating at these short timescales. We have to wait for more years of RAPID data, which could allow the buoyancy forced variance to dominate this mode or give this mode EOF1 status in the decomposition. Although it may be helpful to do some time filtering which could allow buoyancy 480 forced signal to be found in shorter periods of data, the continuation of the RAPID array would be crucial in order to understand the wind-forced inter-annual variability and also the link between the subpolar North Atlantic and the AMOC at 26N.

Conclusions
In this work we have used NEMO1 OGCM experiments which separate buoyancy and wind forced signals (BUOY and WIND experiments; Polo et al., 2014), together with statistical techniques, to develop methods to extract the Atlantic boundary density 485 profile signatures at 26N most associated with the buoyancy-forced AMOC signals from an experiment with both buoyancy and wind forced variability (CTRL). After finding the "best" vertical profile on the western boundary, we describe the temporal-spatial structures related to this signal. The main findings are summarized as follows: • Using EOF analysis and outputs from OGCM experiments we find that the vertical density structure at both the 490 western (WB) and eastern boundaries (EB) at 26N show characteristic signatures that can be unambiguously linked to buoyancy-forcing in the Labrador Sea.
• The vertical structure associated with the leading EOF mode of density variability on the WB (EOF1-WB) shows positive anomalies from 1500 to 3000m depth that can be related to earlier changes in the North Atlantic deep water 495 formation, and to density anomalies over the Labrador Sea, which are seen to lead PC1-WB by ~30 months. The PC1-WB is found to be very robust in both the CTRL and BUOY experiments, signalling buoyancy-forced AMOC variability on decadal timescales. PC1-WB explains 40% and 70% of the density variance for the CTRL and BUOY experiments respectively.

500
• The PC1-WB is found to lead density anomalies at 1000-1500m on the EB (associated with PC1-EB) by ~7months.
The result of combining both boundaries into a single EOF allows to extract the correlated variance and the optimal https://doi.org/10.5194/os-2020-8 Preprint. Discussion started: 21 February 2020 c Author(s) 2020. CC BY 4.0 License. lag between the boundaries. The combined EOF variance shows maxima when lagging the EB between 7 months and 3 years. This lagged relationship is consistent with density propagation at 2700-3000m.

505
• In the CTRL experiment, density anomalies at 2000-3000m propagate southwards along the WB and eastward along the Equator and then up to the African coast impacting the vertical structure on both boundaries at 26N. The propagation is continuous from the Labrador Sea around the basin and up to the British Isles. This density signal propagates at a speed ~0.3m/s consistent with the propagation speeds in the BUOY experiment.

510
• The same method is applied to the RAPID array data for the common period with NEMO1 simulation. The two leading EOFs for the WB have anomalous densities at 1000m and 3000m respectively and are well simulated by the CTRL experiment (Fig. 8). This inter-annual variability is unequivocally wind-forced. The observational record has to be longer in order to identify the buoyancy-forced vertical anomalies, which are more low-frequency signals. 515 • The same method was able to extract boundary signals from the higher resolution model HadGEM3-GC2. Despite the greater complexity in GC2, the vertical density profiles on the WB at 26N can be clearly related to the geostrophic AMOC although some time-filtering is needed in order to separate the time-scales. 520 • After filtering (periods below and above of 2 years), PC1-WB is found to be more related to AMOC at 1000m (2000-3000m) at high (low) frequency with a EOF profile with positive density anomalies at 1000m and 3000m.
PC2-WB is also related to AMOC at 1000m at low frequency and shows positive (negative) density anomalies at 525 deep 1000m (3000m) level.
• We also show clear density propagation from the Labrador Sea around the basin to the British Isles along a wave track (defined by 3138m bathymetry) at 3000m depth level, which it is also explaining part of the AMOC variability. However, temporal filtering is needed to make this stand out above the noise. 530 • We conclude that the buoyancy-forced signal over the density profile at 26N will be captured in the observations (as well as in coupled models) if the available time period is long enough (>26 years), selecting density profiles with opposite anomalies at 1000m and 3000-3500m, with time filtering of periodicity >2 years, which would help to eliminate high frequency wind-driven signals.