A recruitment forecasting model for the Pacific stock of the Japanese sardine (Sardinops melanostictus) that does not assume density-dependent effects

This study developed a recruitment forecasting model based on a new concept of the stock recruitment relationship. No density-dependent effect in the relationship was assumed in the model, which showed that fluctuations in recruitment and spawning stock biomass of Japanese sardine in the northwestern Pacific can be explained mainly by environmental factors and the effects of fishing. The February Arctic Oscillation (AO) and sea surface temperature over the southern area of the Kuroshio Extension (30 35 ̊N and 145 180 ̊E; KEST) were used as the environmental factors. The recruitment forecasting model is proposed:   ˆ t t


INTRODUCTION
Among Japanese fisheries, the Japanese sardine (Sar-dinops melanostictus) is one of the most important species in the northwestern Pacific.In 1997 the Japanese Government introduced a total allowable catch (TAC) system for seven species including sardine.Sardine numbers began to increase during the 1960s and were very high in the 1980s, but decreased markedly in the early 1990s and have remained low.Since the TAC system was introduced, many fishermen who harvest sardine have opposed the schemes for management of the depleted populations.The government of Japan insists that a reduction in the TAC is necessary to rehabilitate the sardine population, but the fishermen insist that a reduction is not appropriate because the fluctuations in sardine population abundance are mainly caused by environmental factors.
To establish management schemes that are acceptable to fishermen but meet government management objectives, fisheries scientists need new models of sardine population abundance that can explain the extremely large fluctuations in recruitment and the spawning stock biomass.
Many studies have investigated the relationship between population fluctuations and environmental factors [1][2][3][4].For instance, Noto [5], and Noto and Yasuda [6] studied the effect of environmental conditions, as represented by the sea surface temperature in the southern area of the Kuroshio Extension during winter.They raised the possibility that environmental factors determine annual changes in the survival of sardine, and also influence long-term fluctuations in the population and its distribution and migration patterns.Sunami [7] noted that water temperature and variations in the abundance of zooplankton (a sardine food source) influence the early growth stages of sardines (up to 1-year-old), which determines their survival rate.
Although most focus has been on the role of environmental factors, density-dependent effects are also be-lieved to be important in explaining fluctuations in sardine stock abundance.For instance, Chen and Irvine [8] developed a stock-recruitment relationship model that incorporates environmental effects and performs better than the traditional Ricker model, on which it is based.Yatsu et al. [9] investigated the reproductive success of sardine with respect to water temperature, zooplankton density, intra-species competition, and density-dependent effects.Ishida et al. [10] adapted modified Ricker models that incorporate environmental factors, and proposed new stock-recruitment models for the Japanese sardine.However, all these models can be categorized as extended Ricker models because they stress the importance of density-dependent effects as well as environmental factors and fishing in the stock recruitment relationship.
Sakuramoto and Suzuki [11] showed even where there is no density-dependent effect in the stock-recruitment relationship, there is a high likelihood of detecting an apparent density-dependent effect on response to observed and/or process errors; the latter is considered to occur mainly because of environmental changes.They also criticized the a priori assumption of a density-dependent effect on the sardine stock recruitment relationship, because such an assumption can lead to an underestimate of the effects of environmental factors.They concluded that when we analyze the stock recruitment relationship data of the northwestern Pacific sardine stock, no density-dependent effect can be detected, and a proportional model should be used for assessment of the stock recruitment relationship for Japanese sardine.
The purpose of this study was to present a recruitment forecasting model based on a stock recruitment relationship where no density-dependent effect is assumed, and to show that large fluctuations in recruitment and the spawning stock biomass of sardine can be explained only by environmental factors and the effects of fishing.

Data
The catch number, average weight, maturity rate, and the natural mortality coefficient were recorded by age and year for the northwestern Pacific sardine stock from 1976 to 2012 [12].The number of fish and the fishing mortality coefficients by age in each year were also used from the report.[12].Data on the Arctic Oscillation (AO) for February was obtained from the NOAA Climate Prediction Center [13], and the sea surface temperature over the southern area of the Kuroshio Extension for 30 -35˚N and 145 -180˚E (KEST) in February was obtained from Ishii [14].

Recruitment Forecasting Model
The recruitment forecasting models developed com-prised two parts: one concerning the role of the stock recruitment relationship f(S), and the other the role of environmental factors, g(*), as shown in Equation (1): , , , where R t is the recruitment at the beginning of fishing period in year t, S t is the spawning stock biomass at the beginning of the fishing period in year t, and x where, t  is an error term having a log normal distribution with mean 0 and variance 2 t  .The parameters were estimated by regression analysis using Equation ( 3), and the statistically significant partial regression coefficients were chosen.
 
The recruitment number at the beginning of the fishing period in the following year, = , was calculated using Equation ( 4): The number of fish from age 1 to 4 in the following year, , i t N  , was calculated using Equation ( 5), when t = 1976, the initial values, N i,t , (i = 1, 2, , 5+) are used.
The number of fish at age 5 and older in the following year ( 5 , 1 t N    ) was calculated using Equation ( 6): where M i,t , and F i,t are the natural mortality coefficient at age i in year t and the fishing mortality coefficient at age i in year t, respectively.
The spawning stock biomass in the following year (t + 1) was calculated using Equation ( 7): where m i,t+1 and w i,t+1 are the maturity rate at age i in year t + 1 and the mean weight at age i in year t + 1, respectively.This process was repeated until t + 1 = 2012.
The partial regression coefficients estimated by the least square method for Equation (3) were slightly modified to minimize the sum of square value (SS) defined by Equation ( 8): where i and i are the referred [12] and forecast recruitment, respectively.R R 2.4.

Forecasting the Recruitment in 2013
When the values for the AO and the KEST for February 2013 are available it will be possible to predict the recruitment in 2013 using the forecasting model proposed here, provided that the maturity rate by age and the mean weight by age in the year 2013 are the same as in 2012.The recruitment in 2013 can be calculated by VPA after the middle of the following year 2014, when the catch at age data, maturity rate, and mean weight by age for 2013 will have been compiled.Therefore, the recruitment in 2013 can forecast more than one year earlier by using this model.

Parameters Estimated by Regression Analysis
A recruitment forecasting model was chosen for which all the partial regression coefficients were statistically significant.The parameters and the statistically significant probabilities are shown in Table 1.The resulting model is shown in Equation ( 9): The partial regression parameters were then modified to reduce the SS value calculated by Equation (8).More where a i is the partial regression parameters estimated using the least squares method and h i denotes the interval given for each parameter.The modified partial regression coefficients are also shown in Table 1.The SS value was reduced from 5.22  10 23 to 6.85  10 22 using this approach.

Results of Forecasts
Figure 1(a) shows the referred recruitment [12] and that forecast using the model.While large differences were found for 1980, 1983, 1985, 1988, and 1991, generally the pattern of the variation was reproduced well.In particular, the dramatic decline in recruitment from 1986 to 1989 was reproduced well, regardless of the assumption of a density-dependent effect in the stock recruitment relationship.Figure 1(b) shows the case where the scale in the Y-axis was enlarged after 1995.The referred recruitment decreased until 2004, after 2005 there was an apparent slight increase, and since 2010 the level of recruitment has been remarkably high.The forecast pattern was quite different, with recruitment continuing to decrease after 2004.The reasons for these differences are discussed below.
Figure 2(a) shows the spawning stock biomass calculated using Equation (7).Prior to 1988 the forecast spawning stock biomass was lower than the referred biomass.This is because prior to 1988 the forecast recruitment was less than the referred recruitment, except in 1983 and 1984.However, the increasing tendency referred was reproduced well.On the contrary, after 1989 the decreasing tendency referred was reproduced significantly well except in 1993 and 1994.differences are also discussed below.
Figure 3 shows the forecast and referred stock recruitment relationships.The variables on both the x and y axes were estimated, and some of the forecast and referred data points do not coincide well.However, the trajectory showing a loop rotating clockwise was reproduced well.
Figure 4 shows the forecast and referred recruitment per spawning stock biomass (RPS) [12].The fluctuations in RPS were large, but the trends were reproduced well except for the years 1977, 1980, 1982, 1991, 2005, and 2008.In 2010 the forecast value estimated by this model was extremely large relative to the referred one.

Validity of the F Estimates after 2001
This section considers the validity of the value of the      ).This reduction in the recruitment and spawning stock biomass may not have been a consequence of overfishing.Because the forecasted catches in ages 0 and 1 were low after 2000.Figure 6 shows the fishing mortality coefficient trajectories by age [12].Figure 6 shows that after 2001 the fishing mortality coefficients for ages 0 and 1 were extremely large compared with those in the previous years.
One explanation for the extremely low level of recruitment and spawning stock biomass after 2005 is that the values were erroneously forecast in response to the large fishing mortality coefficients used in the simulation.Thus, when the sardine population was extremely low the fishermen tended to target species for which the biomass was greater (in this case mackerel or Japanese anchovy).The VPA method did not incorporate qualitative changes in fishing effort in estimating the fishing mortality coefficient.Therefore, it is highly likely that the fishing mortality coefficient was overestimated during the period when the population abundance was extremely low.and up to age 1 in 2005.If the fishing mortality coefficient after 2001 was much lower, the spawning stock biomass would increase, and thus the recruitment would also be expected to increase.
Simulations were conducted when the smaller fishing mortality coefficient, 100d% of fishing mortality coefficient was used.The parameter d was determined so as to make the sum of square value (SS2) the smallest, which is defined by Equation (10).
SS2 was smallest when d was set at 0.50.In Figure 6, it was also shown that the fishing mortality coefficient after 2001 was half of that referred.Figures 7-9 show  the forecast recruitment, spawning stock biomass, and stock recruitment relationship (respectively) when d was set at 0.50.The forecast values were markedly improved using this approach.
Figures 10(a) and (b) show a comparison of the forecast and referred recruitment and spawning stock biomass.The determination coefficients were 0.755 and 0.942, respectively.

Forecasting the Recruitment in 2013
When values of the AO and the KEST for February 2013 are available it will be possible to forecast the re- cruitment in 2013 using the forecasting model proposed here, provided that the maturity rate by age and mean weight by age for 2013 are the same as in 2012.The recruitment in 2013 was forecast to be 61.3  10 8 , which is 35.6% less than in 2012.Corresponding to this value, the spawning stock biomass in 2013 was estimated at 346,000 tons, which is 32.3% less than in 2012.However, the former will vary depending on the maturity rate and mean weight by age in 2013.Although these are estimated values it is important to estimate values that can be confirmed in the following year, as this can aid assessment of the validity of the forecasting model, and provide important information for fisheries resource managers.

DISCUSSION
A noteworthy finding of the present study was that large fluctuations in recruitment and the spawning stock biomass were reproduced remarkably well by the proposed model, which did not assume a density-dependent effect on the stock recruitment relationship.As shown in Figure 3, a large loop emerged in the stock recruitment relationship for the northwestern Japanese sardine stock; however, it is difficult to explain how the loop emerged.Although, this loop is never produced with density-dependent effects, the recruitment forecasting model proposed in this study can reproduce a similar loop in the stock-recruitment relationship in response to the fluctuations of environmental factors.
Shimoyama [16] forecasted the recruitment and spawning stock biomass using the models, f S , she used not only a proportional model but also the models that had the density-dependent effects, and showed the similar results regardless of what kinds of the stock-recruitment models assumed.However, she only used the two independent variables, AO t and KEST t in the function g(AO t , KEST t ).Therefore, the forecasting did not reproduce the recruitment well.Especially the steep decline from 1986 to 1989 could not reproduce well.In this model, the performance was improved much better than the results of Shimoyama [16].
In this study, I did not compare the models between the cases when density-dependent effects were assumed or not.Because, if the results for the case when the density-dependent effects were assumed was little better than the case when they did not assume, as Sakuramoto and Suzuki [11] showed that it did not necessarily indicate the validity of density-dependent models.It just indicated that other environmental factors that would control the recruitments should be incorporated.
Sakuramoto [17] showed that the clockwise loop shown in the stock-recruitment relationship of the Japanese sardine could be explained by the mechanism, which constructed two assumptions.One is the recruitment is proportional to the spawning stock biomass and another is the spawning stock biomass fluctuates cyclically in response to environmental fluctuations.The simulation shown in this study supported the validity of the assumptions.That is, essentially the recruitment is proportional to the spawning stock biomass, however, the recruitments were changed in response to environmental fluctuations, and the spawning stock biomass fluctuates cyclically, and then the clockwise loop emerged, which is never explained by the density-dependent effect on the stock-recruitment relationship.
Shimoyama et al. [18], Sakuramoto et al. [15], and Sakuramoto and Suzuki [11] noted that no density dependent effect could be detected in data related to the stock recruitment relationship for the Japanese sardine, and suggested that a proportional model be used to describe the stock recruitment relationship for this species.
However, the criticism always noted for this conclusion is that in this case the population will increase infinitely or deplete to zero when the environmental conditions are constant.Therefore, assuming a proportional model for the stock recruitment relationship of Japanese sardine appears inappropriate.However, this criticism is not logical because environmental conditions never remain constant in the real world.Therefore, even if an unrealistic conclusion is derived from an unrealistic assumption (e.g., that environmental conditions are constant), this does not imply that a proportional model is inappropriate.As shown in this study, as environmental conditions fluctuate the population never increased infinitely or depleted to zero, even when the stock recruitment relationship was expressed using a proportional model.
Data for the Pacific decadal oscillation index (PDO) were also used as an environmental condition, but in this model the partial regression coefficients related to the PDO were not statistically significant.In future, if other environmental conditions are found that explain the fluctuations in recruitment more accurately, the forecasting model can be further improved.
The sensitivity tests for the F values after 2001 indicated that there was a possibility that F values were overestimated when the population abundance is extremely low.The actual F values may be almost the half of those of estimated [12].The sensitivity tests also indicated that fishing is an important factor controlling the population trajectories, especially when the population biomass is low and environmental conditions are unfavorable for the sardine population.
The earthquake and tsunami after the earthquake occurred on March 11, 2011 must have changed extremely the environmental conditions along the Pacific coast of Japan.These changes might affect the survival of the fish recruited after 2012.The forecasting is done using only two environmental factors, AO and KEST in February, therefore, if other environmental conditions, which affect the recruitment, have extremely changed, of course, the performance of the forecasting will not be well.

Figure 2 (Figure 1 .
Figure 1(a)shows the referred recruitment[12] and that forecast using the model.While large differences were found for1980, 1983, 1985, 1988, and 1991, generally the pattern of the variation was reproduced well.In particular, the dramatic decline in recruitment from 1986 to 1989 was reproduced well, regardless of the assumption of a density-dependent effect in the stock recruitment relationship.Figure1(b)shows the case where the scale in the Y-axis was enlarged after 1995.The referred recruitment decreased until 2004, after 2005 there was an apparent slight increase, and since 2010 the level of recruitment has been remarkably high.The forecast pattern was quite different, with recruitment continuing to decrease after 2004.The reasons for these differences are discussed below.Figure2(a) shows the spawning stock biomass calculated using Equation(7).Prior to 1988 the forecast spawning stock biomass was lower than the referred biomass.This is because prior to 1988 the forecast recruitment was less than the referred recruitment, except in 1983 and 1984.However, the increasing tendency referred was reproduced well.On the contrary, after 1989 the decreasing tendency referred was reproduced significantly well except in 1993 and 1994.Figure2(b) shows the case where the scale in the Y-axis was enlarged after 1995.The referred spawning stock biomass decreased until 2005, was stable until 2008, and subsequently increased markedly.The forecast pattern was substantially different, with the spawning stock biomass continuing to decrease after 2005, following the same pattern as recruitment.The reasons for these

Figure 2 .
Figure 2. Trajectories in spawning stock biomass forecasted and referred of Japanese sardine (a); Trajectories shown in enlarged scale (b).Black and red lines indicate the referred and forecasted spawning stock biomass, respectively.

Figure 3 .
Figure 3. Stock-recruitment relationship for Japanese sardine in the northwestern Pacific from 1976 to 2012.Black and red indicate the referred and forecasted relationship, respectively.Figures denote the year.
Figure 3. Stock-recruitment relationship for Japanese sardine in the northwestern Pacific from 1976 to 2012.Black and red indicate the referred and forecasted relationship, respectively.Figures denote the year.

Figure 4 .
Figure 4. Recruitment per spawning stock biomass (RPS) for Japanese sardine in the northwestern Pacific from 1976 to 2012.Black and red indicate the referred and forecasted RPS, respectively.fishing mortality coefficient, F a,t .The results shown in Figures 1(b) and 2(b) indicate that the forecast and referred values and trends were significantly different.Figure 5 shows the referred and forecast catch trajectories by age.For ages 0 after 2000, the forecast catches were extremely low relative to the referred catches, and the forecast recruitment and spawning stock biomass after 2005 were much less than the referred values (Figures 1(b) and 2(b)).This reduction in the recruitment and spawning stock biomass may not have been a consequence of overfishing.Because the forecasted catches in ages 0 and 1 were low after 2000.Figure6shows the fishing mortality coefficient trajectories by age[12].Figure6shows that after 2001 the fishing mortality coefficients for ages 0 and 1 were extremely large compared with those in the previous years.One explanation for the extremely low level of recruitment and spawning stock biomass after 2005 is that the values were erroneously forecast in response to the large fishing mortality coefficients used in the simulation.Thus, when the sardine population was extremely low the fishermen tended to target species for which the biomass was greater (in this case mackerel or Japanese anchovy).The VPA method did not incorporate qualitative changes in fishing effort in estimating the fishing mortality coefficient.Therefore, it is highly likely that the fishing mortality coefficient was overestimated during the period when the population abundance was extremely low.Figure 1(b) indicates a large difference between forecasted and referred recruitment trends since 2005.The spawning stock biomass was based on mature sardine (1 -5 years old).Therefore, the spawning stock biomass in 2005 was determined based on sardines of age 1 in 2001,

Figure 5
shows the referred and forecast catch trajectories by age.For ages 0 after 2000, the forecast catches were extremely low relative to the referred catches, and the forecast recruitment and spawning stock biomass after 2005 were much less than the referred values (Figures 1(b) and 2(b)

Figure 1 (
b) indicates a large difference between forecasted and referred recruitment trends since 2005.The spawning stock biomass was based on mature sardine (1 -5 years old).Therefore, the spawning stock biomass in 2005 was determined based on sardines of age 1 in 2001,

Figure 5 .
Figure 5. Trajectories of catch by age for Japanese sardine in the northwestern Pacific.Black and red indicate the referred and forecasted, respectively.

Figure 6 .
Figure 6.Trajectories of fishing mortality coefficients by age for Japanese sardine in the northwestern Pacific.Black and red are the values referred and the half of the referred values from 2001, respectively.

Figure 7 .
Figure 7. Trajectories in recruitment forecasted (red line) and referred (black line) when the fishing mortality coefficient by age after 2001 was set at half of the referred values.

Figure 8 .
Figure 8. Trajectories in spawning stock biomass forecasted (red line) and referred (black line) when the fishing mortality coefficient by age after 2001 was set at half of the referred values.

Figure 9 .Figure 10 .
Figure 9. Stock-recruitment relationship for Japanese sardine in the northwestern Pacific from 1996 to 2012.Black and red lines indicate the referred and forecasted when the fishing mortality coefficient by age after 2001 were set at half of the referred values.Figure denotes the year.

Table 1 .
Estimated partial regression coefficeints and modified parameters used in the simulations.