Blackmouth croaker is predominantly harvested by coastal trawl and gillnet fisheries in the northern and southwestern regions of Taiwan, with the trawl fishery having a higher catch proportion in recent years (> 60% in general; Fisheries Agency, 2022). Two distinct stocks were identified in the north and southwest (Hwang and Chen, 1984), and this study concentrates on the southwestern stock due to its higher prevalence in the overall catch. The species is captured alongside nearly 200 demersal species throughout the year by the trawl fishery (refer to Section 2.2.1.2), with a moderate higher in catch during the spawning season from January to May (Hsiao et al., 2017). Mitochondrial cytochrome b gene analyses on 59 samples from five locations in the southwestern waters confirm that blackmouth croakers in this region belong to the same stock (unpublished project report by author SK Chang). The majority of blackmouth croaker caught by the trawl fishery in the southwestern waters (> 90%) is landed at Keziliao Fish Market in Zihquan, Kaohsiung City, the primary market for vessels operating in this region, particularly those engaged in coastal trawl fishing (Huang et al., 2022).

Daily vessel-specific landing data, including species-specific weights spanning from 2011 to 2021, were obtained from the market. The data were from all fisheries (no separation of fishing gears) and included a total of 268,176 landing events (annual average = 24,380) for all species and 33,907 landing events (annual average = 3,082) for blackmouth. Interviews with fishers indicated that, given the species’ high market value in auction markets for whole fish, at-sea discards or dumping are unlikely. This observation supports the reliability of using individual vessels’ landing data as a representative source for their catch information (Chang et al., 2019).

Most small-scale offshore fishing vessels in Taiwan have installed simplified VDR to evaluate fuel subsidies from the active moving hours of the vessels at sea (Chang, 2016). Akin to the “Black Box” on airplanes, the VDR is an equipment-fitted onboard ships that record the various data on a ship which can be used for reconstruction of the voyage details and vital information during an accident investigation (Bhattacharjee, 2019). It was required to be installed on passenger ships and ships other than passenger ships of 3,000 gross tonnage and upwards constructed on or after 1 July 2002, under International Maritime Organization (IMO) regulations that entered into force on 1 July 2002 (IMO, n.d).. The simplified device provides high-resolution temporal position and speed data like that of VMSs ( Figure 1 ), not in real time but with lower device cost, no data transmission fee, and higher resolution at 3-min intervals. This information incidentally can be used to estimate fishing effort through a general criteria of 2–5 knots for Taiwanese trawl fishery that was developed by a five-step procedure for reviewing fishing patterns and defining speed criteria (Chang, 2016). Constrained by the restriction of the VDR data provision policy of the Fisheries Agency (FA), VDR data of vessels with a 5-year average of blackmouth catch exceeding 1.5 mt were obtained. A total of 29 coastal trawlers (about 25% of the vessels that have landed blackmouth for the past 5 years) were used for this study (the sample vessels), ranging from 20 to 100 gross registered tonnage, whose blackmouth catch during 2011–2021 was 381,743 mt, 62% of the total blackmouth catch, from 20,756 trips. Annual number of sample vessels ranged from 22 to 28.

Based on the work of Chang (2016) and interviews with fishermen, the following rules were established and applied to the VDR data to estimate fishing effort through an automated algorithm: (1) the VDR records that leaving from and later returning to the fishing port were marked as a trip. The number of VDR records between leave and return must be over 40 (about 2h). The VDR records located continuously in port were excluded. (2) The records with a speed of 2–5 knots for at least ten consecutive records (approximately 30 min) were considered as fishing. The short period of slow speed leaving and returning the port was considered as navigating. (3) The location in a resolution of 0.1° square that the vessel stayed longest during the fishing day was defined as the representative fishing location. (4) The records with speed > 5 knots were assumed to be navigating (e.g., transiting to fishing grounds).

The spatial distributions of VDR records identified as fishing according to rule (2) are illustrated in Figure 2 . It is important to note that the VDR record distribution may cover a larger area than the CPUE distribution used for standardization analyses, as only one location in a trip was chosen to represent fishing according to rule (3).

Following the outlined procedures, a comprehensive dataset of trip records was collected from the sampled vessels. Multiple sets were deployed during each trip, as illustrated in Figure 1 . However, the fishing hours per trip exhibited variability based on fishing conditions and weather. Therefore, fishing effort was calculated as fishing hours. Subsequently, the trip-basis effort data (measured in fishing hours) were integrated with the corresponding market landing data of the trip for the sample vessels, to construct a trip-basis logbook-like data with time and location information for the calculation of nominal CPUE (kg per hour). On average, these sample vessels undertook approximately 1,887 trips with 3.5 mt of blackmouth annually, each spanning an average of 5.3 fishing hours per trip (excluding transition time). Since the multispecies feature of the trawl fishery and the relatively low abundance of blackmouth stock, about 60% of the trips contained zero blackmouth catch.

This study constructed a model structure with 15 model designs (M1–M15) to compare the influences of the three statistical models (GLM, GAM, and VAST) and the four key factors: environment factor (with or without), vessel factor (treated as fixed or random effect), spatial factor (with or without), and target factors (addressed by k-means or PCA methods) ( Table 1 ). A flowchart showing the compositions of the models investigated is shown in Figure 3 . Each model has a corresponding model to compare with; for example, M1 is to compare with M0 to see the effect of environmental factor, and M2 is to compare with M1 to see the effect of treatment of target factor using k-means versus PCA (the last column of Table 1 ). The study performed four sets of comparisons: the first set (S1) compared M1–M5 for environment factor and treatment of target factor. Including environment factor in the rest models, the second set (S2) compared M6–M7 for vessel factor and target factor, without consideration of spatial factor, and the third set (S3) compared M8–M11 for vessel factor and target factor under the consideration of spatial factor. With fixed vessel factor and PCA-treatment of target factor, the fourth set (S4) compared M12–M15 for statistical models. QQ plots of all the models are provided in Supplementary Material and suggested that they all conformed with the lognormal distribution assumption, affirming the appropriateness of the assumption regarding the error distribution for the CPUE standardization.

The study examined the predictive capability of two categories of methods introduced by Hinton and Maunder (2004) for CPUE model evaluation. The first category applied pseudo-R ² (higher value preferred), a log-likelihood–based measurement representing the improvement in model likelihood over a null model (Faraway, 2016; Chang et al., 2019).

The second category was the cross-validation (CV) method. The study applied a stratified random sampling approach using year and vessel as strata in the fivefold CV procedure. After repeated 20 times, the CV-R ² that determines the overall correlation between the actual and predicted value (Li et al., 2011; Chang et al., 2017; Zhang et al., 2023) was calculated and averaged. CV-R ² was used to evaluate model’s EP.

Although the Akaike information criterion (AIC), falling into the second category, is also a valuable and widely employed tool for model comparison, it is crucial to recognize that its values are impacted by the sample size, because they are derived from a likelihood function. Therefore, as elucidated in Chang et al. (2019) and the citation, caution is warranted when comparing one-stage GLMs and two-stage delta-GLMs with different sample sizes using AIC. Additionally, specifying the variance parameters of the likelihood function in two-stage delta-GLM or delta-GAM cases is a complex undertaking. Therefore, this study used CV-R ² when evaluating model EP.

Figure 1 provides an example of a VDR data track for effort estimation; an explanatory plot modified from Chang (2016) is provided in Supplementary Material . The 2011–2021 data obtained from the FA of Taiwan contains 29 vessels and 45,000 records in the study area. Using the procedures described in Section 2.1, about 43,530 trips with 210 tons of blackmouth catch and 230,809h of fishing effort were constructed from the data. Figure 2 shows the distribution of the fishing efforts, which were used for the following analyses.

Two R ² (pseudo-R ² and CV-R ²) of the 15 models (M1–M15) are shown in Tables 2 (for model sets S1–S3) and Table 3 (for S4). The environmental effect was evaluated through S1 comparison under the consideration of vessel factor as a fixed effect: the R ² between M2 (including environment factors) and M1 (without environment factors) are not much different. The same observation was noted for M3/M4 that used PCA for the target factor, compared to M1/M2 that used k-means. Another pair of tests (M5 and M12) using GAM instead of GLM showed that the model with environment factors (M12) has better EP (i.e., higher R ²) than without (M5), although the differences were not substantial. The EP of M12 were much better than M2 and M4, which might suggest that the nonlinear features of environment factors could be better explained by the nonlinear GAM model.

The effect of the method in addressing target factor was explored by comparing M1/M2 with M3/M4 and M8 with M9 when considering vessel factor as a fixed effect and by comparing M6 with M7 and M10 with M11 when considering vessel factor as a random effect ( Table 2 ). All these comparisons suggested that the PCA method performed better than the k-means method.

The effect of the treatment of vessel factor was evaluated by comparing of M2 with M6, M4 with M7, M8 with M10, and M9 with M11 ( Table 2 ). These comparisons suggested that models perform slightly better when the vessel factor is considered as a fixed effect instead of a random effect.

The spatial effect and EP of different statistical models were evaluated through comparisons with models of M4 and M12–14 that included environment factors, treated vessel factor as a fixed effect, and addressed target factor by PCA ( Table 3 ). The results suggested that considering spatial random effect (delta-GLMM or delta-GAMM, M13–M14) resulted smaller R ² than compared to those without (delta-GLM and delta-GAM, M4 and M12), and the M12 using GAM type model (delta-GAM) without consideration of spatial effect has the second-best EP (highest R ²: pseudo-R ² = 0.604 and CV-R ² = 0.594).

The above comparisons using traditional GLM/GAM type statistical models showed that spatial effect might not be important in the standardization. However, when applying VAST and including spatial and spatiotemporal effects, the EP measures of M9 and M15 were much higher than those of M13 and M14, which applied GLM/GAM (linear and nonlinear) approaches. It might imply that spatial effect was important but needed a better approach to handle. The nonlinear-featured VAST (M15) demonstrated to have the best EP among all model designs (pseudo-R ² = 0.613 and CV-R² = 0.601, Table 3 ).

Table 4 was created for explanatory purposes on the significance of the factors. Table 4 indicates that all the main effects were significantly different from zero at a significance level of α = 0.05 for M4 and M12, except for the chlorophyll factor that was not significant for the PCM sub-model of delta-GLM. Sea temperature and salinity at a depth of 109 m were excluded from PCM but were included and statistically significant in ZPM in the process of model selection. Figure 8 shows the estimates of fixed-effect factors of the year, quarter, and vessel from M12 (GAM), along with the smoothers for the PCA target factor and environment factor. The figure suggested that the highest CPUE occurred in the third and fourth quarters and that the variation among vessels is high. The effect of environmental factors seemed not very strong ( Table 4 ) and without clear patterns ( Figure 7 ). Most of the variation could be explained by the targeting effect, especially the PC1, PC2, and PC4 ( Table 4 ): blackmouth CPUE has general positive (but nonlinear) correlation with PC1 ( Figure 7 ) which was negatively correlated with “others” fishes ( Figure 6 ) (i.e., higher PC1 value has lower “others” fishes composition and higher Sciaenidae composition), and was nonlinearly negative correlation with PC2 and PC4 ( Figure 7 ), which was negatively correlated with Sciaenidae ( Figure 6 ).

The standardized relative CPUE indices are shown for all the 15 models in Figure 9A and a subset of specifically selected eight models in Figure 9B . The general trend shows a decrease to 2013, an increase in 2014 and 2015 that was coincident with a weak to very strong El Niño event (Golden Gate Weather Services, 2023), and a continuous decline since then (although some models, especially the models used VAST, showed a noticeable increase in 2018 and then declined again). In the last year of the analysis in 2021, the indices show an increase to various degrees for all models using GLM- and GAM-type approaches but a continuous decline for models using the VAST approach except for models M8 and M10 that used k-means for targeting effect.

For the specifically selected subset of the CPUE series ( Figure 9B ): (1) M8 and M9 employ different approaches to handle the target factor, with M8 utilizing k-means and M9 employing PCA, resulting in significant differences in the observed trends. M9 has higher EP than M8. (2) M9 and M11 differ in the way of treating vessel factors, and they are almost identical. (3) M5 and M12 differ in whether the environmental factor is considered, and the trends were basically similar with some years of variations. (4) M12 and M14 differ in whether a spatial factor is considered, and the results were very similar. (5) M4 and M12 differ in the statistical model used (GLM or GAM), and the indices were very different. Finally, (6) M9, M11, and M15 all used VAST but differ in linear or nonlinear treatment of nonlinear factors (environment and PCs). The index of M15 differs from the M9/M11 only in the beginning period and shows a slightly stronger decline in the last year. All these models exhibit different trends of indices from GLM/GAM type models. M12 and M15 utilized a nonlinear approach with the GAM and VAST methods to address nonlinear environmental factors, and exhibited higher EP compared to M4, M9, and M11, which employed methods featuring linear GLM and VAST.

Logbooks play a crucial role in providing essential catch and effort data for calculating the index of stock abundance. However, numerous coastal and small-scale fisheries are exempted from logbook submissions. In such cases, particularly for high-commercially valued species, which are more susceptible to intense exploitation, landing data, including date information and vessel identity, serves as a viable alternative to represent the “catch”. Additionally, auxiliary data such as VMS, CSRS, VDR, or the automatic identification system (AIS) can be employed to estimate fishing efforts through fishery-specific algorithms (Gerritsen and Lordan, 2011; Chang, 2014, Chang, 2016; Sonderblohm et al., 2014).

Although VMS data have been extensively utilized for scientific purposes, the associated installation and transmission costs are high, making them less feasible for small-scale coastal vessels with limited capital. CSRS has its application limitations, such as the special system requirements to identify fishing vessels from the port and the potential degradation of radar data quality by environmental effects (Chang, 2014).

AIS, introduced by the IMO to enhance maritime safety and prevent ship collisions, offers comparable data to VMS but with improved temporal resolution at a significantly lower cost and can be used to derive fishing efforts and inform fisheries management (Natale et al., 2015; James et al., 2018). This study conducted a preliminary investigation into the utilization of AIS from Taiwanese coastal fishing vessels for scientific purposes, revealing challenges such as the short data period (since the AIS was introduced to fishing vessels late), low data coverage (because installation and operation of AIS are not compulsive for small-scale coastal vessels), and difficulties in identifying fishing vessels (due to the absence of a cross-reference template for AIS databases and fisheries management systems). These challenges are likely prevalent in other coastal countries.

In contrast, the VDR system proves more manageable for the Taiwanese case. The IMO’s Maritime Safety Committee mandated passenger ships and large vessels to carry a standard VDR for accident investigations since 2000, with a 2006 amendment introducing a requirement for a simplified VDR (IMO, 2006). Recognizing the potential for fraud in fuel cost subsidies, the FA of Taiwan mandated the installation of simplified VDR in 2010 (Chang, 2016). This cost-reduced VDR, which records date, time, GPS position, speed, and direction of the vessel at 3-min intervals, not only addresses fraudulent activities but also provides data for scientific use at minimal additional cost. The affordability of installation and transmission makes it feasible for sample vessels to carry VDRs for independent research purposes. This study demonstrates the feasibility of using VDR data to construct nominal CPUE and derive an abundance index for declining fish stocks.

Confounding factors inherent in the nominal CPUE data need to be carefully addressed before the CPUE can represent the abundance. This study designed 15 models to examine the effect of targeting, vessel, and environmental and spatial factors and explore the best model specification for the blackmouth stock off southwestern Taiwan.

This study utilized two R ² as EP metrics. However, it is crucial to emphasize that a higher R ² do not necessarily equate to a more accurate index of relative abundance. The index of relative abundance, a pivotal input for stock assessment models, directly reflects changes in population abundance over time (Francis, 2011). The primary aim of CPUE standardization is to mitigate the influence of all factors except fish density on CPUE, ensuring that the standardized index of relative abundance approximates population abundance (MaunderTan and Punt, 2004). While a higher R ² suggests a better fit and more accurate predictions for nominal CPUE data, they do not confirm whether the model accurately identifies and addresses other factors impacting catchability. Therefore, a careful consideration of the history of the fishery and the stock (e.g., in terms of fleet dynamics, possible changes in fishing gear configurations, economic incentives, shifts in species distribution) is required when selecting the best models to use for CPUE standardization.

The abundance indices derived from the best model (M15) and the second-best model (M12) both exhibited a general declining trend for the southwestern blackmouth stock. This marks the first abundance index developed for demersal fish species in the region. Limited academic studies on the stock status of fish resources in the region exist in the literature due to a scarcity of catch and effort data (Liao et al., 2019). Although abundance indices for other demersal fish species are unavailable, their landings indicate a similar declining trend, suggesting a broad degradation of fisheries resources in the region (based on fisheries statistical yearbooks, Fisheries Agency, 2022; Shao, 2023). Indirect methods, such as trophic models (Liu et al., 2009) or fishermen surveys (Liao et al., 2019), along with numerous gray literature sources, caution that overfishing has occurred to Taiwanese coastal blackmouth resources due to the overcapacity of fishing vessels and illegal fishing (Chang et al., 2010; Liao et al., 2019). Comparisons of growth parameters over a decade also imply high fishing pressures on blackmouth, resulting in elevated growth coefficients and smaller asymptotic lengths (Huang et al., 2022).

In addition to overfishing, environmental pollution is asserted to contribute to the degradation of the fishery environment in the region. Repeated pollution events in rivers flowing into the area have been documented for years without complete remediation (Liu et al., 2015; Wang et al., 2021; Liang et al., 2024). However, further studies are required to establish the impact of contaminations on fisheries resources.