Revealing the global longline fleet with satellite radar

To estimate the total number of non-broadcasting vessels, including those that were not detected by SAR, we: (1) obtained SAR detections of vessels from RADARSAT-2 and the corresponding vessel lengths as estimated from the SAR image; (2) processed a global feed of AIS data to identify every broadcasting vessel that should have appeared in the SAR images at the moment the images were taken; (3) developed a novel technique to determine which vessels in AIS matched to detections in SAR, which AIS vessels were not detected by SAR, and which SAR detections represented non-broadcasting vessels; (4) after matching SAR to AIS, we could then (a) model the relationship between a vessel’s actual length and the length as estimated by the SAR image (Fig. 3b) and (b) model the relationship between the likelihood that a vessel is detected and its length (Fig. 3a); and (5) finally, we combined these relationships to develop an estimate of the number and lengths of non-broadcasting vessels in the region.

SAR imagery and vessel detections

Working with the satellite company Kongsberg Satellite Services (KSAT), we tasked the Canadian Space Agency’s satellite RADARSAT-2 to acquire SAR images from its ship detection mode (DVWF mode, GRD product), with a pixel size of about 40 m and a swath width over 400 km (19). These images were processed following standard procedures for GRD products (e.g. applying radiometric calibration and geometric corrections)^29,30. Vessel locations were extracted from the images with the widely used ship detection algorithms, which discriminates objects at sea based on the backscatter difference (pixel values) between the sea clutter and the targets³¹. Vessel lengths were estimated by measuring distances directly on the images with the aid of a graphical user interface tool³¹.

Identifying Vessels using AIS

In each region, AIS data, obtained from satellite providers ORBCOMM and Spire, were processed using Global Fishing Watch’s data pipeline¹. The identities and lengths of all AIS devices that operated near the SAR scenes in both space and time were first obtained using Global Fishing Watch’s database¹. To be sure vessels were identified correctly, two analysts reviewed the tracks of every AIS device in each region.

In both regions, it is common practice for fishers to put AIS beacons on their longlines, likely to aid in retrieving them, meaning that many AIS devices were longline gear and not vessels. Because gear outnumbered vessels by several-fold, it was critical to differentiate gear and fishing vessels. In the Indian Ocean, 521 unique AIS devices associated with gear were detected that were likely within the SAR scenes, and 390 unique AIS devices associated with gear in the Pacific that were likely within the SAR scenes. Transponders were determined to be associated with gear by inspecting the name broadcast in the AIS messages (gear frequently broadcasts one of several standard names and/or a voltage reading) and classification using the Global Fishing Watch vessel classification algorithm¹. Most gear also had an MMSI number (unique identifier number for AIS) that started with 1, 8, or 9 or broadcast names that signified gear. We eliminated all gear from the analysis because (1) these gear buoys have reflectors that are only ~ 1 m in size, and they should not be visible in ~ 40 m resolution SAR images, and (2) we found that gear matched to SAR detections only when traveling faster than 2 knots (and thus was on the deck of a boat); of 159 instances of gear in scenes where the gear was traveling slower than two knots, zero matched to a radar detection (Fig. S9).

Generating probability rasters for matching AIS to SAR

Most AIS positions did not correspond to the exact time when the SAR images were taken. Hence, to determine the likelihood that a vessel broadcasting AIS corresponded to a specific SAR detection, we first developed probability rasters of where a vessel was likely to be minutes before or after a GPS position was recorded (Figs. S1,S2). We mined one year of global AIS data, including roughly 10 billion GPS positions, and computed these rasters for six different vessel classes (trawlers, purse seines, tug, cargo or tanker, drifting longlines, and others) and considered six different speeds (1, 3, 5, 7, 9, and 12.5 knots) and 36 time intervals (− 448, − 320, − 224, − 160, − 112, − 80, − 56, − 40, − 28, − 20, − 14, − 10, − 7, − 5, − 3.5, − 2.5, − 1.5, − 0.5, 0.5, 1.5, 2.5, 3.5, 5, 7, 10, 14, 20, 28, 40, 56, 80, 112, 160, 224, 320, and 448 min).

For example, we queried a year of AIS data to find every example of where a tugboat had two positions that were 10 min apart from one another when the vessel had been traveling at 10 knots at the first position. We then recorded each of these locations relative to the location the vessel would have been if it traveled in a straight line, with x coordinates being in the direction of travel and the y coordinates being perpendicular to the direction of travel. When collected for hundreds of thousands of examples across the AIS dataset, the result is a heatmap of where tug boats are located 10 min after a position when it was traveling at 10 knots. The raster is centered on a point that is the extrapolated position of the vessel based on its speed. For instance, the purse seine raster that corresponds to a vessel traveling between 6 and 8 knots between 96 and 128 min after the most recent position is centered at a point that is 13.1 km (7 knots × 112 min) straight ahead of the direction the vessel was traveling. Figure S1 shows samples of these rasters for different vessels.

We built rasters of 1000 by 1000 pixels for each vessel class and time interval, with the area covered by the raster dependent on the time interval (longer time intervals imply longer traveled distances, covering more area). The scale of each pixel was given by:

where Δm is the time interval in minutes, and pixel width is measured in km. Thus, if the Δm is under one minute, the entire raster is one kilometer wide with each pixel one meter by one meter. If the time is 10 min, then each pixel is 10 m wide, and the entire raster is 10 km by 10 km.

Since the pixel width varies between rasters, the units of the rasters are probability per km², thus summing the area of each pixel times its value equals one. Six vessel classes with 36 time intervals for each and six speeds led to 1296 different rasters. This probability raster approach could be seen as a utilization distribution³²—for each vessel class, speed and time interval—where the space is relative to the position of the individual.

Combining probability rasters to produce a matching score

For a few vessels (~ 4%) there was only one AIS position available before or after the scene. This resulted from a long gap in the AIS data due to poor reception, a weak AIS device, or cases where the vessels disabled their AIS. For these vessels, we used the raster values for a single position. For the vast majority of vessels, however, there was a GPS position right before and after the scene, and thus two probability rasters. We used two methods to combine these probability rasters to obtain information about the most likely location:

Multiply and renormalize the rasters

To multiply the rasters, we interpolated the raster values, using bilinear interpolation, to a constant grid at the highest resolution between the before and after rasters. Then, we multiplied the values at each point and renormalized the resulting raster (Fig. S2):

where p_i is the probability in vessel density per km² at location i, p_ai is the value of the raster before the image, p_bi is the value of the raster after the image. The denominator is the sum of all multiplied values across the raster, scaled by the area of each cell, da.

Weight and average the rasters

For this method, we weighted the raster by the squared value of the probabilities of that scene. This has the effect of giving the concentrated raster a higher weight, thus weighting higher the raster that is closer in time to the image:

and the weighted average at location i is:

where w_a is the weight for raster a, w_b the weight for raster b (calculation analogous to w_a’s in Eq. 3), p_i is the probability in vessel density per km² at location i.

To determine whether we should multiply (Eq. 2) or average (Eq. 4) the probabilities, we compared the performance of these two metrics against a direct inspection of the detections. We found that at short intervals, multiplying the rasters and renormalizing often made probability values extremely small (< 1 × 10^–4) despite direct inspection by an expert analyst suggesting there was a high likelihood of match (i.e. that the probability of a vessel being in location i should be high). Overly small values of multiplied probabilities could happen, for example, when there is a position one minute before the image and another one minute after, in which case both rasters contain almost the same information and, therefore, are virtually identical. This is analogous to squaring the probability values of one of the rasters (Fig. S3). On the other hand, at longer time ranges the averaging scheme tended to overestimate the matching rates of some highly unlikely matches. We established an ad hoc rule based on our data: when the closest position was within 10 min of the image, we used a weighted average (Eq. 4), and multiplied (Eq. 3) at longer time ranges.

Ranking and matching potential SAR to AIS pairs

We computed a matrix of scores of potential matches between SAR and AIS. We then greedily assigned matches, i.e. the AIS-SAR pair with the highest score (highest p_i from Eq. 4) is selected as a match and the corresponding row and column are removed from the matrix. Then the next highest score is assigned and so on until all pairs have been assigned. In cases where a single AIS vessel presented an equally high matching score to multiple SAR detections, or vice versa, we manually reviewed the matching pairs (mostly for cases where the scores were within a factor of 100 of one another).

Choosing a threshold for accepting SAR to AIS matches

A major challenge is deciding when the probability score is too low to accept a match. We reviewed matches with two independent analysts, classifying the AIS to detection pairs as: “likely matches” when there was a high likelihood of matching, “potential matches” when it was not possible to determine, and “unlikely matches” when there was a low probability of matching. Comparing these classifications with different score thresholds, we found that a matching score between 5 × 10^–6 and 1 × 10^–4 vessel km⁻² provided a reasonable threshold for accepting/rejecting a match. Only a few pairs had scores between this range; most potential matches had higher or lower scores, and thus our results are not very sensitive to the different thresholds. Only nine likely matches (2% of total) between SAR detections and vessels broadcasting AIS were ambiguous, meaning that a vessel broadcasting AIS could have matched to multiple SAR detections or, conversely, a SAR detection could have matched to multiple vessels broadcasting AIS.

This analysts’ judgment is roughly in line with a theoretical estimate. Given a single SAR detection and AIS vessel, there are three possible options: (1) the detection represents the AIS vessel, (2) the AIS vessel was not detected by SAR and the detection represents a non-broadcasting vessel, or (3) the AIS vessel was not detected by SAR and the detection represents a false positive. We should match the SAR detection to the AIS vessel if the probability of (1) is greater than the probability of (2) and (3):

where ({p}_{v}) is the probability density of the vessel presence at the location of the SAR detection (the score listed above), ({p}_{d}) is the probability that the vessel is detected by SAR, ({d}_{d}) is the density of non-broadcasting vessels in the region, and ({p}_{f}) is the density of false detections in the scene. The greater ({p}_{d}), the more dark vessels there are in a scene, and the more likely it is that any given detection is a dark vessel instead of a vessel broadcasting AIS. The right-hand side of the equation ({d}_{d}cdot {p}_{d} + {p}_{f}) should roughly equal the number of detections per unit area that do not match to AIS in the region. In other words, the probability of the vessel with AIS being at that specific location and detected by SAR (left side of the equation) should be greater than the probability of a dark vessel or a false detection at that location (right side of the equation).

The total number of unmatched vessels in each studied region normalized by total area covered gives a density of non-broadcasting vessels of 2.6–2.8 × 10^–5 vessels km^-2 (Indian Ocean) and 6.8–7.2 × 10^–6 vessels km⁻² (Pacific Ocean), similar to the thresholds estimated by analysts. For the most likely number of matched vessels, we use a threshold that is halfway between the higher and lower bound of the analyst (5 × 10^–5 to 1 × 10^–4), 2.5 × 10^–5 which is also roughly equal to the theoretical estimate of the Indian Ocean.

This threshold approach performed significantly better than a metric based on the distance between the SAR detection and the most likely location of the vessel, where the likely location is based on extrapolating speed and course of the position closest in time to the image (Fig. S4).

Determining whether a vessel with AIS was within a scene

Vessel positions from AIS are usually available before and/or after the SAR images, and sometimes it is unclear if a vessel should have been within the scene footprint at the time of the image.

To estimate the probability that a vessel (with AIS) was within a scene, we used the multiplied probability raster, summing the values inside the scene boundaries. This provides an estimate of the likelihood that the vessel was within the scene footprint at the time of the image. We applied this to every vessel that had at least one AIS position within 12 h and 200 nautical miles of the scene footprint. The vast majority of vessels were either very likely inside or outside the scene footprints, with 516 vessels having a probability of > 95% and only 16 having a probability between 5 and 95%. We filtered out all vessels that were definitely outside of the image footprint before matching.

Estimating the likelihood of detecting a vessel with SAR

The AIS data show that not all vessels broadcasting AIS were captured by the RADARSAT-2 images (Fig. 3a). Using the known lengths of detected vessels with AIS, we estimated the likelihood of detecting a vessel with SAR as a function of vessel length (Fig. 3a). For vessels shorter than 60 m, we approximated the detection rate as a linear function. Treating each vessel as an individual detection, we fitted the 50th percentile using quantile regression to approximate the detection rate. For vessels above 60 m, we assumed a constant detection rate as very few vessels above this length did now show up in the SAR images. Of the 46 unique vessels larger than 62 m, 42 were detected, implying a detection rate of ~ 91%. Given that it is highly likely that large vessels will be captured by medium-resolution SAR imagery, we manually reviewed these cases to confirm that they were (almost surely) inside the scene footprints at the time the images were taken.

We should note that the probability of detecting a vessel in SAR also depends on the sea state, incidence angle, polarization, material of the vessel, and orientation of the vessel. We are unable, however, to measure these effects directly so we cannot explicitly model these effects.With sufficient scenes, these effects should be randomly distributed across our scenes, so they likely account for some of the variability in detectability and the inaccuracy in our length estimates from SAR.

Estimating the number and length of non-broadcasting vessels

Because SAR does not detect all vessels, and because the length as estimated by SAR can be incorrect, there are many possible distributions of actual non-broadcasting vessels that could have produced the distribution of unmatched SAR detections that we found in the scenes. To estimate the most likely such distribution, we built a model to combine the two key relationships—between vessel length and likelihood of detection, and between vessel length and the length as estimated by SAR. This model allowed us to estimate, based on the number and distribution of SAR vessels, the likely number and distribution of actual vessels present (Fig. 3c,d).

We binned the likelihood of vessel detection as a function of length into 1 m intervals, yielding a vector (alpha) of length 400. We also binned into 1 m intervals the population of lengths of all detected vessels ((ell_{D})) as reported by AIS (i.e. number of vessels at each length bin), the population of expected SAR lengths ((ell_{E})), and the population of lengths of all vessels ((ell_{A}), the quantity we wish to estimate). Thus, (ell_{D}) can be expressed as the product of (alpha) and (ell_{A}):

where (odot) is the element-wise product. We then estimated a matrix (L_{{}}) that relates (ell_{D}) to (ell_{E}).

where each element (L_{ij}) represents the probability that a vessel with length in bin j would be estimated by SAR to be of length in bin i. We calculated these probabilities as lognormal probability density functions, with one distribution per column. To estimate the scale and shape parameters of these distributions, we first fitted a quantile regression using the (non-binned) lengths from AIS of detected vessels as the predictor for the lengths reported by SAR. Assuming that the predicted ^1/3 and ^2/3 quantiles (as shown in Fig. 3a) represent the quantiles of a lognormal distribution, allow us to calculate the shape and scale parameters. We chose a lognormal distribution because: 1) the variable of interest, length, was always greater than zero, 2) the population of lengths was skewed towards larger values, and 3) there is an explicit and relatively simple relationship between the lognormal quantiles and the shape and scale parameters that simplified the calculations.

Combining Eqs. (6) and (7) provides a relation between (ell_{A}) and (ell_{E}):

To estimate ({mathcal{l}}_{A}) we minimized an objective function (O({mathcal{l}}_{E},{mathcal{l}}_{o})) between the vector of expected counts binned by length (({mathcal{l}}_{E})) and the vector of counts observed in SAR binned by length (({mathcal{l}}_{o})). For this objective function, we chose the sum of the Kolmogorov –Smirnov distance between length distributions and the squared difference of the total numbers of detections. The first term controls the shape of the resulting distribution while the second one controls the magnitude. Specifically:

where:

Assessing the uncertainty in the estimation

To test how accurately our approach predicts the correct number of vessels, we performed a bootstrap simulation. We computed the vector (alpha) and the matrix L from a random subset of vessels with AIS that had a high confidence (> 95%) of appearing within the scenes. We then used our method on the SAR detections that matched the remaining vessels to predict the number of vessels they corresponded to ((ell_{text{A}})). By running 10,000 experiments we found a mean absolute percent error of + − 9% (Figs. S5 and S6). This provides a rough estimate of the uncertainty in our prediction due to the estimation process itself. We used the distribution of these samples to estimate the 90% confidence interval that we report with our estimates. We note that this uncertainty refers to the parametrization of the model and there may be other sources of error, such as the possibility that vessels without AIS have different radar properties (e.g. made out of materials with different reflectiveness), that we did not account for in our model.

Catch and effort data in the overlapping area between WCPFC and IATTC

We downloaded gridded effort and catch data from the WCPFC and IATTC websites, and compared the reported number of hooks and catch from September to December of 2019 for the area between − 140 to − 150 longitude and − 5 to − 15 latitude, a bounding box that contains our study region in the Pacific and which is entirely within both the WCPFC and IATTC convention zones. We found that the reported number of hooks for Korea is three times higher for the IATTC as it is for the WCPFC (Fig. S7), and the numbers of hooks also disagree by more than 10% for most other flag states. Catch is also 2.5 times higher for IATTC than for WCPFC for Korea as well, with catch also differing by more than 10% for most other flag states. This finding suggests that the different RFMOs may not be accounting for the same vessels in the overlap region between the two RFMOs.

Source: Ecology - nature.com