Recent Advances in Japanese Fisheries Science in the Kuroshio-Oyashio Region through Development of the FRA-ROMS Ocean Forecast System : Overview of the Reproducibility of Reanalysis Products

To address various fisheries science problems around Japan, the Japan Fisheries Research and Education Agency (FRA) has developed an ocean forecast system by combining an ocean circulation model based on the Regional Ocean Modeling System (ROMS) with three-dimensional variational analysis schemes. This system, which is called FRA-ROMS, is a basic and essential tool for the systematic conduct of fisheries science. The main aim of FRA-ROMS is to realistically simulate mesoscale variations over the Kuroshio-Oyashio region. Here, in situ oceanographic and satellite data were assimilated into FRA-ROMS using a weekly time window. We first examined the reproducibility through comparison with several oceanographic datasets with an Eulerian reference frame. FRA-ROMS was able to reproduce representative features of mesoscale variations such as the position of the Kuroshio path, variability of the Kuroshio Extension, and southward intrusions of the Oyashio. Second, using a Lagrangian reference frame, we estimated position errors between ocean drifters and How to cite this paper: Kuroda, H., Setou, T., Kakehi, S., Ito, S., Taneda, T., Azumaya, T., Inagake, D., Hiroe, Y., Morinaga, K., Okazaki, M., Yokota, T., Okunishi, T., Aoki, K., Shimizu, Y., Hasegawa, D. and Watanabe, T. (2017) Recent Advances in Japanese Fisheries Science in the Kuroshio-Oyashio Region through Development of the FRA-ROMS Ocean Forecast System: Overview of the Reproducibility of Reanalysis Products. Open Journal of Marine Science, 7, 62-90. http://dx.doi.org/10.4236/ojms.2017.71006 Received: October 2, 2016 Accepted: December 10, 2016 Published: December 13, 2016


Introduction
In the last decade, many ocean forecast systems that include some sort of data assimilation method have been developed to nowcast and forecast ocean states with a particular focus on mesoscale variability; these include HYCOM/NCODA [1], Mercator [2], FOAM [3], MOVE/MRI.COM [4], JCOPE2 [5], and so on.Progress in this area is attributable primarily to the Global Ocean Data Assimilation Experiment project, which, crucially, supports near-real-time collection, processing, and management of oceanographic data [6] [7].In addition, several state-of-the-art community ocean models, which can be coupled with data assimilation schemes, have been developed as open-source code (e.g., ROMS [8] and ECCO [9]) and have contributed to making ocean nowcasting and forecasting relatively easy.Consequently, physical oceanographic studies have made rapid progress in the areas of ocean dynamics and predictability.
Utilization of outputs from ocean forecast systems is not restricted to physical oceanographic studies.Data assimilation models initialize ocean state variables for forecasts by the synthesis of spatio-temporally limited oceanographic data and ocean circulation model.These initialized state variables, which can be considered to have been optimally estimated by combining observation and simulation, are generally referred to as "analysis" or "reanalysis" data.Such assimilated products are needed for many scientific and practical applications where realistic oceanographic conditions are required: in marine navigation systems used for transportation safety (e.g., [10]); for acoustic propagation forecasting (e.g., [11]); and in investigations of the transport of marine debris (e.g., [12] [13]), oil spills (e.g., [14] [15]), and radioactive materials (e.g., [16] [17]).
In recent years, there has also been an intensive demand for assimilated model products in fisheries science.These products have been applied from both Eulerian and Lagrangian viewpoints to fisheries science with various aims.Assimilated data have been analyzed instead of observed data to clarify the linkage between fish abundance and oceanographic conditions [18], and they have been used for stock management applications [19].Survival, growth, migration, and recruitment of important commercial fish have been examined by using particle-tracking models and individual-based models (IBMs) with assimilated oceanographic conditions (e.g., [20] [21] [22] [23]).Several studies have developed high-resolution nested ocean models with data assimilation products (e.g., [24]) in an effort to understand and predict the sudden occurrence of harmful marine organisms such as jellyfish, toxic algae, and viruses.Potential fishing areas have also been identified by examining statistical relationships between fish habitats and simulated oceanographic conditions (e.g., Habitat Suitability Index (HSI) [25] [26]).Additionally, the spatio-temporal variability of nutrients and plankton, which is associated with food availability for fish, has been simulated by an ocean forecast system coupled with a lower-trophic-level ecosystem model (e.g., [27] [28]).
The Japanese Islands are located in the Kuroshio-Oyashio western boundary current region, where mesoscale variations are prominent (Figure 1(a)).In addition, the islands are characterized by a narrow continental shelf (typically less than 50 km wide), with the result that mesoscale variations intermittently affect coastal waters through interactions with submesoscale variations [29] [30].Moreover, mesoscale variability in the Kuroshio-Oyashio region is a key issue for Japanese fisheries and fisheries science, because many pelagic fish, which are important targets of Japanese fisheries, extensively use the Kuroshio-Oyashio region as spawning and nursery grounds.Among these are the Japanese sardine (Sardinops melanostictus) [31], Japanese anchovy (Engraulis japonicus) [32], Pacific saury (Cololabis sairai) [33], Japanese chub mackerel (Scomber japonicus) [34], Japanese common squid (Todarodes pacificus) [35], and various other species [36] [37].For example, Japanese sardine adults migrate southward from the Oyashio to the Kuroshio region for spawning.During winter and spring, eggs, larvae, and juveniles are transported eastward by the Kuroshio, the Kuroshio Extension, or related mesoscale flows, and then juveniles migrate northward to the Oyashio nursery region via the Kuroshio-Oyashio transition zone [38] [39].Japanese fisheries capture sardine larvae and juveniles in the Kuroshio region during the passive transport period, and they capture adults mainly in the Oyashio region during the nursery and migration period.Hence, an exhaustive dataset across the Kuroshio-Oyashio region is required to study relationships between oceanographic conditions and individual life stages of even this single pelagic fish.
The mission of our organization, the Japan Fisheries Research and Education Agency (FRA), which is one of the National Research and Development Agencies of Japan, is to conduct a wide range of research and development activities related to fisheries and to address many fisheries problems through basic research and practical applications.To respond to the demands of multiple research projects on fisheries oceanography, management, and economics, FRA oceanographers have had to prepare an operational ocean forecast system to create a basic and versatile oceanographic dataset.Such an ocean forecast system had been also expected to extensively enable Japanese prefectural fisheries institutes as well as FRA to systematically conduct many kinds of fishery researches and establish many practical applications, which contribute significantly to sustainable fishery and safety of fishery products around Japan.Although the Japanese Oceanographic Society has several operational ocean forecast systems, including MOVE (Multivariate Ocean Variational Estimation) [4], JCOPE2 [5], and DREAMS [40], none of these forecast systems is specialized to fisheries science.Therefore, the FRA has developed an ocean forecast system based on the Regional Ocean Modeling System (ROMS) [41] and specialized to fisheries science called FRA-ROMS.Here, some readers may question why and how FRA-ROMS is specialized for fisheries science.These questions are addressed in Section 4, where ongoing community efforts of the FRA are introduced.
The configurations of FRA-ROMS are described in Section 2, and the reproducibility of its reanalysis products is examined in terms of both Eulerian and Lagrangian reference frames in Section 3. In Section 4 we summarize recent and ongoing fisheries studies that have used FRA-ROMS, and discuss prospects for work in the near future.

Configurations of FRA-ROMS
FRA-ROMS consists of an ocean circulation model based on ROMS and three-dimen-sional variational (3D-Var) objective analysis schemes, which are not included in the original ROMS code.The configurations of FRA-ROMS are briefly explained below.

Ocean Circulation Model
FRA-ROMS consists of 1/2-and 1/10-degree ocean models connected by one-way nesting, in which volume transport across the lateral boundaries is adjusted according to the method proposed by [42].The 1/2-degree parent model covers almost the entire North Pacific region (Figure 1(b)).This model was designed to simulate basin-scale variations associated with El Niño-Southern Oscillation and the Pacific Decadal Oscillation and impose them at the lateral boundaries of the 1/10-degree child model as external forcings.Thus, daily mean outputs from the 1/2-degree model are imposed at the lateral boundaries of the 1/10-degree model.The domain of the 1/10-degree model is limited to the Northwestern Pacific (Figure 1(b)) to simulate mesoscale variations over the Kuroshio-Oyashio region.Vertical positions in both models are defined by the specific S-coordinate function proposed by [43], and the vertical resolution of both models is 48 levels.Since a sea-ice model is not coupled with the current version of FRA-ROMS, simulated temperature and salinity (TS) values over the Okhotsk Sea and the Bering Sea in the two models are restored to the monthly mean climatology obtained from World Ocean Atlas 2001 (WOA2001) [44] [45] on a timescale of 30 days.
Numerical schemes and parameterizations that have been adapted for the two models are summarized by [29].
The two models were integrated for 20 years from the WOA2001-based initial conditions by using climatological monthly mean heat and momentum fluxes at the sea surface [29].Then, both models were driven by 6-hourly mean atmospheric and radiation elements (i.e., wind speed, air temperature, sea level pressure, specific humidity, precipitation, and net shortwave and downward long wave radiations) from the Japanese 25year Reanalysis (JRA-25) dataset [46], from which heat and momentum fluxes were sequentially estimated by the COARE (Coupled-Ocean Atmosphere Response Experiment) bulk formula [47] during each model run.

Data Assimilation Based on 3D-Var Analysis
To perform reanalysis experiments for the period from 1993 to 2012, 3D-Var analysis schemes were constructed by a technique similar to that used by MOVE, an ocean forecast system developed by the Meteorological Research Institute [4].3D-Var analysis schemes were originally developed for a z-coordinate ocean model, but they can be stably applied to a sigma-coordinate model such as the Princeton Ocean Model [5].
Therefore, this scheme can be expected to work robustly when coupled with S-coordinate ocean models such as ROMS.The 3D-Var analysis schemes applied here should also create reanalysis data more precise than two-dimensional schemes (e.g., Optimal Interpolation).Meanwhile, the 3D-Var schemes may be less perfect than four-dimensional (4D-Var) schemes and ensemble Kalman Filter, because background errors in the 3D-Var schemes are implicitly assumed to be stationary.In addition, heat and salt budget in reanalysis data based on the 3D-Var schemes are not completely conserved, as explained below.Nevertheless, we selected the 3D-Var analysis schemes, which had an advantage that computational cost was much smaller and were more speedy/practical for fisheries science, where additional heavy computations such as PTMs are needed after data assimilation.
where x , f x , and o x are column vectors whose elements are the state variables, the first guess, and the observed values, respectively.B is the horizontal background correlation matrix, R is the observed error covariance matrix, and H is an operator that picks up values of temperature and/or salinity at the observation position.h denotes the dynamic topography estimated from x , o h is the altimetry-derived sea- level anomaly (SLA), and 2 h r is the standard SLA observation error.A vertically coupled TS empirical orthogonal function (EOF) modal decomposition was adapted for B .
The 10 dominant modes of the coupled TS EOF were estimated prior to the objective analysis from six sub-regions of the 1/2-degree model and nine sub-regions of the 1/10-degree model [4].
The fourth term on the right-hand side of Equation ( 1) is a constraint that was proposed by [51] to avoid extremely low temperature estimates.We confirmed that without this term extremely low temperatures appeared frequently in the reanalysis data, particularly in the Kuroshio-Oyashio transition zone (not shown).The low temperatures appear because without this term the cost function basically assumes a Gaussian distribution of TS errors, even when the temperature is around 0˚C and close to the freezing point of seawater (i.e., about −1.7˚C).
Optimal TS values estimated with Equation ( 1) are referred to as analysis data in this paper.The analysis data were assimilated into the ocean model by incremental analysis updating [52] on a 7-day cycle as follows.The ocean model was freely integrated from an initial value for 3.5 days, and the simulated TS values at the end of the 3.5 days were used as first guess values in Equation ( 1).After optimal TS values (i.e., analysis data) were obtained, the difference in TS between the first guess and the analysis data was estimated.This difference was divided by the number of internal-mode time steps, each corresponding to 7 days, and the result was added to the TS governing equations as an increment.The ocean model was again integrated from the same initial value for 7 days with the increment as forcing (hereafter, forced run).These procedures were repeated sequentially on a 7-day cycle.Assimilated data were sea surface temperature (SST), altimetry-derived SLA, and in situ TS profiles (Table 1).In this study, reanalysis experi-

Reproducibility of the Reanalysis Product
We examined the reproducibility of the daily reanalysis data for the 20 years from 1993 to 2012 by comparing the simulated results with previous oceanographic knowledge and observed data over the Kuroshio-Oyashio region in both Eulerian and Lagrangian reference frames.

Analyses from an Eulerian Viewpoint
We examined the main flows in the Kuroshio-Oyashio region by mapping stream func-  Comparison of the standard deviations between the simulated current velocity and altimetry-derived geostrophic velocity at the sea surface (Figure 4) showed that the spatial patterns of current variability were very consistent.The simulated velocity used insufficiently removed from the AVISO altimetry data [62].In addition, in the subarctic region north of 40˚N, the current variability was greater in the reanalysis data.This discrepancy suggests that variations with a timescale of less than a week (e.g., in windinduced currents within the Ekman layer [63], which are associated with ageostrophic components), contributed more significantly to the current velocity at the sea surface in the subarctic region.
Daily SST based on reanalysis data was compared with a daily gridded SST dataset with a horizontal resolution of 0.25˚ [64], which was created by NOAA by using optimal interpolation to merge satellite, ship, and buoy observations (product name: AVHRR_OI).The comparison was performed for the area from 20˚N to 50˚N and 120˚E to 160˚E, but excluding the Okhotsk Sea.Root-mean-square differences (RMSDs), mean differences, and correlation coefficients between the two SST datasets were estimated at monthly intervals (Figure 5(a)).Correlation coefficients (not shown) were greater than .These large RMSDs were partly because in the reanalysis data the SST fronts tended to be sharper than they were in the AVHRR_OI data, which were spatially smoothed by the optimal interpolation.In contrast, positive biases were apparent in summer (Figure 5(a)), reflecting the fact that summertime SSTs in the reanalysis data tended to be higher than those in the AVHRR_OI data (Figure 5(c)).Relatively large RMSDs, greater than 2˚C, were related to this positive bias and showed a scattered distribution along the coasts of the Yellow Sea and along the southern coasts of the Kuril and Four Islands northeast of Hokkaido (Figure 5(e)).These regions are characterized by strong tidal currents and mixing, which were not parameterized or simulated by FRA-ROMS.The positions of the axes of the Kuroshio and Kuroshio Extension were specified from the reanalysis data by combining the methods of [65] and [66].Used alone, the former method could not successfully specify the axis position of the Kuroshio Extension, where mesoscale eddies interact intermittently with the mainstream axis.Therefore, the Kuroshio axis position within the 123˚E -140˚E interval was first determined by the former method.Then, following the method of Qiu and Chen [66], sea surface heights within the 136˚E -140˚E interval were averaged along the pre-determined axis position.The position of the sea surface height isoline with the mean along-stream value was extracted and used as the axis position from 136˚E to 160˚E.Finally, the Kuroshio axis positions in the 136˚E -140˚E interval estimated by the two methods were averaged.
The 20-year mean axis position of the Kuroshio and its standard deviation were compared between the reanalysis data and observation-based QBOC (Quick Bulletin of Ocean Conditions) data digitized by the Marine Information Research Center (MIRC) (Figure 6(a)).Although some regional differences could be detected, the mean and standard deviation of the reanalysis data seemed to be very consistent with those of the observation-based data.The axis-position difference between the two datasets in the interval from 126.5˚E to 140.5˚E was defined following [67] [68] as the ratio of Q to L, where Q denotes the area between the two axes and L is the axis length in the MIRC observation data.The difference fluctuated temporally but was relatively stable overall throughout the period analyzed (Figure 6(b)).The mean difference was 20.2 km, which is almost the same as the width of two grid cells in the 1/10-degree model, the minimum was 7.5, the maximum was 50.9 km, and the standard deviation was 5.5 km.
Year-to-year variations in the simulated axis position of the Kuroshio Extension from 1993 to 2012 are illustrated in Figure 7.In both the simulation and observation [66], the Kuroshio Extension exhibited stable and unstable modes.We compared the path length and the meridional position of the axes in the simulation with the observation-based results of [69] (Figure 8).The simulated and observed values were similar, suggesting that the bimodal variation of the Kuroshio Extension was successfully simulated in the reanalysis data.The path of the Kuroshio Extension seemed to show a higher frequency of variation in the reanalysis data compared with the observation; this difference can be attributed to the difference in time resolution between the weekly data utilized by [69] and the daily reanalysis data.
We also examined the reproducibility in the reanalysis data of the water-mass distribution in the Kuroshio-Oyashio transition zone, where the Kuroshio and Oyashio waters are distributed and intermixed in a complicated way [70].The Kuroshio-derived waters are mainly supplied to the transition zone by the Kuroshio Extension and the Tsugaru Warm Current (Figure 1(a)).Oceanographic conditions in the transition zone are conventionally described by using spatial patterns of isotherms of an index temperature at depths of 100 or 200 m [70].We used this method to examine the reproducibility of the southernmost latitudes of the first and second Oyashio intrusions and of the easternmost longitude of the Tsugaru Warm Current (Figure 1(a)).The southernmost latitudes of the Oyashio intrusions can be used to represent the Oyashio water distribution (Figure 1(a)) and are closely related to the formation of fishing grounds for pelagic fishes (e.g., [71] [72] [73]).The easternmost longitude of the Tsugaru Warm Current is useful not only for representing the water-mass distribution but also for assessing the capability of a mesoscale model, because from summer to autumn the Tsugaru Warm Current forms a clockwise gyre with a maximum diameter of 100 km [74], which is about the minimum magnitude of variation that mesoscale models with a grid size of about 10 km are generally able to reproduce.
The southernmost latitudes and easternmost longitude were determined by using isotherms of 5˚C and greater than 6˚C, respectively, at a depth of 100 m, as explained in the caption of Figure 9.The southernmost latitude of the first Oyashio intrusion was very consistent between the reanalysis and observed values (Figure 9(a)).Both seasonal and interannual variations seemed to be reproduced well in the reanalysis data.The southernmost latitude of the second intrusion, however, was not precisely reproduced in the reanalysis data (Figure 9(b)), mainly because the observation data, which were obtained by an in situ ship survey, were frequently limited to coastal waters near Honshu and Hokkaido Islands; thus, the estimation accuracy of the observed temperature map used to determine the southernmost latitude of the second Oyashio intrusion was probably temporally inhomogeneous and dependent on the observed temperature distribution.
The easternmost longitude of the Tsugaru Warm Current showed good agreement between the reanalysis and observation data (Figure 9(c)), though there were some slight discrepancies.In the reanalysis data, the easternmost longitude exhibited an

Analyses from a Lagrangian Viewpoint
Particle-tracking methods, including use of IBMs and super-IBMs, are among the most powerful techniques used in applications of reanalysis data to fisheries science (e.g., [75]- [81]).However, particle-tracking errors in a Lagrangian reference frame have hardly been evaluated, mainly because no well-established error evaluation method has been shared among individual modelers developing ocean forecast systems.Hence, here we propose a method, which we call the best performance particle (BPP) method, in which particle-tracking errors of mesoscale models are estimated by comparing the positions of numerically tracked particles with those of actual surface drifters.For validation data, we used Surface Velocity Program (SVP) drifter trajectory data that had passed a sophisticated quality check performed by NOAA'S Atlantic Oceanographic and Meteorological Laboratory.Only data for drifters with a drogue were utilized.SVP drifters are designed to measure mixed layer currents in the upper ocean; thus, the center position of the drogue is set at 15 m below the sea surface, and its vertical length is about 6.44 m.Daily mean positions of the drifters were estimated from original 6-hourly data.Data from drifters with a tracking period (=N) longer than 60 days and without any missing data were selected.The data from each drifter were divided into subsets, each consisting of 60 days of daily observations.Thus, the number of the subsets for a single drifter was N-59.A total of 57,608 subsets were available for comparison with particle-tracking results in the Kuroshio-Oyashio region (Figure 10(a)).
A particle-tracking experiment was performed with FRA-ROMS for 60 days, and the results were compared with the daily drifter positions in each data subset.The 60-day period was chosen because for the first few months after spawning, the swimming capacity of eggs, larvae, and juveniles of pelagic fish is generally negligible (e.g., [21] [23] [78] [82]).A particle-tracking model (LTRANS) specialized for ROMS output [83] was modified and applied in this study.Particles were transported by the simulated current velocity at a fixed depth without any diffusion.A group of particles was initially centered around the actual drifter position on the beginning day of each data subset.The particles were evenly spaced at intervals of 0.0267˚ (~3 km) of latitude and longitude within ±0.4˚ (~40 km) of latitude and longitude of the drifter position and at depths of 5, 10, 15, and 20 m.A total of 3844 particles (961 at each depth were tracked) for each data subset.RMSEs were estimated from difference in horizontal position between the drifter and the individual particles at daily intervals for 60 days.Finally, the single particle with the minimum RMSE value, that is, the BPP, was selected and its statistical properties were analyzed.Multiple particles were tracked, and the BPP among them was selected both to reduce the effect of potential bias in the reanalysis data and because some conditions that would affect drifter motion were not appropriately considered in the particle-tracking experiment, as explained in greater detail in Appendix 1.
If, under idealized conditions, the simulated current velocity includes a constant unidirectional error of 2 cm/s (3 cm/s), the position error and RMSE after particle tracking for 60 days would be 104 and 60 km (156 and 90 km), respectively.Therefore, we considered RMSEs of less than 100 km to indicate relatively accurate reproducibility of the reanalysis data.In fact, as a criterion for evaluating tracking errors, the RMSE may be rather strict, because it was often relatively large when the trajectories of the particle and drifter after tracking for 60 days were similar but their daily positions were different (not shown).
The RMSE of the BPP was weakly correlated with the trajectory length of the drifter integrated over 60 days (r = 0.45; Figure 11), the mean RMSE and mean integrated trajectory length were 149 km and 1430 km, respectively, and their ratio was about 0.1.
We examined the spatial distribution of integrated trajectory lengths (Figure 10 regardless of the initial particle position, the ratios were mostly less than 0.1 across the Kuroshio-Oyashio region, although ratios exceeding 0.2 were scattered in the region (Figure 10(d)).This result suggests that the particle-tracking error of the BPP can be roughly estimated by multiplying the integrated trajectory length of the drifter at the particle's initial position by 0.1.

Ongoing Studies and Near-Future Plans
The enterprises are also being developed to clarify processes leading to year-to-year variations in recruitment and stocks, and to improve stock assessment accuracy, of important commercial fish such as Japanese sardine (Sardinops melanostictus) [84], Japanese chub mackerel (Scomber japonicus) [84], Japanese common squid (Todarodes pacificus) [85], Japanese walleye pollock (Theragra chalcogramma) [84], and chum salmon (Oncorhynchus keta) [86].Moreover, sometimes in summer, the giant jellyfish (Nemopilema nomurai) seriously affects the Japanese fishery, so the process and pathway by which this species is transported from the East China Sea to around Japan are also being investigated and predicted by particle-tracking experiments conducted with FRA-ROMS.
Lower-trophic-level ecosystem models lacking the ROMS default options have been tuned and coupled with FRA-ROMS.For example, a simple NPZD (nutrient, phytoplankton, zooplankton, detritus) ecosystem model [87] has been adapted for examining the nutrient supply process to the shelf-slope region facing the Kuroshio south of Japan [88].A more complicated ecosystem model, eNEMURO [89], has also been coupled with FRA-ROMS to project future food conditions for pelagic fish and to evaluate impacts of global warming on Japanese fisheries (Okunishi and Setou, personal comm.).
Many downscaling models that can be linked to the reanalysis and prediction data of To maintain and enhance the reproducibility of FRA-ROMS, an oceanographic data acquisition system has been developed by oceanographers at the FRA.TS profiles obtained by regular monitoring by Japanese prefectural fisheries institutes have been since 2007 distributed to the Global Temperature-Salinity Profile Program (GTSPP) [93].
Recently, new types of observation instruments, such as glider and underway CTDs, have been adapted to measure fishery and oceanographic conditions.These instruments have made it possible to obtain high-resolution measurements able to resolve mesoscale and submesoscale variability.The FRA has developed a system to transfer these high-resolution measurements to FRA-ROMS in the TESAC format [94].Such efforts cannot only improve reproducibility of the reanalysis data (i.e., initial state variables for forecasts) but also enhance predictability of the operational version of FRA-ROMS.
FRA-ROMS is one of the few operational ocean forecast systems in the world that is specialized and optimized for fisheries science.In particular, FRA-ROMS embeds a versatile subsystem with multiple functions so that model outputs can be efficiently applied to fisheries science (Okunishi, personal comm.).Using this subsystem, even a scientist without any special knowledge of physical oceanography and numerical techniques can easily extract, analyze, and illustrate necessary elements from the model outputs.Moreover, a function to conduct particle-tracking experiments is being developed.This subsystem is available to not only all FRA scientists but also registered scientists affiliated with research institutes of the Japanese prefectural fisheries.In the near future, several more useful functions will be developed and the user-friendliness of the subsystem will be enhanced.
Plans exist to upgrade FRA-ROMS.One modification will extend the model domains.The eastern and southern boundaries of the 1/10-degree model will be extended to the North American continent and to 5˚S -20˚S, respectively.The 1/2-degree model domain will also be optimally expanded.Further, the 1/2-and 1/10-degree models will be connected by applying a two-way nesting technique.In addition, several effects (e.g., tidal mixing, sea ice, and freshwater discharge from land) not considered in the current version, will be incorporated into the model.Moreover, a high-resolution (1/50 degree) coastal model is being developed.This model, which is an extended version of the model used by Kuroda et al. [29] [30] [76], is designed to cover the entire Japanese coastline, not only in the Kuroshio-Oyashio region but also along the Japan Sea.The use of the FRA-ROMS platform and its associated models, combined with a more userfriendly subsystem, is expected to lead to dramatic progress in Japanese fisheries science in the near future.

Figure 1 .
Figure 1.(a) Schematic view of the Kuroshio-Oyashio system.The main currents are denoted by arrows, and geographical names are enclosed in boxes.1st and 2nd indicate the first and second Oyashio intrusions, respectively.Contour lines with an interval of 10 cm indicate 20-year mean sea surface height estimated from the AVISO maps of absolute dynamic topography during 1993-2012 (i.e., MDT plus 20-year mean SLA).(b) Model domains and bathymetry (contour).The entire illustrated area is the 1/2-degree model domain, and the region within the dashed rectangle is the 1/10-degree model domain.
ments were conducted for 1993-2012, and daily mean outputs during the forced run (i.e., reanalysis data) were analyzed.It should be borne in mind that FRA-ROMS has two different versions, an operational version and a long-term reanalysis experimental version.The former has been used operationally since May 2012 for nowcasts and two-month forecasts at intervals of tions of the simulated 20-year mean volume transport from the sea surface to a depth of 1000 m or the sea bottom (Figure2).The transport through the Tokara Strait, which is a good indicator of the Kuroshio transport through the East China Sea, was estimated as 26 Sv (=10 6 m 3 •s −1 ), a value comparable to observed transports of 21.9 Sv[53],23.4Sv[54], 25 Sv[55], and 28 Sv[56].At the eastern mouth of the strait, this transport joins with the northward transport of 18 Sv associated with the Ryukyu Current.The

Figure 2 .Figure 3 .
Figure 2. Stream function of the simulated 20-year mean transport from the sea surface to a depth of 1000 m or the sea bottom.Major transports referred to in the main text are emphasized by colored arrows with numbers showing the associated volume transport (Sv = 10 6 m 3 •s −1 ): green arrow, through the Tokara Strait; small red arrow, Ryukyu Current; large red arrow, Kuroshio; orange arrow, Kuroshio recirculation; gray arrow, Oyashio.

Figure 4 .
Figure 4. Standard deviations of current velocity at the sea surface in (a) reanalysis data and (b) altimetry-derived geostrophic velocity data.The standard deviations, which were separately estimated from the east and north velocity components and then superposed, are depicted by the color scale.The vectors show the 20-year mean velocity at the sea surface (as in Figure3).

Figure 5 .
Figure 5. (a) Monthly time series of the RMSD (thick line) and mean difference of SST (thin line) between the reanalysis data and AVHRR_OI.The difference was estimated by subtracting AVHRR_OI values from the reanalysis data.The target area was limited to the area from 20˚N to 50˚N and 120˚E to 160˚E, excluding the Okhotsk Sea.Scatter plots of the two SST datasets are depicted as a frequency distribution (color scale (unit: number of data points for each class)) for (b) January and (c) July.Spatial maps of the RMSD of the SSTs during 20 years (contours; contour interval, 1˚C) are shown for (d) January and (e) July; RMSDs exceeding 0.5˚C are emphasized by color.

Figure 6 .
Figure 6.(a) 20-year mean (thick curves) and meridional standard deviation (thin vertical bars) of the Kuroshio axis position based on the reanalysis data (red) and observation-based MIRC data (blue).(b) Time series of the difference in the Kuroshio axis position in the 126.5˚E -140.5˚Einterval between the reanalysis and MIRC (QBOC) data.

Figure 7 .
Figure 7. Year-to-year variability of the simulated Kuroshio Extension axis from 1993 to 2012.

Figure 8 .
Figure 8.(a) Path length (141˚E -153˚E) and (b) mean meridional position (141˚E -158˚E) of the Kuroshio Extension (KE).Blue lines indicate the estimates of [69], and red lines show the reanalysis data.

Figure 9 .
Figure 9. Comparisons of the southernmost latitudes of the (a) first and (b) second Oyashio intrusion and (c) the easternmost longitude of the Tsugaru Warm Current water (TWCW) between reanalysis (red lines) and observation data (blue lines).These indices were determined from about 10-day (reanalysis data) and monthly (observation data interpolated by Gaussian filters) mean temperature maps at a depth of 100 m.The Oyashio domain was initially specified on a 100-m temperature map as a cold water area surrounded by 5˚C isotherms and distributed continuously off the eastern coast of Hokkaido.The boundary of the Oyashio domain typically exhibits tongue-shaped southward intrusions.The intrusions closest and second closest to the Japanese coast are defined as the first and second Oyashio intrusion, respectively.The easternmost longitude of the Tsugaru Warm Current was determined as the easternmost position of isotherms greater than 6˚C at a depth of 100 m that pass through the Tsugaru Strait (Figure 1) off southern Hokkaido into the North Pacific.

Figure 10 .
Figure 10.(a) Number of 60-day drifter data subsets in each grid cell (resolution, 1/2 degree × 1/2 degree).The data subsets were counted in the grid cell in which the initial drifter position of each data subset was located.(b) Drifter trajectory lengths integrated over 60 days.Lengths greater than 2000 km are shown by red squares, and those less than 1000 km are shown by blue squares.(c) RMSEs estimated from the position difference between the drifter and the BPP.Red squares indicate RMSEs greater than 200 km, and blue squares indicate those less than 100 km.(d) Ratios of RMSE (c) to integrated trajectory length (b).Red squares indicate ratios greater than 0.2, and blue squares indicate those less than 0.1.In (b) to (d), all values are averaged over the grid cell where the initial drifter position of each data subset was located.

Figure 11 .
Figure 11.Scatter plots of the RMSE between the drifter and BPP positions versus the trajectory distance of the drifters integrated over 60 days.Scatter plots were displayed as frequency (unit of percentage) in each class (i.e., data-point number in each class divided by total data number).

FRA-ROMS have
been developed for simulating coastal waters around Japan.A resolution of 1/50 degree has been selected for studying coastal-offshore interactions in the shelf-slope region[29] [30][75].The simulated outputs have been applied, for example, in scientific projects aimed at promoting recovery from the unparalleled damages to fisheries caused by the 11 March 2011 tsunami.For use in small, semi-enclosed bays such as Mikawa Bay and the Yatsushiro Sea in Japan, super-high-resolution ocean models with a grid size of less than 300 m have been developed[75] [90] [91][92] to simulate the transport of jellyfish, including vertical movements, and the occurrence of red tides by using FRA-ROMS in combination with non-passive particle-tracking models.Hasegawa et al. (personal comm.) have also developed a 1/90-degree wave-ocean-sediment coupled model that uses FRA-ROMS reanalysis data to improve understanding of processes affecting the transport of the radioactive materials that were accidentally leaked from the Fukushima Nuclear Power station after the 2011 mega-earthquake and attached to sediments on the sea bottom.

Forestry
and Fisheries Research Information Technology Center (AFFRIT).FRA-ROMS was developed by the physical oceanographers of the FRA with partial support from the Fisheries Agency, but the contents of this study do not necessarily reflect the views of the Fisheries Agency.This work was also financially supported by Japan Society for the Promotion of Science KAKENHI Grant Number 23740365, a Kofu-kin research project by the FRA, and by the "The Study of Kuroshio Ecosystem Dynamics for Sustainable Fisheries (SKED)" project of the Ministry of Education, Culture, Sports, Science and Technology.

Table 1 .
List of assimilation data.Gridded data by optimal interpolation to merge satellite and in situ SST by NOAA https://podaac.jpl.nasa.gov/dataset/NCDC-L4LRblend-GLOB-AVHRR_OI2) OSTIA (daily, 0.05 degree, after 2007) Gridded data by optimal interpolation to merge satellite and in situ SST by UK Met Office http://ghrsst-pp.metoffice.com/pages/latest_analysis/ostia.html 7 days by assimilating available near-real-time data.The output is updated every week and posted on the FRA-ROMS website (http://fm.dc.affrc.go.jp/fra-roms/index.html).Meanwhile, the latter is used for reanalysis experiments aimed at accurate reproduction of past oceanographic conditions by assimilation of high-quality delayed-mode data, which have been massively accumulated.The experimental version can be used more effectively for examination of linkages between oceanographic and fisheries/fish conditions in past time, whereas the operational version is more effective for examining linkages in near real time.The reanalysis data presented here were produced by the experimental version and are expected to be more accurate than the products of the operational version.
FRA has employed FRA-ROMS reanalysis and prediction products in multiple projects, mainly ones sponsored by the Fisheries Agency.Reanalysis and prediction data are routinely utilized for forecasting fishing and oceanographic conditions in the Kuroshio-Oyashio region (http://abchan.fra.go.jp/index2.html).Predictability of an operational version of FRA-ROMS will be demonstrated in near future work.A variety of IBMs based on the reanalysis data produced by the Fisheries Agency's research