Left panel: Surrogate data (FT-method) made from the raw geopotential height data series at the 700 hPa level averaged over 70°–90° S cross-correlated with the real IMF B_y index over the time period 1968–2020 for 20,000 iterations. The result of these iterations makes up a distribution of correlation values for every lead lag. The red area illustrates where only 5% of correlation values land, while the green shaded area illustrates where the remaining 95% of correlation values land. The black line is added as an example and shows a typical iteration when the surrogate data and B_y are cross-correlated. Right panel: After constructing the distribution of correlation values shown in the left panel, the surrogate data and B_y are again cross-correlated 20,000 times. Then, all data points at specific lead lags are measured against the distribution (in left panel) at the specific lead lag to obtain the p-values. To accurately estimate the appropriate α_FDR in a way that takes into account the autocorrelation present in our data, all p-values generated from each iteration of the simulation are processed through the FDR method, where 5 different values for α_FDR are tested. By doing so, we can obtain the specific α_FDR value which ensures that any signal determined to be statistically significant occurs globally in only 5% of cases when the null hypothesis is assumed to be true for our data. As can be seen, when α_FDR is set to 0.09, only 5% of the 20,000 iterations produce a response that passes the global FDR limit.