Skipjack movement and fisheries interaction in the western Pacific

John R. Sibert
Pelagic Fisheries Research Program
Joint Institute of Marine and Atmospheric Research
University of Hawaii
Honolulu, Hawaii 96822, USA

John Hampton
Oceanic Fisheries Programme
South Pacific Commission
B.P. D5, Noumea Cedex, New Caledonia

David A. Fournier
Otter Research Ltd.
P. O. Box 265, Station A
Nanaimo, B.C. V9R 5K9, Canada

ABSTRACT

An advection-diffusion model in which movement parameters are variable in both time and space is applied to the analysis of skipjack tag and recapture data. Two data sets spanning the periods 1977-1982 and 1989-1992 are analysed. The results indicate that skipjack movement is highly variable in time and space with very little coherent pattern. The estimated movement pattern mortality and catchability coefficients, and observed patterns of fishing effort for the period 1989-1992, are used in a general model of the skipjack fisheries in the western Pacific. The model predicts that the pole and line yield may have declined by as much as 20% due to the activities of the purse seine fishery. We conclude, however, that a change of this magnitude would be difficult to detect in actual recorded catches.

1. INTRODUCTION

A quantitative method to describe large scale movement of fish is necessary to predict the effects of a fishery operating in one place on the yield of a fishery operating in another place. Such descriptions comprise movement patterns. The purpose of this paper is to examine movement patterns of skipjack in the western Pacific ocean and to attempt to determine how these patterns may change over time.

A model of the large scale net movements of tagged skipjack is used to estimate movement patterns. In this model, movement is partitioned into directed and dispersive components. The directed component refers to net displacement in a unique direction. The dispersive component, also referred to as diffusion, is an essentially random process. Description of fish movement as "random" poses difficulties for behaviorists who claim with justification that the behavior of individual animals is not random. Movements of individual skipjack on a small scale may have well determined causes such as pursuit of prey, avoidance of predators, or preference for water with certain properties. Although individuals may never move at random on a small scale, there may be considerable variation in the population, and descriptions of large scale movements of populations should incorporate this variability. The relative importance of directed and dispersive movement in determining the distribution of skipjack populations at different scales will also be analyzed. Finally the potential impact of expanding purse seine fisheries on pole and line fisheries for skipjack will be discussed.

2. METHODS: ESTIMATION OF MOVEMENT PATTERN

Movement patterns were estimated from tag and recapture data using a population dynamics model that explicitly incorporates movement. A preliminary version of the model was presented to the First FAO Expert Consultation on Interactions of Pacific Tuna Fisheries (Sibert and Fournier, 1994). Additional model structure, including spatially variable diffusion, monthly changes in fishing mortality, and fleet-specific tag reporting rates, has been added since that presentation. A mathematical summary of the model is presented in the Appendix.

The model partitions the change in population density of tagged skipjack at a point in the ocean between movement and mortality. Movement is modeled using the generalized advection-diffusion partial differential equation; see equation (2) in the Appendix. While this equation is usually considered to describe transport of a solute in a moving fluid, it is also the limiting case of a biased random walk and can be used to describe movements of animals (Okubo, 1980). This formulation recognizes two types of behavior: directed movement represented by the vector (u, v) (also called advection) and random movement represented by D (also called diffusion). Mortality is modeled as the sum of mortality due to fishing and natural mortality. Mortality due to fishing is proportional to the nominal fishing effort reported by each fleet operating in the model area and thus varies over time and space. Natural mortality is assumed to be constant throughout the model area and is properly interpreted as mortality not attributable to either fishing or to movement. The area to which the model is applied is closed so that there is no possibility that emigration from the model area can be confused with mortality. In the absence of mortality, the total number of tagged skipjack in the model area remains constant however large the variation may be in the number at any single point.

The southern boundary of the model area is 25° S latitude, the northern boundary is 10° N latitude, the western boundary is 130° E longitude, and the eastern boundary is 170° W longitude. This model area covers approximately 7.5 x 10⁶ nm². Regions within the model area are defined as geographical subdivisions of the model area, and Seasons are defined periods of time over which the movement parameters (u, v, D) are constant. The model area and the ten regions within it used to analyze the data are shown in Figure 1. Skipjack movements are considered to be seasonal in parts of their range. Movements for two seasons were estimated, Season 1 (from April through September) and Season 2 (from October through March). Selection of seasons is largely arbitrary. Anecdotal evidence suggests that skipjack movements may be seasonal and the two seasons selected correspond to the southern hemisphere cold and warm seasons, respectively.

The borders of the model area and are arbitrary. Simulation results reported by Sibert and Fournier (1994) indicate that parameter estimates can be inaccurate for regions in which no tags were recovered. Therefore, the model area and regions shown in Figure 1 were established so that there were recaptures in all regions during each season.

The advection-diffusion equation treats a quantity as continuous in time and space. Numerical solutions of the equation on a computer require that it be discretized in time and space. Although, in principle, the resolution of the numerical model is limited only by the resolution of the data, in practice, spatial resolution is constrained by computer memory limitations and temporal resolution is constrained by the time one is willing to wait for a sufficiently accurate solution. The spatial resolution used in the present analysis is 60 nm and the temporal resolution of the data is one month.

Figure 1. Map of southwest Pacific ocean showing the model area and the locations of regions within it. The model resolution of the coastlines of Australia and Papua New Guinea is also indicated.

The model requires two types of data, tag release and recapture data, and fishing effort data. All data are aggregated to one degree square by one month strata. A tag cohort is defined as the sum of all tagged skipjack released in a one degree geographic square during a month. Aggregation of the tag release data to this level means that all tags released during one month are assumed to mix instantaneously throughout the one degree geographic square in which they were released.

Tag and recapture data from two different South Pacific Commission (SPC) studies were used, the Skipjack Survey and Assessment Programme (SSAP), 1977-1982, and the Regional Tuna Tagging Programme (RTTP), 1989-1992. The SSAP released 87,941 tagged skipjack in 102 cohorts in the model area, and the RTTP released 77,962 tagged skipjack in 237 cohorts in the same area. Recaptured tags are only useful in this analysis when reported with recapture position, date and fleet. There were 4,018 recaptured tags reported to the SSAP and 6,258 to the RTTP within the model area with useable recapture information.

Fishing effort data were obtained from a variety of sources including daily fishing logs reported to Pacific island states and forwarded to the SPC, published data summaries, and data supplied to the SPC by fishing nations. Fishing effort data for the model are nominal fishing effort expressed in boat-days fished per month in each 1° x 1° geographic square. Gear types included are pole and line and purse seine, the two principal gear types that exploit skipjack in the western Pacific. Fishing nations include Fiji, Japan, Kiribati, Korea, Papua New Guinea, Philippines, Solomon Islands, Taiwan, and USA.

The proportion of tags recaptured in the fishery and returned to the SPC with useful information is defined as the reporting rate (see Appendix). In the absence of fleet specific estimates, a reporting rate of 1.0 is assumed for all fleets returning tags to the SSAP. Reporting rates for the RTTP were estimated for each fleet on the basis of tag seeding results and other information (SPC. 1994). Reporting rates were multiplied by (i) the proportions of tags returned by each fleet that had accompanying data of sufficient resolution to be included in the analysis and by (ii) the global value, 0.923, the proportion of all tag returns in the model area accounted for by the fleets included in this analysis. The resulting overall fleet-specific reporting rates are given in Table 1.

Table 1. Reporting rates used for estimation of movement patterns.

Fleet	SSAP	RTTP
Fiji pole and line	1.0	0.826
Japan pole and line	1.0	0.376
Japan purse seine	1.0	0.360
Kiribiti pole and line	-	0.646
Korea purse seine	-	0.026
Papua New Guinea pole and line	1.0	-
Philipines purse seine	-	0.669
Solomon Islands pole and line	1.0	0.660
Solomon Islands purse seine	-	0.713
Taiwan purse seine	-	0.016
United States purse seine	-	0.509

3. RESULTS: MOVEMENT PATTERN

The movement pattern estimated from the SSAP data is shown in Figure 2A. The most prominent feature of this pattern is the variability between regions and seasons in both directed and dispersive components of movement. In spite of this variability there are some common features. There is a consistent eastward tendency in the directed component in the regions between 5° and 10° N, most prominent in Region 1. Rates of roth directed and dispersive movements estimated for Region 10 are rather low, particularly in Season 2.

The persistence of the movement pattern in the SSAP data was tested by dividing the tag releases into two time periods: (A) October 1977 through September 1978 (47 cohorts released) and (B) October 1978 through August 1980 (55 cohorts released). The mortality parameters (M,Q_i) were fixed at the values estimated for the overall SSAP data set and the movement parameters (u, v, D) were estimated separately for the two time periods. The value of the negative log likelihood function (equation 8 in Appendix A) for the overall data set is 11,696.03 and there are 60 estimated movement parameters. The combined negative log likelihood for the two time periods is 11,620.49 and there are 120 movement parameters (60 for each time period). The model with separate movement parameters thus fits the data better (lower negative log likelihood) than the overall model. The significance of this improvement in fit can be tested using a likelihood ratio test (Brownlee, 1965) in which twice the difference between the negative log likelihood for the two models is interpreted as a x² variable with degrees of freedom equal to the difference in numbers of parameters in the two models. In this case the x² value is 75.55, the degrees of freedom is 60, and the improvement in fit of the 120 parameter two period model is significant with 0.05 < a < 0.10. Therefore there are significant differences between the movement patterns for the two periods.

The role of dispersive movements was tested by fitting the SSAP data set to a model in which dispersive movements were constrained to be zero: that is, all movement was assumed to be directed. This model produced a negative log likelihood value of 12,759.00 with 26 parameters compared to the complete model with a negative log likelihood values of 11,696.03 with 66 parameters. The improvement in fit of the complete model is significant with a < 0.001. Similarly a model was fit in which the directed movements were constrained to be zero; that is, all movement is assumed to be dispersive. Again the complete model fits the data significantly better with a < 0.001. These tests indicate that both dispersive and directed movements are required to describe the movements of tagged skipjack during the SSAP.

Figure 2. Skipjack movement patterns estimated from SSAP (A.) and RTTP (B.) tag returns. The number in the upper left corner of each region is the region number. The number (in parentheses) in the lower right comer of each region is the number of tags returned from that region during a season. The area of the circles is proportional to the rate of dispersive movement. The length of the arrows is proportional to the directed movement rate in the direction indicated. The radius of the circles and the length of the arrows are drawn to the same scale as indicated at the right.

A. SSAP - Season 1 (4..9)

A. SSAP - Season 2 (10..3)

B. RTTP - Season 1 (4..9)

B. RTTP - Season 2 (10..3)

The overall movement pattern estimated for the RTTP is shown in Figure 2B. Although the details are quite different, the degree of regional and seasonal variability is similar to that estimated for the SSAP. The eastward directed movement appears to have shifted to the band between the equator and 5° N (Region 3, Season 1) and the lack of movement in Region 10 during Season 2 is still detectable. The year to year variability was tested by dividing the RTTP tag releases into two groups: (A) those released between July 1989 and October 1990 (86 cohorts) and (B) those released between February 1991 and December 1992 (151 cohorts). The difference in movement pattern between the two years is significant with a < 0.001.

The relative importance of dispersive and directed movement components can be expressed by a dimensionless ichthyokinematic ratio, R, computed as

where is the magnitude of the resultant vector formed by the u, v components of directed movement and It. is a characteristic movement length. Values of R greater than one indicate that directed movements are more important than dispersive movements over distances specified by the characteristic length. This ratio is similar to the Peclet number (Huyakorn and Pinder, 1983; Boudreau, 1986) and can be considered as a Reynolds number for the movement of fish. Values of R less than one indicate that dispersive movements are more important than directed movements over a distance specified by the characteristic length. Table 2 shows values of R computed for three values of l characteristic of the model resolution (60 nm), a scale typical of an Exclusive Economic Zone (200 nm), and a longer distance (1,000 nm). At the model scale, dispersive movements are dominant; as the scale of movement (l) is increased, directed movements become dominant.

Movement and mortality are integral components of population dynamics. The estimated catchability and natural mortality coefficients along with estimates of their standard deviations are given in Table 3. Natural mortality was estimated to be approximately 0.07^-mo for the SSAP and near zero for the RTTP. Catchability coefficients vary more between gear types than between the SSAP and the RTTP data sets. Catchabilities for purse seine fleets range from 0.001 to 0.05 while catchabilities for pole and line fleets are approximately an order of magnitude lower.

4. RESULTS: FISHERY INTERACTION

The impact of expanding purse seine fisheries on pole and line fisheries is one of the major preoccupations in the analysis of fishery interaction. Analysis of the RTTP data provides estimates of fishing mortality by four pole and line, six purse seine fleets and other parameters necessary to model regionally expanding fisheries. The population model described in the Appendix was used to simulate the skipjack population and fisheries that exploit it. The average monthly distribution of fishing effort for each of the 10 fleets used in the RTTP analysis was calculated by averaging over years from 1989 through 1993. These averages were used to simulate the current distribution of fishing by these fleets.

Table 2. Ichthyokinematic ratios for SSAP and RTTP analyses. Movement parameters u, v (in nautical miles per day) and D (in square nautical miles per month) are shown for each region and season. Dimensionless dispersion indices are calculated for three different scales of movement: 60, 200 and 1,000 nautical miles.

Season	Region	Parameters			Ichthyokinematic Ratio
		u	V	D	60	200	1000
		SSAP Analysis
1	1	4.67	-0.15	7097	1.19	3.95	19.76
	2	0.25	-4.16	14811	0.51	1.69	8.43
	3	-5.00	-0.36	84576	0.11	0.36	1.78
	4	-0.46	-0.41	15666	0.07	0.24	1.19
	5	-1.00	0.28	250	7.45	24.82	124.11
	6	-1.09	2.13	5431	0.79	2.65	13.23
	7	-0.26	-0.91	0	>>1	>>1	>>1
	8	-17.77	1.00	56448	0.57	1.89	9.46
	9	-0.03	0.79	11658	0.12	0.41	2.03
	10	0.95	2.39	35	131.55	438.49	2192.45
2	1	2.27	-2.07	19207	0.29	0.96	4.79
	2	4.68	0.98	3006	2.86	9.54	47.69
	3	6.82	10.05	148720	0.15	0.49	2.45
	4	-3.65	-0.72	32531	0.21	0.69	3.43
	5	-2.50	-6.47	8895	1.40	4.68	23.38
	6	6.47	1.05	8773	1.34	4.48	22.41
	7	1.82	0.03	8786	0.37	1.24	6.20
	8	1.23	-4.12	2470	3.14	10.45	52.25
	9	0.29	2.85	0	>>1	>>1	>>1
	10	-0.49	-1.25	321	7.51	25.04	125.19
		RTTP Analysis
1	1	0.03	-7.56	12269	1.11	3.70	18.50
	2	0.03	-7.85	9404	1.50	5.01	25.03
	3	19.02	2.11	16976	2.03	6.76	33.82
	4	2.50	-10.24	206440	0.09	0.31	1.53
	5	0.45	3.13	9564	0.59	1.98	9.91
	6	0.04	-3.63	22508	0.29	0.97	4.84
	7	-0.19	-0.69	1861	0.69	2.30	11.48
	8	22.09	-13.48	10	4735.31	15784.38	78921.90
	9	-1.39	-1.69	0	>>1	>>1	>>1
	10	13.51	-26.22	27588	1.92	6.42	32.08
2	1	-0.33	1.49	1	2275.39	7584.64	37923.22
	2	-5.71	4.07	0	>>1	>>1	>>1
	3	-1.95	3.43	117402	0.06	0.20	1.01
	4	0.65	3.00	98	56.13	187.09	935.47
	5	-3.11	1.04	9591	0.62	2.05	10.26
	6	2.46	3.04	18853	0.37	1.24	6.22
	7	-0.07	-1.23	1163	1.90	6.34	31.68
	8	1.75	-0.71	174947	0.02	0.06	0.32
	9	-1.93	1.78	4536	1.04	3.47	17.34
	10	-0.37	-0.01	409	1.62	5.41	27.06

Table 3. Estimates of global natural mortality (M) and catchability (Q_f) by fleet for the SSAP and RTTP analyses. Numbers set in small type to the left and right of the estimates are the estimated lower and upper 95% confidence limits of the parameter estimates. A dash indicates that the fleet was inactive.

M	SSAP			RTTP
	0.067962	0.074339	0.080950	0	0.000000	0.001713
Qf
Fiji P&L	0.004630	0.004947	0.005274	0.003932	0.004575	0.005264
Japan P&L	0.001720	0.001991	0.002280	0.000507	0.000624	0.000753
Kiribati P&L		-		0.001062	0.001385	0.001749
Papua New Guinea P&L	0.000522	0.000569	0.000617		-
Solomon Islands P&L	0.000525	0.000594	0.000668	0.000537	0.000588	0.000642
Japan PS	0.005294	0.006737	0.008340	0.015577	0.016898	0.018260
Korea PS		-		0.040285	0.049578	0.058885
Philipines PS		-		0.008769	0.009477	0.010209
Solomon Islands PS		-		0.012241	0.013382	0.014566
Taiwan PS		-		0.009570	0.012664	0.016121
United States PS		-		0.001062	0.001385	0.001749

The population was initialized at zero and allowed to grow by the process of constant monthly recruitment without fishing for 120 months until a steady state was approached (see equation (10) of the Appendix). At that point, simulated pole and line fisheries began to operate at current levels of effort. Simulated purse seine effort was increased in ten steps from zero to two times the current effort levels. Effort was held constant for 24 months after each increase to allow the system to approach equilibrium. The histories of purse seine yield and pole and line yield summed over all fleets starting in month 120 is shown in Figure 3. The stepwise growth of the purse seine yield and decline in pole and line yield can be seen. The pole and line yield relative to the pole and line yield at zero purse seine effort was computed for each fleet and for the sum of all pole and line fleets. The results are shown graphically in Figure 4. It would appear that the current pole and line yield is approximately 80% of what might have been attained at zero purse seine effort. The decline is greatest for those pole and line fleets operating nearest the large scale purse seine fleets (Japan, Kiribati) and least for those pole and line fleets operating farther away (Fiji).

Estimates of natural mortality are known to influence estimates of interaction potential (Sibert, 1984). The interaction analysis was repeated at three different levels of natural mortality by setting the numerical value of M to a constant and estimating the movement and catchability parameters. The resulting parameter estimates were used to estimate theoretical reductions in pole and line yields at current purse seine effort levels. The results are presented in Table 4. The potential for interaction generally decreases as natural mortality increases, although the trend in not monotonic, as in the case of Solomon Islands.

Figure 3. Simulated growth of the purse seine (PS) fishery from zero to two times its current effort levels and the hypothetical effect on the pole and line (PL) fishery. The duration of the simulation is 22 years.

5. DISCUSSION AND CONCLUSION

One of the striking features of these results is the large variability in both directed and dispersive components of movement between seasons and between regions (Table 2, Figure 2). The variability in the dispersive components (D) seems particularly large in relation to the directed components (u, v), but it must be remembered that the units of D contain the square term, nm², whereas the units of u, v are linear. The relative radii of the circles in Figure 2 are a better indicator of the relative magnitude of the variability in D than simple comparisons of the numerical point estimates. The parameterization of the advection-diffusion PDE used in this analysis associates the movement parameters with specific region-season strata. It is therefore difficult and probably not meaningful to attempt to interpret the parameter estimates as functions of environmental variables such as water temperature, circulation or topography.

Consistent large scale directed movement is one of the characteristics of animal migration. Skipjack are classified as "highly migratory" species, yet there appears to be little coherent pattern or temporal consistency in skipjack movement at any scale. Nevertheless, the hypothesis of totally random movement (i.e., (u, v) = 0) has been rejected. The extent to which skipjack movements can be regarded as random depends on the scale of the movement. Although it is doubtful that individual fish ever move randomly, random movement generally appears to be the dominant mode over shorter distances. Directed movements are generally dominant at larger scales, although in some region-Season strata, dispersive movements are dominant. One of the few consistencies is the eastward movement north of the equator during the SSAP period. One obvious explanation for this movement is transport in the North Equatorial Counter Current. It would be relatively simple to incorporate such water movements into the model and to test the significance of physical transport of fish by ocean currents in determining the distribution of skipjack.

Figure 4. Hypothetical decline in relative pole and line (PL) yield in response to relative increases in purse seine (PS) effort. Lines labeled FJ, SB, JP, and KI indicate yield from the pole and line fisheries of Fiji, Solomon Islands, Japan, and Kiribati, respectively. Crosses indicate the average PL yield.

Table 4. Potential percentage reduction in yield of pole and line fisheries at different levels of natural mortality (M) at current purse seine effort levels relative to pole and line yield without purse seine fisheries.

M	Potential Reduction by Fleet
	Fiji	Solomon Islands	Japan	Kiribati
0.0	10.6 %	18.4 %	28.5 %	31.7 %
0.07	5.3 %	16.4 %	23.3 %	25.5 %
0.15	2.4 %	28.4 %	21.4 %	20.5 %

The models used in this analysis are spatially disaggregated population dynamics models with a resolution of 1° square. They attempt to describe population dynamics on a local scale in contrast to the more conventional approach of attempting to describe the dynamics of the population as a whole. No assumptions are made about the integrity of the stock and whether it can be considered a discrete unit. The only assumptions about immigration and emigration are those made at the boundaries of the model area, and these boundaries can, in principle, be placed to encompass the entire population.

In a spatially disaggregated model, mortality refers to the instantaneous rate of mortality of the population at a point, not to the rate of mortality of the population as a whole. Mortality due to fishing varies over time and space from zero at times and places where no fleets are operating to over 0.5^-mo at times and places where fishing is intense. In contrast, natural mortality and fishing mortality in a unit stock model refer to the instantaneous rate of mortality of the population as a whole. Natural mortality estimated using the disaggregate model appears to be low in comparison to estimates from unit stock models. Higher rates of natural mortality estimated from unit stock models may be confounded with "mortality" due to emigration from the unit stock. In this sense, natural mortality is a measure of the departure from one of the major assumptions of unit stock models.

The analysis of fishery interaction concludes that purse seine fishing may have an adverse impact on pole and line fisheries for skipjack. This conclusion is based on the insupportable assumption that exploited fish populations approach equilibrium. Furthermore, all interannual and monthly variability has been averaged out, and the decrease is detectable only after a two year equilibration period. Whether changes on the order of 20% could be separated from interannual variability in real fisheries is unlikely. Verification of these hypothetical interactions using real data is essential before considering major changes in fishing policy.

Estimates of potential fishery interaction depend on the rate of natural mortality. In general, adverse effects of fisheries operating at a distance are less at higher levels of natural mortality. The relationship is not simple, however, because the estimated movement pattern also depends on natural mortality.

The advection-diffusion model appears to be a useful description of large-scale population movement. The parameters of the model can be estimated from existing data and the resulting empirical description of movement can be applied to predict the potential interactions between fisheries operating in different areas. Several improvements are desirable, however. The current parameterization in terms of regions and Seasons is too inflexible and ties behavior to specific times and places rather than to features of the environment. A more appropriate parameterization would express the movement parameters as functions of observable characteristics of the environment such as water temperature, depth, or oxygen concentration. Spatial gradients of these variables may also provide additional information from which to calculate model parameters.

The central role of natural mortality is obvious from this analysis. In almost all fisheries models, natural mortality is like a "residual", the quantity left over when mortality due to other causes (i.e., fishing) has been counted. Further work on the causes and age dependence of mortality in tunas (and fish in general) is warranted.

The numerical approximations used to solve equation (2) employ a technique known as "upwind differencing". This approximation stabilizes the advective terms by introducing "numerical viscosity" into the approximation. The effect of this artificial diffusion seems to be low estimates of the parameter (D) in certain season-region strata (see Table 2). Preliminary attempts using alternative approximations show promising results.

The results of this analysis of large scale skipjack movement suggest several topics for future research. These topics range in generality from research on the management of specific fisheries to matters of fundamental importance in population dynamics.

The analysis of interaction suggest that pole and line fisheries for skipjack may have declined by as much as 20% in response to the growth of purse seine fisheries. This conclusion is based on the assumption that exploited populations reach equilibrium after a twenty year period of exploitation at a constant level of fishing mortality. Alternative models of population growth, movement, and exploitation that incorporate the high level of interannual and interdecadal variability, both in oceanic production systems and in social and economic factors that drive fishing activities, are required. Preliminary efforts to create such models are currently being explored at the University of Hawaii and at the South Pacific Commission, but this challenging task will require more than a short term research and development effort.

The predicted reduction in pole and line yield must be considered to be a theoretical potential until empirical verification of the reduction is obtained. Attempts at verification through analysis of the historical record of catch and effort statistics have not yielded convincing evidence. Alternative methods of verification should be attempted. One alternative method is the use of what has become known as "adaptive management". Certain areas could be closed to purse seine fishing, but kept open to modest levels of pole and line fishing. Changes in pole and line catch rate inside the closed area could be compared with catch rates outside the closed area. An actual adaptive management experiment would require careful planning and a long-term commitment to catch rate monitoring in order to separate environmental effects from the effects of management policy.

Natural mortality emerges from the analysis of fisheries interaction as a critical process. Most fisheries models define natural mortality as mortality not attributable to fishing. Furthermore, natural mortality is usually assumed to be constant over the (exploited) life of the fish and throughout its range. In fact, we know very little about how fish, particularly tunas, die. It is generally recognized that predation is an important source of mortality. In some cases, emigration from the "stock" is recognized as a source of mortality. Senescence may also be important, but there is almost no information on the effects of ageng in fish. Research is required to improve our understanding of mortality in fish. A suite of model age-dependent mortality schedules, similar to those developed for human populations (Coale and Demeny, 1983), need to be developed for use in age-structured population models. A modest physiological research program on change in condition of tunas with age would greatly assist understanding senescence and could form a basis for developing age dependent mortality schedules.

Population models that explicitly include spatial variability and movement are clearly important in tuna fisheries. Several models of this type were discussed in the process of the Expert Consultation. In the future, these models will also be applied to other fisheries, for instance, to cod in the Atlantic ocean. However, more work is required to convince fisheries managers of the importance of spatially disaggregated models and to make such models more flexible and easily applied to different situations.

6. ACKNOWLEDGEMENTS

This project has been a complex undertaking involving several institutions. The work described was supported by FAO's project on "Cooperative Research on Interaction in Pacific Tuna Fisheries". The Joint Institue for Marine and Atmoshpheric Research at the University of Hawaii provided computing facilities and assistance in preparation of the manuscript under Cooperative Agreement Number NA37RJ0199 from the National Oceanic and Atmospheric Administration. The authors would like to thank Dr. Peter Bills for comments and suggestions regarding finite difference methods and Mr. Richard Bailey for implementing DOS graphics in X windows.

The views expressed herein are those of the authors and do not necessarily reflect the views of FAO or of NOAA or any of its sub-agencies.

Bard, Y. 1974. Nonlinear parameter estimation. Academic Press. N.Y. 341 p.

Brownlee, K.A. 1965. Statistical theory and methodology. John Wiley & Sons, N.Y. 590 p.

Boudreau, B. P. 1986. Mathematics of tracer mixing in sediments. I. Spatially-dependent, diffusive mixing. Am. J. Sci. 286:161-198.

Carnahan, B., H.A. Luther and J.O. Wilkes. 1969. Applied numerical methods. John Wiley & Sons, N.Y. 604 p.

Coale, A.J., and P. Demeny. 1983. Regional model life tables and stable populations. 2nd Ed. Academic Press. New York. 496 p.

Feller, W. 1968. An introduction to probability theory and its applications. John Wiley & Sons, N.Y. 509 p.

Griewank, A., and G.F. Corliss. 1991. Automatic differentiation of algorithms: Theory, practice and application. SIAM, Philadelphia, 353p.

Huyakron, P. S., and G. F. Pinder. 1983. Computational methods in subsurface flow. Academic Press. London. 473 p.

Okubo, A. 1980. Diffusion and ecological problems: mathematical models. Springer, N.Y. 254 p.

Press, W.H., B.P. Flannery, S.A. Teukolsky and W.T. Vetterling. 1988. Numerical recipes in C. Cambridge Univ. Press, N.Y. 735 p.

Sibert, J.R. 1984. A two-fishery tag attrition model for the analysis of mortality, recruitment, and fishery interaction. Tech. Rep. Tuna Billfish Assess. Programme, S. Pac. Comm., 13:27p.

Sibert, J.R., and D.A. Fournier. 1994. Evaluation of advection-diffusion equations for estimation of movement patterns from tag recapture data. In: Shomura, R.S., J. Majkowski, and S. Langi (eds.). Proceedings of the First FAO Expert Consultation on Interactions of Pacific Ocean Tuna Fisheries; 3-11 December 1991, Noumea, New Caldeonia; Vol. 1: Summary report and papers on interaction. FAO Fish.Tech.Pap. (336/1):108-121.

SPC (South Pacific Commission). 1994. Oceanic Fisheries Programme: Work programme review 1993-94 and work plan 1994-95. Paper presented at SPC Seventh Standing Committee on Tuna and Billfish, Koror, Palau, July 1994, S. Pac. Comm., 66 p.

A. APPENDIX - MODEL DESCRIPTION

A.1 Integrated Model of Movement and Population Dynamics

The model used in this analysis is an extended version of the model initially proposed by Sibert and Fournier (1994). Movement patterns were estimated by fitting a population dynamics model that explicitly incorporates movement to tag and recapture data. Let N_xytr symbolize the density of tagged skipjack at point (x, y) in the ocean at time t that comprised tag release cohort r. The aggregate density of tags from all cohorts released up to time t is given by

(1)

The aggregated tag density is assumed to satisfy the following partial differential equation

(2)

Equation (2) partitions the rate of change of tag density into "dispersive" movements, "directed" movements, and mortality. The first two terms on the right hand side of equation (2) characterize dispersive movements in terms of a "diffusion" parameter D. The next two terms in equation (2) characterize "directed" movements in terms of two "advection" parameters (u, v) that describe eastward and northward directed movements. At any time t the movement parameters (u, v, D) may vary in space. The final term in equation (2) describes the loss of tagged fish at point (x, y) due to mortality. In this analysis the movement parameters (u, v, D) are interpreted to characterize the behavior of tagged skipjack at a point in space and time rather than to be features of the moving water in which the fish live. The form of the model is such that water movement could be easily added to the behavior parameters.

Mortality, Z, is separated into two components by

(3)

where F_xytf is the mortality due to fishing by fishing fleet f operating at point (x, y) during time t, and M is mortality due to other causes or "natural" mortality. Natural mortality is assumed be constant at all times and places. Fishing mortality is assumed to be a simple function of observed fishing effort

(4)

where E_xytf is the observed fishing effort of fleet f operating at point (x, y) during time t, and Q_f is a fleet specific proportionality constant or catchability coefficient.

Boundary conditions and initial conditions are required to complete the specification of the model described by equation (2). Closed (reflective) boundaries are imposed:

x = 0, x = X

y = 0, y = Y

The obvious initial conditions for Ñ are

where N_xy01 is the number of tagged skipjack released at point (x_r, y_r) at time 0 in tag cohort 1.

The predicted number of skipjack tags returned during 1 month is given by

(5)

where b_f is the reporting rate, i.e. the proportion of tags captured by fleet f returned with useable recapture information, Dt = 1 month, and Ñ satisfies equation (2) for month t.

A.2 Solution of PDE

Equation (2) is solved on a grid with a spatial resolution of 60 Nmi with a discrete time step. The finite differencing scheme and the solution methods were selected on the basis of numerical stability, convergence to a solution of equation (2) and robustness to coarse step sizes in both time and space. The second derivatives in the dispersive terms are approximated by three point central finite differences. The first derivatives in the directed terms are approximated by two point finite differences computed in the direction opposite the sign of the parameters (u, v), a scheme known as "upwind differencing" (Press et al., 1988). The following finite difference approximations are used:

In the subsequent discussion, the continuous subscripts x, y, t are replaced by their discrete equivalents i, j, k.

The "alternating direction implicit" method was used to solve the equations. This method is very robust and converges to a solution unconditionally for all steps sizes (Carnahan et al., 1969; Press et al., 1988). This method was selected because equation (2) is solved hundreds of times under the control of a numerical function minimiser. A method that is tolerant of large step sizes and possibly unrealistic parameters is therefore required. It was found in practice that a time step of ¼ months was sufficiently accurate.

A.3 Movement Pattern

Solution of equation (2) requires parameters to be specified at each grid point and at every time step. Regions and seasons are introduced to reduce the number of parameters that must be specified. A region is defined as a subdivision of the model area over which the movement parameters (u, v, D) are constant, and a season is defined as the period of time during with these parameters are constant. Let G_ij be a matrix that contains the region number for each model cell indexed by (i,j) and S_k be a vector that contains the season for each time step indexed by k. The model parameters are specified at each grid point by the following equations:

(6)

, , and are matrices of parameters to be specified. For a model with 10 regions and 2 seasons, there are 60 movement parameters so , , and . are 2 x 10 matrices. These 60 parameters are used to specify values for u, v, and D in the 2100 points in the 60 x 35 model area used in the present analysis.

A.4 Parameter Estimation

Observed numbers of tag returns, C_xytf, are related to predicted numbers of tag returns, xytf, by a Poisson likelihood function.

(7)

This function assumes that the predicted number of tag returns in each cell during one month is the expected value of a random variable with a Poisson distribution. This distribution is appropriate for an observation of a rare event such as the radioactive decay of an atom (Feller, 1968) or the recapture and return of a tagged skipjack.

Maximum likelihood parameter estimates are obtained by finding the values of the parameters that maximize equation (7). The maximization is accomplished by minimizing the negative log of equation (7) using a quasi-Newton numerical function minimiser which in turn depends on the gradient of partial derivative computed using adjoint functions (Griewank and Corliss, 1991).

(8)

The standard deviations of the parameter estimates are computed from the approximate Hessian matrix, the elements of which are

where q_i and q_j are any two model parameters. The variance of the parameter estimates are the diagonal elements of the inverse of the Hessian matrix S = H^-1 (Bard, 1974).

A.5 Fishery Simulation

For the purposes of analyzing fisheries interaction, the skipjack population density, P, is assumed to satisfy a PDE similar to that used to describe the dynamics of tagged skipjack

(9)

where R is the rate of recruitment of skipjack density at a point. The initial condition is P = 0 everywhere and the boundaries are closed. Equation (9) is solved by the same methods used to solve equation (2). Under some recruitment assumptions, the total skipjack density in the model area tends to a steady state. The simplest recruitment assumption is R = 1^-mo everywhere. Under this assumption, a steady state is approached after approximately 120 months using parameters estimated from the RTTP data.

The predicted monthly skipjack catch at a point is given by

(10)

where is the catch by fleet f and P satisfies equation (9). The units of P and are arbitrary and determined by the units of R. Thus population density and yield are expressed in arbitrary skipjack recruitment units. In the interaction analysis, yield is normalized to remove the ambiguity imposed by the skipjack recruitment units.