Subtropical Western Boundary Currents (WBCs) are narrow, rapidly flowing warm water currents that are important for ecosystems, climate, weather and cross-shelf exchange (e.g. Lambaerts et al., 2013; Malan et al., 2020; Oliver et al., 2021; Li et al., 2022a). As fast-flowing WBC jets become unstable, they shed O(100) km-wide mesoscale eddies. The formation, structure and evolution of these eddies and associated structures are important due to the impact they have on the transport of heat and salt (Abernathey and Haller, 2018), weather (Frenger et al., 2013), mixing (Klocker and Abernathey, 2014), and the delivery of nutrients (Everett et al., 2012).

Given their location adjacent to populous coastlines, WBCs have a pivotal role in coastal fisheries and other blue economies (e.g. Young et al., 2011; Brieva et al., 2015), weather and climate, and search-and-rescue and navigation. An impediment to manyof these end-users is the limitation in model predictability resulting from the dynamically changing eddy field. Mesoscale dynamics are inherently sensitive, where divergent evolution results from small differences in initial conditions. This leads to a timescale limit on predictability for techniques such as search and rescue, navigation and other methods that require accurate forecasts of ocean weather.

Accurate estimates and predictions of WBCs and eddy-rich regions are generally sought through data-assimilating models (Oke et al., 2013). The technique of data assimilation combines observations with a model forecast, often in an iterative process, to produce an analysis or optimal estimate of the ocean state. Hence, the evolution of eddies can be continually updated within model forecasts, providing a best estimate of eddy structure, timing and location (Oke et al., 2013).

One of the key ways in which operational forecasting systems differ from non-operational, research-focussed assimilation systems is the types of observations that can be assimilated. For example, operational forecast systems typically assimilate sea surface height (SSH), sea surface temperature (SST), a smaller number of subsurface observations, such as vertical temperature profiles from expendable bathythermograph (XBT) probes and Argo floats, as well as occasionally wind stress and sea surface salinity observations (e.g. Brassington et al., 2007).

In contrast, hindcast reanalyses with a research focus (e.g. Kerry et al., 2016; Siripatana et al., 2020), can augment these traditional observation types with more novel, often delayed-mode observations, e.g. high frequency radar-inferred surface currents, hydrography from autonomous gliders, and measurements from subsurface moorings. The constraint limiting operational forecast systems from assimilating the full gamut of non-traditional and subsurface observations, is the ability to have data prepared and available in near-real time for the next operational window — which can often not be met when quality control or other data preparation must be conducted with human supervision, or obviously if there is a delay in data recovery. As a result, operational systems will often be limited to just surface observations combined with a small number of near-real time subsurface measurements.

The East Australian Current (EAC; see Figure 1 for region) is one such WBC with routine subsurface sampling. The EAC flows southwards from 27.5°S as a coherent jet, before beginning to meander and lose coherency between 31°–33°S and then feeding an ‘evolving’ field of cyclonic and anticyclonic eddies in the Tasman sea. Like many other WBCs, an important source of subsurface real-time measurements in the EAC are repeated XBT transects and Argo floats. A long-term program of XBT deployments from ships of opportunity has been operated along repeat transects in the southern Pacific Ocean since the late 1980s (and 1991 for the two transects through the EAC region, named ‘PX30’ and ‘PX34’) by Scripps Institution of Oceanography, the Australian Commonwealth Scientific and Industrial Research Organisation and the New Zealand National Institute for Water and Atmospheric Research. The original use of this data was for estimating boundary current heat budgets (e.g. Roemmich and Cornuelle, 1990; Morris et al., 1996; Roemmich et al., 2005). However the near-real time data delivery and consistent transects lends itself well to data assimilation into ocean models. As the XBT data is delivered to the global telecommunications system it is readily available in near realtime for ingestion into operational modelling systems and reanalyses (e.g. Carton et al., 2000).

Argo floats are a second source of subsurface observations that have also been used in operational forecasting. Argo floats return vast amounts of deep (to 2000 m) vertical temperature and salinity profiles over a much broader area of the ocean and have revolutionised understanding of the ocean (Wijffels et al., 2016; Wong et al., 2020). However, they still have relatively low spatial distribution and are not measuring systematically e.g. along repeat transects at routine time and space scales. Thismeans the data cannot easily be used for closed box heat budgets, their Lagrangian paths could lead to sample aliasing, nor can we systematically assess observation impact. For these reasons we do not consider Argo data in this analysis.

It has been shown with data assimilation experiments that weekly subsurface temperature (XBT-like) observations have a significant impact on representation of the EAC: improving mean surface and subsurface circulation patterns, upper ocean heat content estimates (Gwyther et al., 2022), as well as baroclinic mode structure and eddy representation (Gwyther et al., 2023). However, the actual EAC XBT observing system only employs an approximately quarterly transect repeat time (or less). Hence, there is strong motivation to assess how impactful the existing XBT system is on representation of the EAC and its eddy field. Further it is useful to explore how representation of the EAC System is improved by increasing the observation sampling frequency.

This assessment is conducted with Observing System Simulation Experiments (OSSEs), which are a method of assessing observation impact, where a free-running simulation is taken as the true ocean or reference state, from which synthetic observations are extracted (e.g. Halliwell et al., 2014; Halliwell et al., 2015; Gasparin et al., 2019; Moore et al., 2020). These observations are then assimilated into a simulation that had an initial perturbation applied, and the resulting estimate can be compared to the ref state (see Figure 2 of Gwyther et al., 2022 for a schematic of the full OSSE procedure). We are thus able to assess the impact of the observing strategy on the representation of key ocean properties. In this study, we compare how several temporal sampling frequencies impact representation of subsurface temperature and eddy kinetic energy. This approach has the advantage of being able to assess the impact of a range of sampling configurations on ocean state estimates, without the time or cost of obtaining the ocean observations.

We present a series of model experiments that assimilate synthetic XBT observations with increased temporal resolution approximately matching the existing XBT transect network in the EAC System. In particular, we focus on assessing the impact of different sampling frequency on subsurface temperature fields and eddy kinetic energy, which have previously been shown as challenging to represent accurately (Gwyther et al., 2022; Gwyther et al., 2023). We characterise the ocean variability spectrum in order to demonstrate the potential for undersampling and aliasing by low frequency observations. Lastly, we explicitly separate the systemic error that is introduced by the data assimilation system from the endemic error resulting from inadequate representation ofthe mesoscale dynamics.

Numerical simulations are conducted with the Regional Ocean Modeling System (ROMS; Shchepetkin and McWilliams, 2005) which is a finite-difference primitive equation model with a terrain-following s-coordinate. The model domain extends from 27°S to 38°S and over 700 km offshore (with meridional grid resolution of 2.5 km linearly increasing to 6 km towards the east), and with constant meridional resolution of 5 km). There are 30 model layers in the vertical, with the model s-coordinate configured for more resolution in the surface boundary layer. This discretisation leads to cell thicknesses in the EAC of 1–3 m immediately below the surface, ~50–100 m thick cells in several hundred metres of water and ~300 m thick cells in the deep ocean below 3000 m. The model grid is rotated by 20° clockwise so as to approximately align the model coordinates with the along-shore and across-shelf directions (Model grid shown in Figure 1A . The bathymetry is sourced from the Geoscience Australia 50 m multibeam survey (Whiteway, 2009). The model domain is identical to that used in several other studies that explore the EAC velocity variability (e.g. Kerry and Roughan, 2020), biogeochemistry (e.g. Rocha et al., 2019), eddy dynamics (e.g. Li et al., 2021; Li et al., 2022b) and observation impact (e.g. Kerry et al., 2018; Gwyther et al., 2022; Gwyther et al., 2023).

We use two model configurations: a free-running and a data assimilating configuration. The free-running model uses lateral boundary forcing (currents, temperature and salinity conditions) from BRAN2020 (Chamberlain et al., 2021) and surface forcing conditions from the Bureau of Meteorology Atmospheric high-resolution Regional Reanalysis for Australia (BARRA-R; Su et al., 2019). More details are given in Gwyther et al. (2022) and Gwyther et al. (2023).

The data assimilating configuration used for the OSSEs is based on the model setup developed by Kerry et al. (2016), and uses an Incremental Strong Constraint 4-Dimensional Variational scheme (IS4D-VAR; e.g. Moore et al., 2011). This scheme calculates the differences between a free-running ‘forecast’ and observations over a chosen assimilation window (5 days in our case), where model and observations have associated error fields. The data assimilation scheme then generates new initial and boundary conditions such that a new ‘analysis’ simulation running with these adjusted initial, boundary and surface forcing conditions has minimised differences (in a least-squares sense) from the observations. The cycle then increments forward, using the previous analysis as the initial conditions for the new forecast. This data assimilating configuration has also been used in previous studies (e.g. Kerry et al., 2016; Kerry et al., 2018; Kerry et al., 2020; Siripatana et al., 2020) including for OSSEs (Gwyther et al., 2022; Gwyther et al., 2023). The lateral forcing conditions are from BRAN2020, while atmospheric conditions are sourced from the Australian Bureau of Meteorology’s ACCESS reanalysis (Puri et al., 2013). Both free-running and data assimilating configurations use a bulk flux parameterisation (Fairall et al., 1996) for calculating surface fluxes. This difference in surface forcing conditions between configurations is a necessary requirement for ‘fraternal twin’-type OSSEs. As summarised by Halliwell et al. (2014), a balance must be sought between slightly different configurations (forcing conditions, mixing parameters or parameterisations) that introduce error and realistically test the data assimilation system, and a long-term bias that the assimilation cannot correct for. An analysis of long-term mean bias resulting from the different forcing conditions showed that bias is minimal and constrained to the surface ocean, where it can be readily corrected by assimilating SST (see Gwyther et al., 2023).

We use the free-running configuration as the reference state (referred to as the ‘ref state’), to which a series of data-assimilating simulations (the OSSEs) are compared. Values are extracted from the ref state and a representative level of error is added as a normally distributed perturbation with standard deviation equal to the observation error. The observation error is set at 0.04 m for SSH, 0.5°C for SST, and with a depth-dependent profile for XBT observations, ranging from 0.6°C to 0.12°C (more information is given in Gwyther et al., 2022). These values are then taken as the synthetic observations, and are assimilated into a perturbed data-assimilating simulation (the OSSE). Following the procedure of (Gwyther et al., 2022), we generate a perturbed simulation by initialising the free-running simulation with an 8-day offset. This perturbation is enough to cause a free-running simulation to diverge significantly from the ref state, and is thusan effective test of the performance of the data assimilation system in assimilating observations into a realistically-diverged background state. Several different perturbations were tested, including a 1-day offset, 1-month, 1-year and a climatological state, but all were found to eventually generate similarly high levels of error (not shown). We thus chose a 8-day offset as it relatively quickly diverged, but it still had mesoscale features in the approximately correct locations, which is analogous to initialisation error in a true data assimilation system. For a more detailed description of the OSSE method, Gwyther et al. (2022) gives further information and includes a schematic of the procedure (their Figure 2).

Previous work has showed that the free-running model produces accurate seasonal and interannual representations of the EAC, including the eddy field, the separation latitude, eddy kinetic energy and volume transport (e.g. Kerry et al., 2016; Kerry and Roughan, 2020; Li et al., 2021). This demonstrates that the free-running ref state is representative of the true ocean. Hence, by analogy, the impact of assimilating the synthetic observations should translate to the true ocean.

This study assesses the performance of four different hypothetical XBT observing strategies through OSSEs, by comparing these against the ref state. All OSSEs assimilate the same synthetic surface observations representing SST and SSH by extracting observations from the ref state simulation with the appropriate location and timing, as per Gwyther et al. (2022; 2023). Synthetic SSH observations are based on the spatial and temporal coverage of the global ocean along-track multi-mission sea level altimeter data. Synthetic SST observations are based on the spatial and temporal coverage of the near-real time Himawari-8 satellite product. Each OSSE also assimilates subsurface XBT-like temperature observations, also extracted from the ref state,along two transects: at ~28.5°S and ~34°S, with a horizontal observation spacing of ~12.5 km at the continental shelf break. However for each experiment the temporal repeat time of the XBT transects is different, ranging from weekly (the ‘W-12.5’ OSSE), fortnightly (the ‘2W-12.5’ OSSE), monthly (the ‘M-12.5’ OSSE) to quarterly (the ‘Q-12.5’ OSSE). Experiment names, and temporal and spatial resolutions are shown in Table 1 .

Subsurface observations consist of vertical temperature profiles down to 900 m, with 10 m vertical resolution (or one observation per model layer for model vertical layer spacing greater than 10 m). These observations are taken along two transects, one in the north (referred to as XBT-N; ~28.5°S) and one in the south (referred to as XBT-S; ~34°S) of the domain. The two XBT lines are chosen to represent the approximate location and vertical sampling rates of the PX30 and PX34 lines ( Figure 1A shows transect locations and Figure 1B shows vertical distribution of observations along XBT-S). For ease of implementation, in our experiments, the XBT lines are oriented along grid coordinates (i.e. normal to the shore), whereas in reality the PX30 and PX34 XBT lines are along shipping routes from Brisbane to Nouméa, New Caledonia and Sydney to Wellington, New Zealand and thus at an angle to our model grid. The unique niche of the XBT dataset is the relatively fine spacing (10–100 km) between vertical profiles along defined repeat transects, allowing the resolution of a broad spectrum of processes, from eddies to basin-scale circulation (Smith et al., 1999).

All OSSEs assimilate surface observations (SSH and SST; see above) and subsurface temperature observations along both the northern transect and the southern transects, with horizontal XBT spacing of 12.5 km to 30 km (every 5 model cells) and 5-days to sample the transects. Each OSSE has different XBT transect repeat times: The W-12.5 OSSE has a transect repeat every week, the 2W-12.5 OSSE has transect repeats every two weeks, the M-12.5 has transect repeats every month (30 days) and the Q-12.5 has transect repeats every quarter year (90 days). The Q-12.5 OSSE represents the true Scripps XBT lines most closely in spatial resolution and temporal repeat time. The other OSSEs (W-12.5, 2W-12.5 and M-12.5) represent how a higher-frequency sampling scheme will impact simulated representation of the ocean.

There are different modes of variability present in the ocean, from fast (e.g. tidal) to slow (e.g. climate modes) timescales and from small (e.g. eddy cascade) to large (e.g. gyre circulation) spatial scales. When designing an observation strategy, a choice must be made as to what portion of this period-wavelength phase spectrum should be sampled. Processes that are outside of the sampled portion of the spectrum are not measured enough (in time or space) to discern the process and may be aliased into thesampled portion of the spectrum. These processes either require more rapid or finer spaced sampling (i.e. for discerning small or fast scale processes), or, longer or broader sampling (i.e. for capturing large scale or slow processes). Likewise, model resolution must be fine enough that any small scale processes that are observed can actually be resolved by the model simulation – though there could be inherent limits to predictability at very fine, submesoscale resolutions (Sandery and Sakov, 2017).

We have shown that in the EAC (using our ref state simulation), ocean surface processes display a wide variety of scales, from weekly changes to seasonal variability in SST ( Figure 2A ). In the subsurface ocean at 500 m, the most obvious processes are eddy dynamics, which have small scale and fast changes in temperature and EKE ( Figures 2C, D ). Indeed, inspection of the variability period-wavelength spectrum of EKE at 500 m shows two key features: higher power at the large scales (>100 km and 30–180 days), and a band of increased power stretching from ~350 km/monthly to 30-40 km/daily wavelengths and periods, which likely represents the increased variability associated with the Rossby-mesoscale-submesoscale energy pathway ( Figure 3A ). Sampling any less frequently than monthly will truncate a large portion of this mesoscale energy pathway. Given that we can quantify the amount of variability with different spatial scales, we can directly describe the portion of the total variability that is captured by repeat sampling at chosen frequencies. In particular, for EKE at 500 m, quarterly sampling captures ~20% of total variability, while increasing to monthly sampling captures ~50%. More rapid sampling rates at fortnightly and weekly frequencies captures ~70% and ~90% of the total EKE variability ( Figure 4 ). Note that since data assimilation systems will not exactly match observations (for e.g. due to observation error), we would only expect these patterns to be approached in the long-term.

While there is no obvious reduction in transect-mean error when sampling is increased from quarterly to fortnightly, reduced RMS is noticeable with weekly sampling ( Figures 5A–C ). This shows that infrequent observations (i.e. quarterly through to fortnightly subsurface observation repeats) have limited impact on constraining mean subsurface conditions, and highlights the importance of regularly assimilating subsurface observations to maintain an accurate estimate.

This is further confirmed in the vertical transects of temperature RMS, where the highest error region (~250 m) is represented with similar error in quarterly, monthly and fortnightly sampling, but is improved with weekly sampling (cf. Figures 6A–D ). Horizontal fields of temperature RMS show that the greatest improvement in representation is indeed in the 250 m depth range, and particularly in the high eddy energy region. This confirms that repeat sampling of at least weekly frequency is required to improve the mean representation of heat in the critical 250 m depth range in the northern upstream region ( Figure 7A ), separation region ( Figure 7E ) and southern high eddy region ( Figure 7I ). Indeed, this confirms the importance of higher repeat frequency observations as suggested by the variability analysis, where sampling at fortnightly–weekly frequencies is required to capture 60–80% of variability in the higher EKE region ( Figures 4A, B ).

Systemic error is introduced by the data assimilation system and can only be mitigated through improvements to the data assimilation system. This error could present as overly deep eddies (Siripatana et al., 2020), incorrect vertical distribution of temperature and heat content (Gwyther et al., 2022), or baroclinic mode structure (Gwyther et al., 2023). Endemic error results from the resolution of mesoscale processes and is improved by faster or higher resolution sampling. Maps of bias ( Figure 8 ) suggest that the data assimilation process is producing a temperature structure that is too cold between 250–500 m and too warm at 1000 m. The location of the strongest bias suggests that eddies suffer disproportionately from this systemic error. The slight decrease in bias that is observed in the W-12.5 OSSE, particularly for the eddy region, indicates that weekly subsurface temperature observations are playing an important role in repeatedly correcting the vertical structure. The bias-corrected temperature RMS shows that the largest improvement in representation is in the 250 m depth range, which is achieved by taking XBT measurements at faster than monthly frequency. Improvement at 500 m is smaller, though, as most error is concentrated in the 250 m range, this may be acceptable. The improvement in increasing XBT sampling time from quarterly to monthly (see reduction in RMS in the separation zone; Figure 9D ) may represent improvements in capturing seasonal-scale processes, such as mean separation latitude (Ypma et al., 2016; Oke et al., 2019), jet core velocity and EKE (Archer et al., 2017). The comparatively larger reduction in RMS observed when moving to fortnightly or weekly sampling likely represents better representation of mesoscale dynamics. This is supported by the frequency analysis (Section 3.1) and by the largest RMS reduction being in the eddy region ( Figures 9A, B ). In Figure 7 , the greatest improvement in RMS error, particularly in the eddy region, is achieved by weekly sampling. In contrast, Figure 9 shows that fortnightly sampling achieves sufficient reduction of bias-corrected RMS error in the eddy region. This indicates that the gain achieved from moving from fortnightly to weekly sampling is through reduction of the bias – likely due to the presence of subsurface observations in each assimilation window.

The large temperature bias ( Figure 8 ) likely results from the data assimilation process itself. This could result from the specification of background error covariance matrix, which controls the influence of observations in the horizontal and vertical (Bannister, 2008) and the weighting of the model forecast (Lee and Huang, 2020). Correctly specifying the background error covariance matrix is a major challenge to accurate data assimilation simulations (see review in Moore et al., 2019). Indeed, several studies have shown that 3-dimensional structure suffers in data assimilation, even with the inclusion of the limited subsurface observations (see for e.g. Zavala-Garay et al., 2012; Pilo et al., 2018; Siripatana et al., 2020; de Paula et al., 2021; Gwyther et al., 2022; Gwyther et al., 2023).

In any case, the result is that subsurface representation needs to be constrained frequently (i.e. weekly), otherwise subsurface structure drifts too far from the truth and any improvement from observations is minimal, though this result may change for different assimilation techniques. This suggests that improvements to assimilation schemes could improve representation of the subsurface structure even in the absence of observations. Furthermore, existing observations (i.e. quarterly subsurface measurements) could have more impact than they currently do, and potentially, sampling on the mesoscale timescale would be sufficient.

Operational and research-focussed data assimilation systems benefit greatly from the high spatial resolution, coverage and relatively rapid temporal repeat times of satellite-inferred sea-surface measurements. As a result, these systems can provide accurate and representative estimates of sea surface temperature and height (and hence geostrophic currents). However, subsurface conditions such as temperature and velocities (Gwyther et al., 2022), isothermal slopes and baroclinic mode structure (Gwyther et al., 2023) are not represented with the same accuracy. These simulation inaccuracies are further compounded when trying to simulate dynamically complex 3-dimensional structures such as eddies (Gwyther et al., 2023). As a result, the assimilation of subsurface observations is critical to improving representation. While some subsurface observations are routinely assimilated into operational ocean models, many of these are sampling sparsely (e.g. Argo floats), or were designed for climatological studies(e.g. the Scripps high resolution XBT program). So, there is some justification for the re-purposing of existing subsurface sampling frameworks, at potentially faster repeat times, so as to be of more use for operational models that estimate and forecastocean conditions.

Our results have showed that frequencies of variability of interest need to be considered when assimilating subsurface observations. We have explored the benefit of assimilating data from an XBT observing network designed to observe climatological-scale processes i.e. heat fluxes through ocean basins over interannual periods. Our results show that assimilating XBT data at fortnightly repeat times would give a better representation of higher frequency processes such as mesoscale eddies. However, we also show that systemic errors in the data assimilation process itself limit the ability to represent accurate vertical structure. Improvements to the assimilation scheme to reduce systemic error and biases would therefore increase the positive benefit obtained from current and future subsurface observation systems.

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

DG, MR, CK, SK conceived the study and designed and experiments. DG and CK configured the data assimilating model. DG performed the analyses. DG and MR wrote the first draft of the manuscript. All authors contributed to the manuscript editing and revisions. All authors contributed to the article and approved the submitted version.