Spatial distribution of northern shrimp Pandalus eous in the Sea of Japan in relation to the seafloor environment

Abstract

To explore relationships between spatial distributions of northern shrimp Pandalus eous and environmental factors along the coastline of Honshu Island in the Sea of Japan, we built delta-type two-step generalized additive models (delta-GAM) based on bottom trawl surveys conducted in 2013–2022. The models provide the first quantitative analysis of the species’ habitat, showing that its distribution is associated with bottom sediment type, depth, slope and topographic position index (TPI) alongside the effects of year and region. From the delta-GAM response plots, species habitat preferences in the Sea of Japan were estimated as follows: seafloor deeper than 283 m; muddy rather than sandy bottoms; gently sloped to flat bottoms (<0.8° in slope); and valleys rather than ridges (TPI < 0.9). These results were reviewed in detail along with previously reported distribution records of northern shrimp. Standardized density (s-density) per fishing grid cell (10′ square latitude–longitude mesh) estimated from the delta-GAM results indicated that this species is widely distributed on the continental slope along the coast of Honshu Island. To test plausibility of the s-density analysis, we compared s-density per fishing grid cell with nominal CPUE (kg per haul) per fishing grid cell based on official logbook data from large offshore trawl fisheries. The two were generally positively correlated by year, and the delta-GAM results were assumed to be fairly robust. Finally, potential fishing grounds were explored based on the delta-GAM results.

Introduction

The spatial distribution of marine organisms is affected by various environmental factors, including abiotic (e.g. water temperature, depth) and biotic factors (e.g. intra- and interspecific interactions, including predation) (e.g. Lalli and Parsons, 1997). Species distribution models (SDMs) incorporate data on the abiotic and biotic environment to predict species’ geographic and environmental ranges (Robinson et al., 2017), which in turn can be used to predict the potential fishing grounds of fishery species and improve the precision of stock assessments (e.g. Zainuddin et al., 2008; Cao et al., 2017). SDMs for marine species have been built with various data sources such as field surveys conducted by researchers and global and regional databases of biomass and environmental conditions (Melo-Merino et al., 2020). Research vessel surveys that can simultaneously sample marine organisms (i.e. collect catch data) and conduct oceanographic observations (i.e. collect environmental data) are particularly effective for building SDMs.

The northern shrimp, Pandalus eous Makarov, 1935 is a cold-water shrimp distributed on the continental shelves and slopes of the North Pacific at depths of 16‒1380 m from the Sea of Japan to Puget Sound, as well as in the Chukchi Sea (Komai, 1999). Members of the genus Pandalus are protandrous hermaphrodites, meaning that individuals function as males early in life then change sex and reproduce as females (Bergström, 2000). Northern shrimp populations in the Sea of Japan are notable for their ontogenetic migration: post-larval shrimp settle at a depth of around 300‒350 m after the planktonic larval stage, then move deeper (400‒600 m) with growth in the juvenile stage; ovigerous females enter shallower waters (around 200 m) to hatch and then return to deep water (around 400 m depth) after hatching in winter to spring (Sadakata, 2000b). Until the early 1990s, northern shrimp was treated as a local population of Pandalus borealis Krøyer, 1838 (e.g. Kubo, 1965; Butler, 1980), but Squires (1992) distinguished northern shrimp from P. borealis by several morphological characteristics. We followed Squires (1992) and assumed that previous reports of P. borealis in the north Pacific were of northern shrimp.

Northern shrimp is a high-value fishery resource and is caught across much of the Sea of Japan. Off the Sea of Japan coast of Honshu Island (the largest island of the Japanese archipelago), it is primarily caught at depths of 200–500 m around Sado Island, off the northern Oki Islands, in Wakasa Bay, around Noto Peninsula, and on the Yamato Bank (Figure 1), and the majority of the catch in this region is accounted for by large offshore trawlers (Danish seiners with gross tonnage >15) operating west of 136.5°E (Figure 2).

Figure 1. Map of annual trawl survey locations (red dots). Dotted and dashed lines show 200 and 500 m isobaths, respectively. Uppercase letters (A–I) show administrative offshore trawl fishing areas, with borders shown by black lines.

Figure 2. Timeseries of northern shrimp catch (bars) along the coast of Honshu Island in the Sea of Japan and standardized CPUE (solid line) calculated from official logbook data for offshore trawl fisheries (after Naito et al., 2024). Total catch is divided into catch by large offshore trawlers operating west of 136.5°E (regions A‒F; red), large offshore trawlers operating east of 136.5°E (regions G‒I; blue) and other fisheries (green).

The relationship between the distribution of northern shrimp and environmental factors such as water temperature, sediment type, salinity and depth has long been discussed (Shumway et al., 1985). The species has generally been recorded in cold waters (usually 0‒5°C) on seafloors with gently sloped to flat bathymetry, and mud to sand/mud bottoms (e.g. Ivanov, 1970; Kitano and Yorita, 1978). However, such information on the distribution of the species, bottom sediment and seafloor topography has been mostly qualitative. Importantly, the factors for which there is a lack of quantitative knowledge on the species’ preferences could be critical to understanding its distribution (Kurata, 1981; Bergström, 2000), and so comprehensive analyses of these preferences along with the other environmental variables are urgently needed. Ship-based survey data could fill this knowledge gap and be used to build quantitative SDMs, which could identify the environmental factors that affect northern shrimp distributions in the Sea of Japan.

In the present study, we clarify the distribution of northern shrimp along the coastline of Honshu Island in the Sea of Japan. To fully understand the habitat of the species, we built a spatially explicit SDM by using a generalized additive model (GAM) framework based on bottom trawl surveys conducted at more than 200 stations over a 10-year period. The length of the monitoring period allows us to identify correlations between species distributions and environmental factors while reducing uncertainty stemming from interannual variability. Finally, we discuss the characteristics of the spatial distribution of northern shrimp along the Sea of Japan coast of Honshu Island based on the output of the SDMs and potential fishing grounds. Our results will hopefully contribute to improved stock assessment and more precise management of this species in this region.

Materials and methods

Study area

The study was conducted on continental slopes and shelves in the Sea of Japan at longitudes of 130°E to 140°E and at depths of 60–500 m along the coastline of Honshu Island in the Japanese Archipelago (Figure 1). This study area includes the main local habitat of northern shrimp, which extends from depths of 200 to 550 m (Ito, 1976). To compare the various fishing grounds used by large offshore trawl fisheries in this area, we divided the study area into nine regions (Figure 1A‒I) according to the conventional fishing areas of Japanese offshore trawlers. The regions were established by the Japanese Ministry of Agriculture, Forestry and Fisheries, and take into account the characteristics of benthic fauna, the stock structures of some benthic organisms, and differences in trawl operations and fishing strategies among regions.

Bottom trawl surveys

Bottom trawl surveys were conducted to sample the density distribution of northern shrimp and obtain environmental parameters in Japanese coastal areas of the Sea of Japan. In the western and central Sea of Japan (from 130°E to 138°E), an otter trawl (net mouth dimensions approximately 17 × 3.5 m, cod end with 20 mm stretch mesh) was deployed from T/S Tanshu-Maru (Hyogo Prefectural Kasumi High School, Japan). In the northern Sea of Japan (from 138°E to 140°E), otter trawls were unsuitable because of the prevalence of steep and rough seafloor bathymetry. Instead, a beam trawl (net mouth dimensions 6.8 × 1.0 m, cod end with 20 mm stretch mesh) was towed from R/V Tenyo-Maru (Japan Fisheries Research and Education Agency, Japan). A summary of the two bottom trawl surveys is shown in Table 1.

Table 1. Summary of bottom trawl surveys

At each site, the collected shrimps were counted and weighed onboard. Latitude, longitude and depth at the beginning and ending of each tow were recorded with a GPS (GP-500, Furuno, Hyogo, Japan) and fish finder (FCV-1500L, Furuno) aboard T/S Tanshu-Maru and with the SCANMAR system (SCANMATE6, Scanmar AS, Vestfold, Norway) aboard R/V Tenyo-Maru. Towed distance was calculated from the GPS data. The net-mouth opening of the otter trawl was monitored with a combination of the SCANMAR system and a Net Recorder (CN-2220-40L-P, Furuno). Bottom temperature was measured at each site with a compact salinity–temperature–depth profiler (RINKO-Profiler ASTD-150, JFE Advantech, Hyogo, Japan) and the values from the Net Recorder were used when the profiler was not available (0.4% of all records). The area swept was calculated as a product of the towed distance and the average net-mouth width (for the otter trawl) or the width of the beam trawl (6.8 m), and the density of northern shrimp (kg km^–2) was calculated as catch per area swept.

Details of SDM dataset

The dataset was composed of the survey year, region (A‒I, Figure 1), presence/absence and density of northern shrimp, environmental variables (depth [obtained from the field surveys], bottom-water temperature [from the field surveys], bottom sediment type [mentioned below] and bathymetry [mentioned below]). Bottom sediment type was based on the 15-arc-minute gridded seafloor sediment data for the NW Pacific provided by the Marine Information Research Center, Japan Hydrographic Association (Tokyo, Japan) (hereinafter the ‘SS data’). The SS data are compositional data consisting of the per cent coverages of 12 sediment types [land, rock, gravel, coarse sand, sand (medium sand), fine sand, sand and mud (weight ratio of sand 50%), sand and mud (weight ratio of sand 25%), mud, clay, coral and shells] linked to a 15′-square mesh grid. We recalculated the bottom sediment percentages of each mesh grid after excluding the land percentage (e.g. the sediment composition of a grid that is 50% land, 40% sand and 10% mud is recalculated to 80% sand and 20% mud) and summarized the SS data into a single principal component matrix to facilitate interpretation using principal components analysis. Negative values of the first principal component (PC1) were associated with a high per cent coverage of mud, and positive values with a high per cent coverage of sand (Figure S1). Because PC1 explained 60.8% of the total variance (Figure S1), we assumed that it sufficiently reflects the characteristics of the bottom sediment and included it as a parameter in the subsequent analyses. The slope and topographic position index (TPI; a terrain classification index for which a negative value indicates a valley and a positive value indicates a ridge) of each 0.25′-square mesh were calculated based on the depth of each grid and its eight neighbouring grids from the General Bathymetric Chart of the Oceans (GEBCO, British Oceanographic Data Centre) by using the ‘terrain’ function in the package ‘raster’ (Hijmans, 2022) in R (R Core Team, 2021). For each site in the trawl surveys, the PC1 of the nearest mesh grid was calculated as a proxy for bottom sediment, slope and TPI based on the position at the beginning of a trawl haul by using the function ‘st_within’ in the package ‘sf’ (Pebesma, 2018) in R.

Generalized additive models

To explore relationships between the spatial distribution of northern shrimp and environmental factors along the coastline of Honshu Island in the Sea of Japan, we applied a delta-type two-step model (Lo et al., 1992) incorporated with GAMs (Hastie and Tibshirani, 1986). GAMs can account for non-linearity in the spatial distribution of northern shrimp, which is abundant at a narrow range of depths and bottom temperatures (e.g. Kojima et al., 1969; Ito, 1976). In addition, a delta-type two-step model can separately identify environmental factors that influence presence/absence and density, and enables estimation of a standardized density considering zero catches. Prior to the model construction, we examined correlations among the potential explanatory variables other than year and region, to avoid issues stemming from multicollinearity. Because depth and bottom temperature were highly correlated with each other (r = −0.61, Figure S2), we excluded bottom temperature from both models, following the collinearity threshold identified by Booth et al. (1994) and Rhodes et al. (2009). The full delta-type two-step GAMs (delta-GAMs) were as follows:

$$\eqalign{& {\rm logit\ }( P ) {\rm \sim }\alpha _ 1{\rm} + {\rm Year} + {\rm Region} + s( {{\rm PC1}} ) {\rm} \cr & + s( {{\rm Depth}} ) {\rm} + s( {{\rm Slope}} ) {\rm} + s( {{\rm TPI}} ) {\rm} + {\rm Region\ast }s( {{\rm Depth}} ) \cr & ( {{\rm First}\hbox{-}{\rm step\ model}} ) } $$

$$\eqalign{& {\rm log\ }( {{\rm Density}} ) {\rm \sim }\alpha _ 2{\rm} + {\rm Year} + {\rm Region} + s( {{\rm PC1}} ) {\rm} \cr & + s( {{\rm Depth}} ) {\rm} + s( {{\rm Slope}} ) {\rm} + s( {{\rm TPI}} ) {\rm} + {\rm Region\ast }s( {{\rm Depth}} ) \cr & ( {{\rm Second}\hbox{-}{\rm step\ model}} ) } $$

where P and density are the probability of presence and density of occurrence of northern shrimp, and the first and second steps are represented by GAMs with binomial and Gaussian errors, respectively. α ₁ and α ₂ represent the intercept in each step model. Year and region were used as categorical explanatory variables with fixed effects, and the environmental factors (PC1, depth, slope and TPI) were used as continuous explanatory variables. s indicates a smoothing spline and is a penalized cubic regression spline which has its penalty shrink towards zero at high smoothing parameter values. The explanatory variables in each step are year (2013‒2022, 10 categories), region (A‒I; 9 categories), PC1 as a proxy for grain size (mud or sand, dimensionless), depth (m), slope (°), TPI (dimensionless) and the interaction of region and depth. Here, year represents annual abundance variability. Region explains putative regional differences caused by variables that could not be included in the model (e.g. settlement variability, regional fishing mortality and ontogenetic migration) but can also include differences in catch efficiencies between otter trawls (used in regions A–G) and beam trawls (used in regions H‒I). Depth reflects the general depth preferences of northern shrimp regardless of region as well as variations in individual weight. This is because northern shrimp undergo ontogenetic migration towards deeper habitats with growth (Sadakata, 2000b) (i.e. northern shrimp show a size-specific distribution in which small shrimp occur at shallower depths than large shrimp). The interaction term was added because the habitat depths and preferred bottom temperatures of northern shrimp differ among regions (regional differences in depth preferences: e.g. Ito, 1976; Sadakata, 2000b) and cannot be explained by individual region and depth effects alone. Moreover, as large (i.e. old) individuals are much heavier than small (young) ones and the response variable is total catch (in kilograms), the interaction of depth and region might not only reflect differences in the depth at which shrimp are abundant in each region but also differences in size distribution in each region.

The optimum degree of smoothing was determined by the generalized cross-validation criterion by using the function ‘gam’ in the package ‘mgcv’ (Wood, 2011) in R. We selected the best models with the lowest Akaike information criterion (AIC) value from the candidate models with all possible combinations of the explanatory variables, and used the selected model in the following analysis. Performance of the best first-step model was evaluated by using receiver operating characteristic (ROC) analysis (Swets, 1988). ROC evaluates the relationship between the sensitivity and specificity of a model and the results are summarized through the metric of the area under the curve (AUC). AUC varies between 0.5 and 1.0 and the model performance is higher if AUC is closer to 1.0. AUC was estimated in the ROC analysis by using the function ‘performance’ in the package ‘ROCR’ (Sing et al., 2005) in R. To confirm that the logarithm of density based on the best second-step model follows a Gaussian distribution, we drew a quantile–quantile (Q–Q) plot (Wilk and Gnanadesikan, 1968) by using the function ‘qq.gam’ in the package ‘mgcv’ (Wood, 2011) in R.

Standardized density (hereinafter s-density) in each fishing grid cell estimated with the best occurrence and density models was used to infer the locations of potential fishing grounds with high northern shrimp density along the Sea of Japan coast of Honshu. Briefly, s-density along a 10′ latitude–longitude square mesh grid (hereafter, the fishing grid) was estimated on the basis of the GEBCO bathymetry for depths of 200–550 m, the main habitat depth of northern shrimp in the Sea of Japan (Ito, 1976). We did not include Yamato Bank, Oki Bank or Hakusanse Bank in the prediction domain because these potential fishing grounds are geologically separated from the continental shelves along the Honshu coast and are not covered by the surveys. The fishing grid, bottom sediment (PC1), slope and TPI data for predicting the distribution of northern shrimp were extracted or calculated from the latitude, longitude and depth of each GEBCO grid, and were applied to all possible combinations of years and regions. The s-density was then calculated as a product of catch probability and a positive catch density, both predicted from the dataset and the best models by using the function ‘predict.gam’ in the package ‘mgcv’ (Wood, 2011) in R. The s-density was finally averaged by year and grid cell.

We also obtained nominal CPUE for each fishing grid cell from the official logbooks of offshore trawl vessels to test the plausibility of the s-density by testing the concordance between the predicted and actual fishing grounds. Furthermore, we compared s-density with nominal CPUE to identify new potential fishing grounds. Because the northern shrimp catch by large offshore trawlers in the northern and central Sea of Japan (>136.5°E) is minuscule, we focused on the region to the west of 136.5°E, namely regions A‒F. The logbooks include the northern shrimp catch (kg), the number of tows, the date of each cruise and the location of the representative fishing grid cell.

We considered data only with a northern shrimp catch >0 kg from May to July in 2013–2022 (the same period as used in the model construction). May to July corresponds to the peak season for northern shrimp catch along the coastline of Honshu in the Sea of Japan (21 March to 30 June and 1 September to 5 November, after and before the snow crab season, respectively). The nominal CPUE per fishing grid cell was then calculated as total catch/total number of tows (kg/tow) by year based on the extracted dataset.

Results

Northern shrimp SDMs

The full SDM including seven variables (year, region, PC1, depth, slope, TPI and the interaction of region and depth) and a partial model excluding year were selected as the best first- and second-step models, respectively, based on AIC (Table 2). The best first- and second-step models had high explanatory power, accounting for 68% and 78% of the deviance, respectively (Table S1). The AUC from the ROC analysis of the first-step binomial model was 0.97 (Figure S3). Swets (1988) divided model performance into three categories of ‘poor’ (0.5‒0.7), ‘medium’ (0.7‒0.9) and ‘high’ (0.9‒1.0) based on AUC. Under this categorization, our first-step model shows ‘high’ performance. A Q–Q plot for deviance residuals of the second-step model showed a distribution similar to a straight line, and the density of northern shrimp was found not to deviate significantly from a Gaussian distribution (Figure S4).

Table 2. Summary of generalized additive models

+ indicates the factors included in each model. PC1, first principal component; TPI, topographic position index; AIC, Akaike information criterion; ΔAIC, change in AIC as explanatory variables are removed from the best model.

Comparison of ΔAIC among candidate models indicated that the interaction of region and depth was the most important explanatory variable for models of occurrence and density, followed by PC1, year, depth, region, TPI and slope in the first-step model, and region, slope, PC1, TPI and depth in the second-step model (Table 2). The effect of each explanatory variable in the presence/absence and density models are shown in Figures 3‒6, respectively. The effect of year on northern shrimp presence was greater during the years 2014−2022 than in 2013 (Figure 3). The effect of region A was higher than that of region E in the first-step model (Figure 3), and higher than that of regions C–G in the second-step model (Figure 5). Estimated averages of the bottom sediment effect were positive at PC1 < 0.0, depth >282 m, slope <0.7° and TPI < 0.0 (but almost flat) in the first-step model (Figure 3), and at PC1 < 0.0 (with minor exceptions), depth >345 m, slope <0.8° and TPI from −17 to 0.9 in the second-step model (Figure 5). The effect of depth by region varied across both models. In the first-step model, the depth effect tended to be higher at greater depth but with some exceptions (e.g. region A) (Figure 4). Also, the trends were sometimes unclear due to large uncertainties (e.g. regions C and D). In the second-step model, complex patterns were observed for each of the regions but generally, stronger effects were observed at deeper depths (Figure 6). Importantly, the strongest depth effect occurred at a different depth in each region.

Figure 3. Response plots for the best generalized additive models (first-step models) showing the effect of year, region and environmental variables on catch probability from 2013 to 2022 along the coast of Honshu Island in the Sea of Japan. In the bar graphs, error bars show standard error. In the line graphs, rug plots along the x-axis indicate the number of observations. Solid lines show response curves and shaded regions show 95% confidence intervals. PC1, first principal component; TPI, topographic position index.

Figure 4. Response plots for the best generalized additive models (first-step models) showing the effect of the interaction between region and depth in regions A–I on catch probability from 2013 to 2022 along the coast of Honshu Island in the Sea of Japan. Rug plots along the x-axis indicate the number of observations. Solid lines show response curves and shaded regions show 95% confidence intervals.

Figure 5. Response plots for the best generalized additive models (second-step models) showing the effect of region and environmental variables on positive catch density from 2013 to 2022 along the coast of Honshu Island in the Sea of Japan. In the bar graph, error bars show standard error. In the line graphs, rug plots along the x-axis indicate the number of observations. Solid lines show response curves and shaded areas show 95% confidence intervals. PC1, first principal component; TPI, topographic position index.

Figure 6. Response plots for the best generalized additive models (second-step models) showing the effect of the interaction between region and depth in regions A–I on positive catch density from 2013 to 2022 along the coast of Honshu Island in the Sea of Japan. Rug plots along the x-axis indicate the number of observations. Solid lines show response curves and shaded areas show 95% confidence intervals.

Northern shrimp spatial distribution

Estimates of s-density indicate that northern shrimp is widely distributed across the continental shelves and slopes along the coast of Honshu Island in the Sea of Japan (Figure 7A), and that interannual differences in s-density trended to be fairly small (Figure S5). s-density was high in fishing grid cells that include the 500 m isobath, and low in shallower grids (Figure 7A). In addition, the fishing grid cells in the northern area (east of Noto Peninsula) with steep slopes and including both the 200 and 500 m isobaths tended to have lower s-density (Figure 7A). Among the nine regions, regions A and C had numerous low s-density grid cells, whereas the other regions included higher s-density grid cells (Figure 7A). In the major offshore trawl fishing areas (regions A‒F), fishing grid cells with higher s-density were generally characterized by higher nominal CPUE (Figure 7A, B).

Figure 7. Density and catch per unit effort (CPUE) of northern shrimp along the Sea of Japan coast of Honshu Island. (A) Standardized density (s-density) in each fishing grid cell estimated with the best occurrence and density models. Grid cells on Yamato Bank, Oki Bank and Hakusanse Bank were excluded from the analysis. (B) Nominal CPUE in each fishing grid cell based on logbook data of offshore trawl fisheries from May to July of 2013 to 2022 (values are averaged across years). The nominal CPUE on Yamato Bank is not shown. Dotted and dashed lines show 200 and 500 m isobaths, respectively. Uppercase letters show administrative regions with borders shown in black lines.

Survey-based s-density was positively correlated with logbook-based nominal CPUE by year and fishing grid cell (r = 0.16‒0.63) in most cases, and correlation coefficients ranged from 0.339 to 0.628 in all years, which corresponds to moderately positive correlation (r = 0.3‒0.7; Ratner, 2009), except for 2014 and 2022 (Figure 8).

Figure 8. Scatter plots of s-density and nominal CPUE in each year for each fishing grid cell. Numbers show Pearson correlation coefficients (r) between s-density and nominal CPUE.

Discussion

Our results demonstrate that the distribution of northern shrimp along the coast of Honshu Island in the Sea of Japan is influenced by a range of environmental parameters (Table 2). We successfully built delta-GAM models based on data from two bottom trawl surveys. Among all explanatory variables in models of both occurrence and density, the interaction of region and depth was the most important variable with the highest ΔAIC (Table 2). Because the depth at which the effect was strongest differed among regions, we assume that regional deviations from the ‘standard’ ontogenetic depth distribution are important in explaining northern shrimp distributions in the Sea of Japan.

Habitat preferences of northern shrimp in the Sea of Japan

The environmental characteristics of northern shrimp habitat as inferred from the present study are largely consistent with previous reports on the species (Table 3). Our findings on preferred seafloor topography generally agreed with previous reports, but the estimated habitat depth and preferred bottom temperature in the Sea of Japan were deeper and colder, respectively, than in other regions, as was previously mentioned by Ito (1976) and Sadakata (2000b) (Table 3). These habitat differences could explain the differences in growth and breeding of this species in the Sea of Japan as compared to other regions – northern shrimp in the Sea of Japan show slower growth rates, longer lifespans, longer ovigerous periods and a biennial spawning pattern (Sadakata, 1999, 2000a, 2000b). Each of the key environmental factors is discussed below.

Table 3. Habitat characteristics of northern shrimp in various regions

Btemp, bottom temperature; Bsed, bottom sediment; TPI, topographic position index; range, range of occurrence; presence, range of highest presence; density, range of highest density; ND, no data.

The estimated depth effect suggests that northern shrimp in the Sea of Japan are abundant on seafloors deeper than 283 m (Figures 3 and 5), where bottom water temperatures are below 2.6 °C. These results are generally consistent with catch records from previous studies in the Sea of Japan (Kojima et al., 1969; Ito, 1976). This depth boundary is roughly consistent with that of Japan Sea Proper Water, a homogeneous water mass that occurs at depths deeper than 200 m and is characterized by extremely low temperature (0‒1 °C, Uda, 1934). Because northern shrimp is cold adapted, this suggests that populations in the Sea of Japan prefer thermal environments with little seasonal change under the strong influence of Japan Sea Proper Water.

Regional variations in the effect of depth on northern shrimp occurrence (Figure 4) and biomass (Figure 6) can partly be explained by ontogenetic migration towards deeper areas (Sadakata, 2000b). For example, differences in the depth of maximum effect on biomass in regions B (ca. 300 m) and C (ca. 400 m) could reflect the greater horizontal distance required for ontogenetic migration in region B, where the continental slope is much broader than in region C (Figure 6). Also, regional differences in the relationship between depth and environmental characteristics (e.g. bottom temperature) could also play a role. Further data collection is needed to test these hypotheses (see also the Limitations section below).

Our results also show that northern shrimp prefers muddy rather than sandy bottoms (Figures 3 and 5). This could reflect the species’ tendency to feed not only on benthic organisms such as crustaceans but also on organic matter contained in the bottom sediment itself (Rice et al., 1980). The density of the congeneric P. borealis is positively correlated with the organic matter content of bottom sediment (e.g. Wigley, 1960; Haynes and Wigley, 1969). Moreover, Saito and Ikehara (1992) have noted that organic matter content is generally positively correlated with mud content in Japanese waters. If organic matter is an important food source for northern shrimp, this would explain its high abundance on muddy sediment containing high levels of organic matter. However, it is also possible that northern shrimp does not feed selectively on bottom sediment and that any consumed sediment is incidentally ingested alongside prey items (Rice et al., 1980). Further study is needed to determine the importance of organic carbon as a food source for this species.

The results also indicate that northern shrimp prefers gently sloped to flat seabeds and flat valley bottoms as habitat (Figures 3 and 5). This is consistent with previous reports. Sadakata (2000b) reported that northern shrimp habitats off Noto Peninsula occur on gentle slopes. Additionally, according to submersible observations, this species is more abundant in natural shallow holes and on ploughed seafloors (plough tracks) to improve seabed environment (e.g. Yasuda, 1987; Kasutani, 1993).

Year was only selected as a variable in the first-step model and was not included in the second-step model (Figures 3 and 5). We hypothesize that northern shrimp density is already reaching the limit of ecological carrying capacity in the study region and that annual fluctuations of stock size reflect expansions or contractions in range, thus only affecting occurrence probabilities. The standardized CPUE of the offshore trawlers have been at a historical high since 2015 (Figure 2; Naito et al., 2024), and our results may reflect the strong stock status of northern shrimp in this area. Further biological information such as annual fluctuations in food resources for this species, competition and predation is needed to prove our hypothesis.

In addition, the importance of some variables differed between the two models. For instance, PC1 (a metric of grain size) was important for northern shrimp presence/absence but had little effect on density (Table 2). We hypothesize that there is a threshold of grain size beyond which northern shrimp cannot inhabit, and that this directly affects the presence or absence of the species. In fact, the effect of PC1 in the first-step model drastically changed around a value of −0.1 (the inconsistency of the trend for PC1 > 0.3 can be attributed to data deficiency) (Figure 3). Bottom sediment type has a strong influence on the distribution of many benthonic marine organisms via direct effects (e.g. physical characteristics), indirect effects (e.g. food availability) or a combination of both (Wigley, 1960). As we detailed above, northern shrimp can ingest organic matter in the bottom sediment (Rice et al., 1980), and the indirect effect of this organic matter ingestion could explain the effect of sediment type on the presence/absence of the species. Further biological information (e.g. predator–prey interactions and the habitat preferences of other marine organisms) is needed to uncover further conflicts between the model variables (i.e. region and slope).

Finally, we note that there is some uncertainty around the species status of northern shrimp in the Sea of Japan. Komai (1999) has reported that shrimp in the Sea of Japan have a different number of dorsolateral spines on the telson than do shrimp in the Bering Sea and Alaska. Future studies involving genetic analysis of the various populations may be helpful in determining whether these differences are intra- or interspecific.

Predicting potential fishing grounds of northern shrimp along the Honshu coastline of the Sea of Japan

Survey-based s-density was generally positively correlated with logbook-based nominal CPUE by year and fishing grid cell (r = 0.16‒0.63), suggesting that our models reflect the empirical distribution of northern shrimp in this area fairly robustly (Figure 8). Furthermore, we identified two areas as possibly underutilized fishing grounds, namely, high s-density fishing grid cells that include 500 m isobaths in region A and similarly high s-density cells including 500 m isobaths north of 37°N in region B (Figure 7A, B). These areas are not heavily used by the large offshore trawl fishery. Because the northern shrimp fishery in region A is operated by trawlers from Tottori and Hyogo prefectures (facing regions C and D, respectively), it may be more cost effective for these vessels to operate in fishing grounds closer to these prefectures (e.g. regions B and C) to adapt to sharply rising fuel costs. Moreover, the potentially underutilized area of region B is part of the provisional measures zone between Japan and South Korea, and trawl operation in this region is challenging because of the presence of numerous abandoned pots and gill nets from South Korean fisheries. However, these areas could become viable new fishing grounds pending the resolution of fuel cost issues and the implementation of cooperative fishery management by Japanese and South Korean fisheries authorities.

Limitations

Although the present study provides robust results supported by empirical fishery data, there are some key limitations to our approach. First, the underlying factors affecting the variables in our models are not fully accounted for. Although we suggest that the interaction between depth and region is associated with regional differences in ontogeny, we are unable to fully evaluate this hypothesis because of the lack of size composition data in our study. Similarly, although we suggest that some differences in the region variable reflect differences in abundance among regions (e.g. regions B and C; Figure 5), some of this difference could also be caused by differences in trawl gear efficiency as described in the Materials and methods. This is because the stronger region effect in regions H and I than in the region A in the second-step models (Figure 5) can be attributed to differences in the mesh sizes of beam and otter-trawl nets (25–33 vs 45–50 mm). The inclusion of size composition data, gear efficiency and other possible environmental variables (e.g. sea surface temperature and productivity) may help to more fully explain the links among our explanatory variables and possible underlying factors.

Finally, the distribution of northern shrimp in the Sea of Japan could vary by month because the species undergoes seasonal vertical migration (Sadakata, 2000b). The effect of this migration on northern shrimp distribution in our study area is expected to be relatively small because there is a little movement of females in particular during the months (May to July) when the trawl surveys were conducted (Sadakata, 2000b). Because of this, we did not consider seasonality in our delta-GAMs, but future trawl surveys conducted during other parts of the year could help improve estimation of habitat preferences in this species.

Conclusion

The distribution of northern shrimp along the coast of Honshu Island in the Sea of Japan was modelled by using delta-GAMs. The results indicate that the species distribution is influenced by environmental factors including bottom sediment composition, depth, slope and TPI along with year and region. Our study provides the first quantitative analysis of the habitat preferences of this commercially important species in this region. Additionally, we reviewed previous reports of the species’ distribution in detail. Finally, we predicted the locations of potential fishing grounds and provided an accurate spatial distribution map. These results could be used to conduct more accurate stock assessments of this species and inform more effective stock management strategies.