Spatial dependency in abundance of Queen conch, Aliger gigas, in the Caribbean, indicates the importance of surveying deep‐water distributions

To evaluate how the spatial distribution of a heavily exploited marine gastropod (i.e. Queen conch) varies in response to a number of known biotic and abiotic variables within and between study areas that vary in environmental conditions.


| INTRODUC TI ON
Queen conch (Aliger gigas) is an economically and culturally important marine gastropod (Appeldoorn, 1994;Brownell & Steveley, 1981) found throughout the wider Caribbean region and the southern Gulf of Mexico. Its biological characteristics combined with high cultural and economic value make the species vulnerable to overfishing (Appeldoorn et al., 2011). The species has been over-exploited throughout large parts of its geographical distribution range (Acosta, 2002), resulting in concerns for the species' future and its listing in Appendix II of the Convention on International Trade of Endangered Species of Wild Fauna and Flora in 1992.
Although conch are most common in depths <25 m (Ehrhardt & Valle-Esquivel, 2008;Weil & Laughlin, 1984), the species can be found to depths of 60 m (Randall, 1964). The depth range of queen conch is believed to be restricted based mostly on light attenuation limiting their photosynthetic food sources (e.g. filamentous alga; Creswell, 1994;Randall, 1964;Ray & Stoner, 1994). Queen conch move to deeper waters with age and size (Randall, 1964;Weil & Laughlin, 1984). Unfished deep-water adult populations are believed to provide significant recruitment to shallow-water stocks and are considered critical spawning stock refugia (Appeldoorn, 1997).
Deep-water populations (>25 m) have been found in the Caribbean (e.g. the Bahamas, Belize, Martinique, Turks and Caicos, Jamaica, Puerto Rico [Appeldoorn, 1997;Berg, 1975;Berg & Olsen, 1989;García-Sais et al., 2012;Reynal et al., 2009;Stoner & Sandt, 1991]), but these populations have rarely been studied in detail. The depth range of conch is affected by fishing: where adult densities in shallow depths are reduced, the remaining conch are often found at greater depths where they are more difficult to reach by fishermen.
Adult conch can be found in a range of habitats (Stoner & Schwarte, 1994) with a preference for sand/algal flats, but they are also found on hard-bottom habitats such as coral rubble (Acosta, 2001;Stoner & Davis, 2010;Torres-Rosado, 1987). They are rarely found on soft bottoms (silt/mud) or in areas with high coral cover (Acosta, 2006), likely due to restrictions in their movement in such habitats (Dujon et al., 2019). Tidal channels with high water flow are of ecological importance to conch (Stoner, 2003), with increasing abundance of conch of all sizes closer to such tidal channels . Due to the positive association of conch abundance with tidal channels, other areas with high water flow, such as around the edge on offshore banks (Hamner & Hauri, 1981), are also expected to have high conch abundance.
Most published studies on the relationships between environmental variables and conch distribution and abundance were based on data collection through dive surveys. However, due to the safety limitations of surveying using scuba, areas below 20 m were rarely surveyed, and areas below 30 m were generally excluded from biomass estimates (MRAG, 2013;Queen Conch Expert Workshop Group Report, 2012). In addition, dive surveys are logistically demanding and relatively expensive, in particular when surveying remote offshore areas (Queen Conch Expert Workshop Group Report, 2012). The rapid technical progress of video systems allows developing new cost-effective sampling tools to study benthic organisms, beyond depths safe for diving (Sheehan et al., 2010;Stevens, 2003). Video systems have shown to produce similar results as dive surveys and may be better at detecting smaller individuals (Cruz-Marrero et al., 2020).
A towed video method capable of accurately determining densities of live adult conch throughout the species' depth range (0-60 m) was developed by Boman et al. (2016) using a belt transect recorded by a laser calibrated towed video sled. This method facilitates abundance estimation and distribution of adult conch at deep, offshore locations. In addition, the method enables studying relationships between environmental variables and conch abundance throughout its depth range, which are, until now, poorly studied.
A complicating factor when studying the relationships between environmental variables and conch abundance is the tendency of conch to aggregate, especially during the reproductive season (Glazer & Kidney, 2004). A patchy distribution pattern with spatial autocorrelation among locations is thus to be expected Vallès & Oxenford, 2012). Therefore, it is important to use statistical tools that account for both autocorrelative spatial patterns and responses to various biotic and abiotic variables (Carroll & Pearson, 2000;Keitt et al., 2002;Zuur et al., 2017). In this study, Bayesian hierarchical spatial models, using integrated nestled Laplace approximation (INLA), were used to account for both aspects.
To predict adult conch abundance throughout the species' entire depth range, surveys were conducted at three study areas in the eastern Caribbean with varying environmental characteristics.
These surveys provide an opportunity to examine how biotic and abiotic factors influence abundance of adult conch. This study aimed to (a) evaluate how the spatial distribution of adult queen conch varies in response to a number of known biotic and abiotic variables within and between study areas, which vary in environmental conditions; (b) evaluate whether the patchy distribution and spatial dependency observed by previous studies is a general pattern in queen conch, and what likely causes this; and (c) test the hypothesis that significant densities of adult conch are a common occurrence in deep areas (>25 m), which would require a depth extension of conch surveys to be used for providing advice on fishing quotas. Maarten. The island is surrounded by a mixture of patchy, barrier and fringing reefs, interspersed with seagrass beds, sand channels and algal flats (Wynne, 2010). Survey transects (N = 132) were made between 31 August 2011 and 11 December 2015 in the waters around Anguilla, covering an area of approximately 940 km 2 in depths of 3 to 54 m, with 48% of transects (N = 63) conducted in the peak breeding season (May-September) ( Table 1; Boman et al., 2018). The area around Anguilla was divided into 5 × 5 km grid cells. Within each grid cell, 5-6 survey transects were dispersed randomly (N = 117).

| Study area and survey design
In addition, within the survey area a total of 15 fixed transects for long-term monitoring of conch abundance in nearshore shallow areas around Anguilla were added.  (Debrot et al., 2014). Rocky reef areas are limited to the southern and southwestern shelf areas, whereas seagrass beds are confined to the west and north (Debrot et al., 2014). Survey transects (N = 167) were made in and in the close vicinity of the SNMP at depths of 6.5 to 45 m, between 11 June 2013 and 5 March 2014, with ca 62% of transects (N = 103) made in the peak breeding season (May-September) ( Table 1; Boman et al., 2018). Sampling locations in the waters surrounding St. Eustatius were selected by a random sampling design (Ehrhardt & Valle-Esquivel, 2008). Along the eastern edge, an actively growing coral reef zone is present Van der Land, 1977). The central, lagoon-like, part of the bank alternates between several habitat types: bare sand, patch reefs, rubble reefs, and sand or pavement (hard substrate) with macroalgae, sponges, gorgonians and/or coral structures (Lundvall, 2008;Meesters, 2010;Toller et al., 2010

| Abundance estimates
Adult conch (defined here as conch with a fully developed lip: [Boman et al., 2018]) abundance was estimated using two comparable methods: towed video (Boman et al., 2016) and standard belt transects (CRFM, 2013). Belt transects, using scuba, were primarily conducted in high-relief habitats (mainly reef habitats) in which the towed video array was not suitable to use (Boman et al., 2016). Total adult conch counts for each transect were calculated for the purpose of modelling conch abundance and distribution patterns. Density (number of conch / ha) estimates for each transect were also calculated but only used as a comparative measure to previous surveys and not used for the models.
Towed video transects at all three study areas were conducted following Boman et al. (2016), with a transect width of 1 m and a transect length between 330 and 806 m. All adult conch inside the transects (i.e. more than 50% of shell inside the transect) were counted, and life status was determined based on visual cues (Boman et al., 2016). Substrate (sand, rubble and reef) and macrobenthos (algae and seagrass) cover were determined for each transect by analysing 20 frames, with an equal spread in time over the transect.
In each of the 20 frames, 10 set dots, in a 2 × 5 pattern, were overlaid, and per dot, the underlying substrate and macrobenthos were identified. If a dot was blocked or the frame was blurry, the next frame in which substrate and macrobenthos for each dot could be determined was used. Substrate and macrobenthos cover for each transect were calculated based on the average total of dot points from the 20 frames. A hand-held GPS system was set to track the position, every 10 s, to follow the transect and calculate transect length accurately. Mean depth (m) was also recorded for each transect and determined by the vessel's depth sounder.
Belt transects at each survey site covered between 200 and 13,900 m 2 , with 94% of the transects being between 200 and 1,500 m 2 , and all adult conch within each transect were counted.
Differences in transect sizes were accounted for in the models by adding it as a log-transformed parametric covariate (sampling effort).
Life status was determined by visual inspection after turning over the animal. The divers estimated the percentage coverage of substrate (sand, rubble and reef) and macrobenthos (algae and seagrass) for each transect. The locations of transects were determined with a hand-held GPS (Garmin GPSMAP 78; Garmin Ltd., Olathe, KS, www. garmin.com). Mean depth was recorded for each transect as determined by the diver's computer. Measured substrate and benthos, as well as depth, distance to ocean, sampling effort, and longitude and latitude, were used as covariates for the models (Table 2).

| Data exploration and model selection
The three study areas had different geomorphology, with Anguilla and St. Eustatius being islands with a narrow shelf and Saba Bank being an offshore submerged bank. Due to these differences in geomorphology, tidal flow patterns, macrohabitat cover and size, relationships between conch counts and covariates were expected to vary among study areas. Thus, each site was analysed separately.
Data exploration was applied following the protocol described in Zuur et al. (2010). Due to the patchy distribution pattern of conch TA B L E 2 Summary and definition of covariates Vallès & Oxenford, 2012), the data were expected to contain high percentages of zeros. This suggests that zeroinflated models will be needed, but this is not necessarily so. If a covariate in the selected model can explain the zeros in the data, applying a zero-inflated model is not necessary. To determine whether a model can explain the zeros in the data, a simulation of data sets from the model (by sampling regression parameters from their posterior distribution) can be made, where the number of zeros for each simulated data set is counted and compared with the observed number of zeros in the original data (Zuur & Ieno, 2016). If these are comparable, then there is no need to extend the Poisson generalized linear models (GLMs) or generalized additive models (GAMs) to zeroinflated models (Zuur & Ieno, 2016).
Data exploration was carried out in order to identify (a) whether any of the covariates contained extremely large or small values, using Cleveland dot plots ( Figure S1), and (b) whether there was collinearity among covariates, using variance inflation factors (VIF) and pairwise scatterplots (  Figure S5). Thus, the covariates "log sampling effort" (total area [m 2 ] of the transect) and "sand cover" were modelled as fixed effects (Equations 2 and 3); the covariate "reef cover" was also modelled as a fixed effect: as categorical for Anguilla and as linear for St. Eustatius. In contrast, nonlinear patterns were found for the covariates "algal cover," "depth" and "distance to open ocean" (DOO) for all three study areas. These covariates were modelled as smoother functions.
Pearson residuals were used to determine the presence of spatial correlation, and GAMs with and without spatial dependency were applied (Zuur & Ieno, 2018) and compared using the Watanabe-Akaike information criterion (Watanabe, 2010), and the model with the lowest (best) WAIC value was selected.
All models were estimated using integrated nested Laplace approximations (INLA; Lindgren & Rue, 2015;Rue et al., 2009;Zuur & Ieno, 2018;Zuur et al., 2017). In INLA, the covariance matrix Ω of the spatially correlated random effects is quantified using the Matérn correlation function, which makes dependency a function of distance and a set of parameters of unknown value. To obtain these parameter values, a series of steps were carried out. First, a mesh was defined using a collection of small triangles. Greater numbers of triangles create a finer mesh and a more accurate INLA solution. At the node of each triangle, a value w k is estimated, using continuous domain stochastic partial differential equations (SPDE). These w k s are directly linked to the spatially correlated random effect and together form the spatial random field, which is visualized by standard interpolation techniques (Zuur et al., 2017).
The spatial correlation represents either real dependency originating from population interactions, spatially structured environmental controls, or missing covariates. Besides the posterior mean of the spatial random field (w k values), the posterior standard deviation for each w was estimated (based on a normal distribution), which can be used to infer which parts of the spatial random field are important (Zuur & Ieno, 2018), with smaller values of the standard deviation indicating higher importance.
One complicating factor for the Anguilla and St. Eustatius data was that the sampling locations were located around an island, and this meant that by default, the spatial correlation crosses land. This is not problematic if the spatial correlation only represents a missing covariate that affects study areas on both sides of the island, but if it represents real spatial dependency, then the model needs to take into account that conch do not cross land. This was done using barrier models (Bakka et al., 2019), which ensure that the spatial correlation matrix does not cross land. Our choice for the barrier was based on biological relevance, and not only on statistical grounds, because barrier models tend to have larger WAIC values than standard models (Zuur & Ieno, 2018). Simulations of 10,000 data sets from the model with the best fit were conducted, to determine whether the models met assumptions and could cope with the number of zeros in the data (Zuur & Ieno, 2016).
Equation (1)  where Conch i is the count of conch at site i, Covariate i represents parametric effects, and f(Covariate i ) is a smoothing function (spline) (Zuur & Ieno, 2018). The term u i in the model is a spatially correlated random intercept that is normally distributed with mean 0 and covariance matrix Ω. In summary, the following steps were taken to arrive at the final models: (1) Data exploration for outliers and collinearity; (2)  (1) fitting Poisson GAMs with and without spatial correlation using INLA (a barrier model was used where appropriate); (9) comparison between models with and without spatial dependency using the Watanabe-Akaike information criterion, and selecting the model with the lowest (best) WAIC value; and (10) Final determination whether the models with the best fit met assumptions and could cope with the number of zeros in the data, through data simulation.
Eustatius and Saba Bank, respectively. Saba Bank had the highest overall mean conch density (126/ha) of the three study areas (Table 3). Conch were found throughout the Saba Bank from 16 m to 50 m depth, except along the eastern edge with high reef cover where no conch were found (Figures 2 and 3). Anguilla had the lowest overall mean density (27/ha), and conch were found around the entire island with higher densities on the western and eastern sides of the island at depths between 4 and 50 m (Figures 2 and   3; Table 3). St. Eustatius had a mean density of 62/ha, and conch were concentrated to the south-western and south-eastern parts of the island ranging from 11 to 45 m depth (Figures 2 and 3; Table 3). The initial Poisson GLMs fitted were overdispersed for all three study areas, and the subsequent scatterplots of conch counts versus the covariates indicated some nonlinear patterns ( Figure S5).
Furthermore, the initial GAMs, without spatial correlation, could not cope with the percentage of zeros in the data and were also overdispersed. Model validation showed that the Pearson residuals were spatially correlated for all sites.
For the three study areas, the models with spatial correlation showed lower WAIC values than the models without spatial correlation (Table 4). Model validation indicated that the spatial Poisson GAMs for Anguilla and St. Eustatius did not contain any remaining spatial correlation in the Pearson residuals, and the models could cope with the zero inflation ( Figure S6). Similarly, the spatial Poisson GAM for Saba Bank did also not contain any remaining spatial correlation in the Pearson residuals and the model could cope with the zero inflation ( Figure S6).

TA B L E 3 Summary statistics of adult conch densities (mean with 95% confidence interval [CI]), at the three study areas
The following GAMs with spatial correlation were therefore In Equations (2) and (3), the terms f 1 () to f 3 () are smoothing functions (splines; Zuur & Ieno, 2018), while the other terms are parametric effects. The term u i in the model is a spatially correlated random intercept that is normally distributed with mean 0 and covariance matrix Ω.  (Table 5). The spatial random field showed two areas with relatively high values and also 1 area with relatively low values of conch ( Figure 2) with spatial correlation up to 1.2 km. Where the value of the spatial random field was around 2.5 (dark red areas), this indicated an exp(2.5) ≈ 12.2 times larger number of conch than the overall mean, due to spatial dependency. Conversely, in the dark blue areas (values ca. −2), the model indicated conch values that were a factor exp(−2) ≈ 7.3 smaller than the overall mean. The posterior standard deviations were around 0.8, which means that ws larger than 1.6 or smaller than -1.6 were deemed important (based on the normal distribution; Figure S7).

St. Eustatius
In St. Eustatius, all smoothers were important (i.e. algae, DOO and depth). The algal cover smoother showed a bell-shaped pattern that caused lower values of conch at low algal cover (0-0.12) and higher values of conch at medium cover (0.25-0.62) (Figure 3). The depth smoother showed a similar bell-shaped pattern, which caused lower values of conch between the 0-and the 14-m interval and marginally higher values of conch in the 18-to 27-m interval. The DOO smoother showed lower values beyond 2.8 km (Figure 3). Sampling effort, sand and reef cover all had an important negative effect on conch counts (Table 5). The spatial random field (Figure 2) showed several areas with relatively high and relatively low values of conch, with a spatial correlation up to 1 km. When the spatial random field was around 3.5 (dark red areas), this indicated an exp(3.5) ≈ 33 times larger number of conch than the overall mean, due to the spatial dependency. Conversely, in the dark blue areas (values ca. −2.5), the model indicated conch values that were a factor exp(−2.5) ≈ 12.2 smaller than the overall mean. The posterior standard deviations were around 1, which means that ws larger than 2 or smaller than −2 can be deemed important ( Figure S7).

Saba Bank
Depth was important and caused marginally higher values of conch between the 22-and 27-m interval (Figure 3). DOO was also important and caused lower values of conch in the 0-to 1-km interval and higher conch values within the 3.5-to 7-km interval (Figure 3). Algal cover had no influence on conch abundance (Figure 3). Neither of the two fixed covariates, sampling effort (i.e. transect length) and sand, were important (Table 5). The spatial random field presented in Figure 2 showed several areas with relatively high (red areas) and relatively low (blue areas) values of conch with spatial correlation of up to 7 km. When the spatial random field was around 3.5 (dark red areas), this indicated an exp(2.5) ≈ 12.2 times larger number of conch than the overall mean, due to the spatial dependency. Conversely, in the dark blue areas (values ca. −2), the model indicated conch values that were a factor exp(−2) ≈ 7.3 smaller than the overall mean. The posterior standard deviations were around 1, which means that ws larger than 2 or smaller than −2 can be deemed important ( Figure S7).

| D ISCUSS I ON
This study showed distinct spatial distribution patterns of adult conch, which occurred in patchy distributions, with areas of high and low abundance. The patchy distribution of conch was caused by spatial dependency (Figure 2) and had a maximum magnitude of range up to 7 km for Saba Bank. The magnitude of range for the spatial dependency for Anguilla and St. Eustatius was smaller, up to 1.5 and 1 km, respectively, which could have been (partly) caused by the smaller survey area and the presence of land barriers at these locations. The patchy pattern appeared at all survey locations, and past studies have shown that including spatial dependency is important, because ignoring it may point to radically different conclusions (Keitt et al., 2002).
The observed spatial dependency could represent either real dependency, originating either from endogenous processes (population biological interactions) or exogenous processes (spatially structured environmental controls) (Planque et al., 2011), or from other missing covariates, such as fishing pressure . A previous study linked conch's patchy distribution pattern to habitat with higher abundances of conch in algal habitats and lower abundance in coral habitats (Vallès & Oxenford, 2012). However, our study did not support this link. Here, spatial dependency from the model was not captured by the habitat covariates. Endogenous processes in conch that influences conch to aggregate have been identified. Stoner and Ray-Culp (2000) demonstrated an Allee effect in conch where reproductive activity begins to decline as densities of conch decrease below approximately 200 conch/ha. Consequently, conch tend to aggregate, which can at least partially explain the spatial dependency observed in conch abundance during the reproductive season (approximately 50% of transects were conducted in the reproductive season). However, due to their limited home range (ca. <1-6.5 ha) and daily movement patterns (ca. 11 m/day) (Delgado & Glazer, 2007;Doerr & Hill, 2013;Glazer et al., 2003;Stieglitz & Dujon, 2017) it is likely that the aggregation persists outside of the breeding season.
Still, currently unknown missing factors cannot be ruled out and the discovery and addition of such factors could potentially shed further light on the cause of spatial dependency of adult conch.
The Bayesian hierarchical spatial models for St. Eustatius also showed fewer conch at shallow depths (0-15 m) and a small increase in conch abundance just below 20 m for St. Eustatius and Saba Bank (Figures 2 and 3). Environmental factors likely influenced the lack of conch in shallow depths (0-15 m) at the three locations.
Fishing has been known to change the depth distribution of conch and can shift conch distributions to greater depths depending on the methods used (Stoner & Schwarte, 1994). Although Saba Bank has had a complete moratorium on the conch fishery since 1994 (Hoetjes & Carpenter, 2010) and can thus be expected to have conch population, which is close to its natural distribution, it also has a shallowest point of ca. 15 m, and thus, the conch population has a natural deeper distribution in comparison with the literature (Ehrhardt & Valle-Esquivel, 2008;Weil & Laughlin, 1984). At St. Eustatius, conch are caught, although at relatively low amounts (3% of total adult population; Meijer zu Schlochtern, 2014). However, the lack of high densities of conch at shallow depths was likely due to shallow habitats being either unsuitable (high-relief areas with corals reef; Debrot et al., 2014) or exposed to the elements that deter conch from settling in the shallow areas around the island. Contrary to St. Eustatius, the shallow areas around Anguilla were often more sheltered and the lack of high abundance of conch in these seemingly suitable areas is likely to be at least partially driven by fishing pressure, which was approximately 6% of total adult population (Kuramae-Izioka, 2016).
In contrast to shallow waters, high densities of conch (> 250 /ha) were found in waters >30 m depth at all three study areas. Although conch are known to be most common in depths <25 m (Ehrhardt & Valle-Esquivel, 2008), high densities (ca. 300/ha) of adult conch during the peak reproductive season have been found in deep-water habitats (García-Saiset al., 2012). Similarly, the current study found densities >300 adult conch/ha at Saba Bank at a depth of 40 m both during and outside the peak reproductive season (Boman et al., 2018) (Acosta, 2006). Sand cover had an important negative effect on conch abundance (Table 5), which is in contradiction with the general notion of sand habitats being recognized as suitable conch substrates (Acosta, 2001;Stoner & Davis, 2010;Torres-Rosado, 1987). In this study, conch were also present in high abundance (>250/ha) in areas with high sand cover (>90%). Therefore, it is not likely that the results indicate a general negative effect with increased sand cover on conch distribution. Instead, within suitable substrates (e.g. sand and rubble) the abundance and distribution of conch is likely more influenced by other factors such as depth and algal cover, and through the natural patchy distribution of conch.
In the current study, high levels of algal cover were mostly associated with low conch numbers. However, algal cover included all macroalgal and cyanobacterial mats, while further in-depth analysis down to species level was not conducted. As the majority of the algal species in the surveys were calcified or possessed chemical defences (e.g. Halimeda spp, Dictiota spp, Caulerpa spp), they are non-palatable for most marine species and thus were not likely to be a significant source of nutrition for conch (Erickson et al., 2006;Hay et al., 1987;Pereira & de Gama, 2008;Pereira et al., 2002). Associations between queen conch abundances and macroalgal cover were thus likely not based on foraging behaviour. However, the results from this study (Figure 3) indicate that there is a maximum threshold for algal cover and possibly an optimal level of macroalgal cover for conch. We hypothesize that macroalgal cover to a certain level may be beneficial for adult conch due to a potential increased in available food sources (macrophyte epiphytes). However, when algal cover reaches a certain level (ca. 0.6-0.7) other sources of nutrition (e.g. benthic diatoms) may be reduced by light competition with macroalgae (Hill, 1996;Yang & Flower, 2012), which cannot be fully compen- which is an actively growing coral reef zone ( Van der Land, 1977). No conch were found in this reef zone, unsurprising given that conch are usually not found in coral reef habitat (Acosta, 2006). As reef cover was not included as a covariate in the model for Saba Bank, this factor could not provide explanatory effect on the low amount of conch in this habitat shown by the DOO covariate. Tidal current flow has been found to govern the distribution of conch in The Bahamas where adult conch were positively associated with tidal channels with high tidal flow and found in higher abundance closer to the tidal channel . A significant inverse association was found between distance away from the mouth of a tidal channel and conch size and age, suggesting that older and larger animals migrate towards higher flow (Kough et al., 2019). Offshore banks (e.g. Saba Bank) have a tidal current flow, which is strongest at the edge of the bank, closer to the open ocean and weaker in the centre of the bank (Hamner & Hauri, 1981). The higher abundances of conch found towards the edge of the Saba Bank in this study are thus likely an effect of tidal flow on the distribution of conch and support previous suggestions of an ecological importance of water flow for conch and their distribution pattern Stoner, 2003). Small islands such as Anguilla and St. Eustatius with a narrow shelf have a more complex water flow pattern, not related to distance to ocean and which is dependent on a range of factors (e.g. shape and size of island mass, direction of water flow, topography of area around the island) creating areas of high and low water flow downstream of the island (Hamner & Hauri, 1981). Therefore, as seen for Anguilla and St. Eustatius, it is not expected that abundance of conch should concur with distance to open ocean, and more complex measurements of water flow patterns will likely be necessary to elucidate whether abundance and distribution of conch are influenced by areas of high and low water flow at such locations.
Spatial dependency and a patchy distribution of conch, independent of the biotic and abiotic factors tested (e.g. substrate, macrophyte cover, distance to ocean), appear to be generally applicable to conch populations as this pattern was visible in areas with different geomorphology, size and habitat homogeneity. The lack of strong generic relationships between biotic and abiotic factors and adult conch abundance and distribution at all study areas is likely partly due to this spatial dependency and different locationspecific factors (e.g. patterns of water flow), which affect different phases of conch life history. One such factor could be larval transport, which shapes the metapopulation of queen conch in the wider Caribbean (Truelove et al., 2017). Also, smaller-scale patterns of larval transport appear to be influencing local population demographics (Kough et al., 2019). Larval connectivity is a biologically relevant covariate, which was not considered in this study, and it is thus unclear whether and to what extent larval connectivity could have played a role in the spatial patterns we found. Furthermore, conch are often unlikely to position themselves optimally in accordance with important factors due to trade-off situation, which was likely the situation seen on the Saba Bank where unsuitable habitat (high-relief reef) likely prevented conch to move to the far most edge of the bank where the highest water flow is found (Hamner & Hauri, 1981).
We demonstrated that substrate and depth can often predict conch abundance and distribution patterns. However, specific percentage coverage of suitable substrate (i.e. sand and rubble) and specific depth in the approximate range of 0-45 m are not reliable factors for predicting adult conch abundance and distribution patterns. The results from the current study showed high densities of conch (>300/ha) at depths >30 m and indicate that the current notion of a most common depth range of conch (<25 m) (Ehrhardt & Valle-Esquivel, 2008;Weil & Laughlin, 1984) may need revision and that deeper areas should not be excluded from conch surveys when present. This general notion of a preferred depth range could in part be due to conventional survey methods being impractical for surveys beyond this depth range and thus be biased towards deeper areas (>25 m). It could also in some areas be a by-product of longterm over-exploitation that has restricted the remaining population to less-accessible deeper habitats .
Another factor that can complicate the prediction of natural abundance and distribution patterns of conch is fishing pressure, which most visibly disrupts the natural depth distribution patterns of conch when shallow more easily caught conch are usually tar- influence is likely to be limited here. A direct comparison between our projections of the spatially random field and subsequently collected data from fishers would support or refute the role of fishing.
In the meantime, we can show where variability occurred to inspire future works, identify relevant factors and track annual variability.
Analysis of spatial dependency gives researchers an additional tool to resolve distribution patterns and to identify potentially missing relevant factors influencing such patterns. Conch populations are managed over many spatial scales throughout their range. Our results identify key components of conch distributions that are within the means of many management entities to consider. We recommend that surveys operate over a range of depths while sampling multiple macrobenthos covers. Results should be analysed with spatial dependency to disentangle local factors from range-wide tendencies. While the management strategy may be different for a small-island park compared with a nation whose jurisdiction includes banks, islands and continental shelf, we show that the same factors are likely to drive conch distributions even if they have different levels of effect. Here, we provide methods to evaluate different statistical models, assess covariates and describe spatial patterns of queen conch abundance.

ACK N OWLED G EM ENTS
We would like to acknowledge the contributions of the staff from the

Department of Marine Resources in Anguilla, St Eustatius National
Parks Foundation and Saba Conservation Foundation, who assisted in the data collection. We thank the three institutions for their support of this work.

CO N FLI C T O F I NTE R E S T
The authors declare that there is no conflict of interest.

PE E R R E V I E W
The peer review history for this article is available at https://publo ns.com/publo n/10.1111/ddi.13392.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data used to derive the distribution patterns of adult conch using