Using Benthic Indices to Assess the Ecological Quality of Sandy Beaches and the Impact of Urbanisation on Sandy Beach Ecosystems

Abstract

As the global population continues to grow, sandy beaches, one of the most valuable ecosystems, have been widely impacted by human activities. Therefore, to develop policies for the conservation and management of sandy beaches, the impact of human activities on sandy beaches must be accurately assessed. We used seven benthic indices to evaluate the ecological quality of sandy beaches in Anmyeon Island, Korea. However, these seven indices were found to either over- or underestimate their ecological quality. Moreover, despite incorporating beach morphodynamics into our study, these indices did not respond to the pressure of urbanisation on beaches. Given the suboptimal performance of benthic indices in reflecting the actual state of Korean beaches, our study indicates that beaches without human interference but with the same morphodynamics must be selected as control groups to further explore the effectiveness of these indices. This is critical for advancing our conservation efforts and managing sandy beach ecosystems under increasing human influence.

1. Introduction

With over 3.5 billion humans, accounting for more than 45% of the world’s population, living within 100 km of the coast, anthropogenic activities modify the coastal environment [1]. In coastal regions, particularly on sandy beaches, the escalating demand for economic development, resource extraction, and recreational activities, fuelled by population growth, is exerting unprecedented pressure on the world’s sandy beach environments [2]. Human pressures have led to changes in the ecosystems and physical environments of sandy beaches [3,4,5,6]. For example, activities such as trampling, the construction of breakwater barriers, and mechanical beach cleaning seriously impact the diversity of macrobenthic organisms on sandy beaches [3]. Rocking, scraping, and using off-road vehicles also affect the distribution of beach vegetation [7]. Moreover, humans are one of the primary contributors to the concentration of Escherichia coli on sandy beaches [8]. Therefore, the impact of human activities on sandy beaches must be accurately assessed to provide adequate policies for their conservation and management [9,10].

Sandy beaches, recognised as valuable ecosystems, have been the focus of considerable research interest in recent decades due to their provision of critical ecosystem services. These services include storm buffering, nutrient cycling, water purification, offering nursery areas for economically important species, and serving as feeding and breeding habitats for key species, such as endangered sea turtles and shorebirds [11,12]. In response to the need for managing and conserving these areas, various indices such as the beach evaluation index (BEI), beach quality index (BQI), conservation value index (CI), and recreation potential index (RI) have been developed [13,14,15]. However, these indices predominantly focus on the physical and social aspects of sandy beaches. Although no specific biotic index has been exclusively developed for sandy beaches, benthic indices based on the evaluation of macrobenthos have been widely applied to assess various coastal ecological environments [16]. To date, scientists have developed over 30 benthic indices based on various macrobenthic characteristics [16]. Most indices consider the tolerance of macrobenthos to organic enrichment or their life history traits for categorisation (e.g., AMBI, BENTIX, and BPI). In contrast, some indices consider the abundance or biomass of macrobenthos (e.g., ISEP and W-value). Additionally, some indices use multivariate approaches (e.g., M-AMBI and M-gAMBI). No single benthic index is the most suitable; the appropriate index must be selected based on the characteristics of the study area [17].

Specific benthic indices have effectively indicated the impacts of various human activities on marine ecosystems. For example, the AZTI marine benthic index (AMBI) and multivariate AZTI marine biotic index (M-AMBI) serve as effective instruments for identifying the effects of fish-farming activities along the coasts of Sardinia on benthic environments [18]; in the coastal waters of Fujian, the benthic polychaetes amphipods (BPA) index responded to dissolved inorganic nitrogen (DIN), copper (Cu), and arsenic (As) [19]. The M-AMBI is suitable for assessing the ecological quality of the eutrophic Yangtze River estuary [20]. Although benthic indices have achieved notable success in assessing the quality of marine ecosystems, only a few studies have examined benthic indices’ effectiveness for sandy beach ecosystems [21,22]. To the best of our knowledge, no study has assessed the efficacy of benthic indices on sandy beaches in Korea.

Anmyeon Island is located south of the Taean Peninsula, bordered by Cheonsu Bay to the east and the Yellow Sea to the west. The sandy beaches of Anmyeon Island have long been popular tourist destinations in Korea. To protect the marine environment of the Taean Peninsula, the South Korean government established the Taeanhaean National Park in 1978 [23]. However, due to reclamation projects and tourism development, Anmyeon Island has experienced substantial changes in its coastal topography [24]. Although some studies have proposed conservation and management plans for the sandy beaches of Anmyeon Island, these are predominantly based on their physical and social aspects [24,25]. Therefore, the conservation and management of beaches on Anmyeon Island must be considered from a biotic aspect.

In our study, we selected seven widely used benthic indices, the AZTI marine benthic index (AMBI), the benthic index (BENTIX), the benthic polychaetes amphipods index (BPA), the multivariate AZTI marine biotic index (M-AMBI), the abundance biomass comparison (W-value) [26,27,28,29,30], and the benthic pollution index (BPI) and Shannon–Wiener evenness proportion (ISEP), developed based on the Yellow Sea ecosystem, to assess the ecological quality of sandy beaches on Anmyeon Island, Korea [31,32]. In addition, we employed the urbanisation index to assess the pressure of human activities on these sandy beaches [33]. The aim of this study was to (1) evaluate the ecological quality and the degree of urbanisation of these sandy beaches, (2) assess the applicability of these seven benthic indices to sandy beaches on Anmyeon Island and benthic index responsiveness to human pressure, and (3) provide vital information for the conservation and management of the sandy beaches of Anmyeon Island to serve as a reference for future research on the ecological quality of sandy beaches in Korea.

2. Materials and Methods

2.1. Study Area and Sampling Design

This study focused on the sandy beaches of Anmyeon Island, Korea. The study area is defined by a latitude range of 36°24′ N to 36°34′ N and a longitude range of 126°19′ E to 126°22′ E (from the Baekssajang sandy beach in the north to the Baramarae sandy beach in the south) (Figure 1). The spring tidal range around Anmyeon Island is approximately 6.9 m, whereas the neap tidal range is approximately 3.1 m. The annual average temperature is 11.8 °C and the yearly average rainfall is 1223 mm on Anmyeon Island [34]. From our observations, human activities such as extensive trampling, beach urbanisation, and shellfish collection were seriously impacting all ten beaches on the island.

We selected the thirty stations on the ten sandy beaches of Anmyeon Island for our study. We observed that most tourists congregate in the middle of these beaches, near the paths leading to the seashore, or around the nearby squares. Consequently, we chose our sampling stations in areas where tourists frequently gathered. At each of the sandy beaches, three stations were selected for sample collection, located in the upper tidal region, middle tidal region, and lower tidal region (for example, for Baekssajang beach, S1 was the upper tidal region, S2 was the middle tidal region, and S3 was the lower tidal region). The sampling stations were surveyed from north to south (from Baekssajang to Baramarae), with three stations investigated on each beach, totalling 30 stations. We used a personal navigation assistant (GPSMAP 64S, Garmin Ltd., Lenexa, KS, USA) to determine the tidal regions. Additionally, macrobenthos were used as a critical indicator for identifying these regions. Notably, the upper tidal region of the beach was predominantly observed to be inhabited by ghost crabs. Conversely, the middle and lower regions exhibited the presence of Umbonium thomasi and burrows created by Chaetopteridae. The sandy beaches of Anmyeon Island attract many visitors annually, particularly in July and August. Therefore, we conducted our sandy beach surveys during the spring tides in June and September of 2023 to compare the impact of tourists on the sandy beaches’ urbanisation index and ecosystems.

2.2. Sample Collection and Processing

Macrobenthos were collected at each station using a 0.25 m × 0.25 m quadrat frame, with four replicates, resulting in a total sampling area of 0.25 m². The collection depth was 0.3 m, and the sampling occurred during the spring tides of June and September 2023. At the field site, samples were sieved through a 1mm mesh to isolate the macrobenthos. Macrobenthos were then preserved in a 4% neutral buffered formalin fixative (Avantor Inc., Radnor, PA, USA) for transportation to the laboratory. Concurrent with the collection of macrobenthos, approximately 200 g of surface sediment was collected at each station using a plastic spoon. Species were identified in the laboratory under a stereomicroscope (Olympus SZX-10, Olympus Co., Ltd., Tokyo, Japan). The macrobenthos were then counted, and their biomass was measured. The size of the sediment samples was measured employing dry sieve analysis, followed by the computation of the mean grain size using GRADISTAT software 8.0 [35]. To determine ignition loss (IL%), 30 g of dried sediment samples was placed in a muffle furnace (HY-800, Hwa Sueng Ind Co., Ltd., Busan, Republic of Korea) and heated at 550 °C for 2 h. The dissolved oxygen (DO), pH, and salinity of the surface seawater were measured using a multiparameter water quality meter (YSI-556MPS, YSI Inc., Yellow Springs, OH, USA) at each station. The soil temperature at each site was measured 13 cm below the surface using a waterproof thermometer (DT400, Summit Co., Ltd., Seoul, Republic of Korea). Based on a previous study, the morphodynamic type of the six beaches was identified [36].

2.3. Urbanisation Index

An urbanisation index based on sandy beach characteristics was employed to assess the human pressure on the sandy beaches of Anmyeon Island [33]. This index integrates seven key indicators to quantify the anthropogenic impact on sandy beaches (Table S1a,b): proximity to urban centres, buildings on the sand, cleaning of the beach, solid waste in the sand, vehicle traffic on the sand, quality of the night sky, and frequency of visitors [32]. To ensure the objectivity of this evaluation, the scoring of each indicator was based on the average values of three observers.

When the urbanisation index (UI) is less than 0.35, the beach experiences low urbanisation pressure; when the UI ranges between 0.35 and 0.6, the urbanisation pressure level is deemed moderate; and when UI exceeds 0.6, the beach is regarded as facing high urbanisation pressure [22].

2.4. Dominant Index

The dominance index (Y) was used to identify dominant species within the sampling area, defined as a value of 0.02 or higher. This index is calculated using the formula

Y = ni/N × fi

where N represents the total number of individuals across all species, ni is the number of individuals of the ith species, and fi denotes the ith species’ frequency [37].

2.5. Benthic Indices

We selected seven benthic indices to assess the ecological quality of the sandy beaches and evaluate their response to urbanisation. The AZTI marine benthic index (AMBI) and the benthic index (BENTIX) were used to categorise the macrobenthos into different ecological groups based on their tolerance to organic enrichment [26,27]. The benthic polychaetes amphipods index (BPA) was calculated using the relative frequencies of polychaetes and amphipods [28]. The benthic pollution index (BPI) categorises macrobenthos into different functional groups, considering their type of feeding and life history [30]. The Shannon–Wiener evenness proportion (ISEP) and the abundance biomass comparison (W-value) were calculated using the abundance and biomass of the macrobenthos [30,32]. Lastly, the multivariate AZTI marine biotic index (M-AMBI) was calculated based on the AMBI, the Shannon index, and the richness of the macrobenthos [29].

We used AMBI software 6.0 (https://ambi.azti.es/, accessed on 1 January 2024) to calculate the AMBI and M-AMBI. The reference conditions for M-AMBI and the functional group assignments for BPI were based on those used in previous studies [38,39]. The BENTIX was calculated based on the AMBI ecological groups (Table S2).

To facilitate the evaluation and comparison of the ecological quality assessments of the seven benthic indices of sandy beaches, we classified their ecological quality into ‘acceptable’ and ‘unacceptable’ categories based on benthic values [40,41]. The formulae for calculating the benthic indices and the threshold values for categorising the ecological quality as ‘acceptable’ or ‘unacceptable’ are referenced in

Table 1. The formulae for calculating the benthic indices and the threshold values for categorising the beaches’ ecological quality as ‘acceptable’ or ‘unacceptable’.

Index	Calculation	Ecological Quality as ‘Acceptable’ or ‘Unacceptable’	References
AZTI Marine Benthic Index (AMBI)	$\begin{array}{l}=\left[\right(0\times \%\mathrm{E}\mathrm{G}\mathrm{I})+(1.5\times \%\mathrm{E}\mathrm{G}\mathrm{I}\mathrm{I})+(3\times \%\mathrm{E}\mathrm{G}\mathrm{I}\mathrm{I}\mathrm{I})\\ +(4.5\times \%\mathrm{E}\mathrm{G}\mathrm{I}\mathrm{V})(6\times \%\mathrm{E}\mathrm{G}\mathrm{V})]/100\end{array}$	Acceptable ≤ 3.3; Unacceptable > 3.3	[26]
Benthic Index (BENTIX)	$=[6\times \mathrm{\%} \mathrm{GI} +2(\mathrm{\%}\mathrm{G}\mathrm{I}\mathrm{I}+\mathrm{\%} \mathrm{GIII}\left)\right]/100$	Acceptable > 3.5; Unacceptable ≤ 3.5	[27]
Benthic Polychaetes Amphipods Index (BPA)	$=\log \left[\right(fP)/(fA+1)+1)]$	Acceptable ≤ 0.211; Unacceptable > 0.211	[28]
Benthic Pollution Index (BPI)	$\begin{array}{l}=[1-(\mathrm{a}\times \mathrm{N}1+\mathrm{b}\times \mathrm{N}2+\mathrm{c}\times \mathrm{N}3+\mathrm{d}\times \mathrm{N}4)/(\mathrm{N}1+\\ \mathrm{N}2+\mathrm{N}3+\mathrm{N}4)/\mathrm{d}]\times 100\end{array}$	Acceptable ≥ 40; Unacceptable < 40	[38]
Shannon–Wiener Evenness Proportion (ISEP)	$={\mathrm{l}\mathrm{o}\mathrm{g}}_{10}(1/[H^{\prime }\left(biomass\right)/H\prime \left(abundance\right)]+1)$	Acceptable ≥ 0.359; Unacceptable < 0.359	[42]
Multivariate AZTI Marine Biotic Index (M-AMBI)	$=K+(\mathrm{a}\times \mathrm{A}\mathrm{M}\mathrm{B}\mathrm{I})+\left(\mathrm{b}\times H^{\prime }\right)+(\mathrm{c}\times \mathrm{S})$	Acceptable > 0.53; Unacceptable ≤ 0.53	[29]
Abundance Biomass Comparison (W-value)	$={\displaystyle \sum _{i=1}^S} \left(B_i-A_i\right)/[50(S-1)]$	Acceptable ≥ 0.15; Unacceptable < 0.15	[30]

Notes: For AMBI, EGI = disturbance-sensitive species; EGII = disturbance-indifferent species; EGIII = disturbance-tolerant species; EGIV = second-order opportunistic species; EGV = first-order opportunistic species. For BENTIX, GI = EGI + EGII; GII = EGIII + EGIV; GIII = EGV. For BPA, fP = polychaetes’ frequency; fA = amphipods’ frequency. For BPI, N1 = filter feeders or large carnivores; N2 = surface deposit feeders or small carnivores; N3 = subterranean deposit feeders; N4 = opportunistic species. For ISEP, H’ = Shannon diversity index. For M-AMBI, S = number of species. For W-value, S = number of species; Bi = biomass of species i; Ai = abundance of species i.

.

2.6. Data Analysis

To assess the seasonal differences among environmental factors, the normality of the data was first evaluated using a quantile–quantile plot. If the data conformed to a normal distribution, the independent samples t-test was applied for analysis. In cases where the data did not follow a normal distribution, the Mann–Whitney U test was used. For Spearman’s rank correlation coefficient analysis, we used data from each station to examine the correlations between benthic indices and environmental factors. Due to the impacts of seasonal changes and morphodynamics on benthic communities, a Spearman correlation analysis was conducted on the data from the same type of beach across different seasons to confirm the correlation between the benthic and urbanisation indices. In addition, kappa analysis was employed to assess the consistency of the benthic index in evaluating the ecological quality of the sandy beaches. We determined the threshold for their kappa values by referencing previous studies [43]. The above analysis was performed using IBM SPSS Statistics 29.0 (International Business Machines Corp., Armonk, NY, USA).

3. Results

3.1. Urbanisation Index and Environmental Factors

The urbanisation index of the sandy beaches ranged from 0.257 to 0.686 in June and from 0.257 to 0.714 in September. Regarding urbanisation pressure levels, seven, two, and one sandy beach experienced low, moderate, and high levels of pressure in June, respectively. In September, four, five, and one sandy beach experienced low, moderate, and high urbanisation pressure levels, respectively (

Table 2. The urbanisation index values of sandy beaches in June and September.

Urbanisation Index Value/Urbanisation Pressure Level	Baekssajang	Sambong	Giji	Anmyeon	Batkkae	Bangpo	Kkotjji	Saetppyol	Jangsampo	Baramarae
June	0.286	0.314	0.257	0.257	0.34	0.457	0.686	0.314	0.343	0.429
	Low	Low	Low	Low	Low	Moderate	High	Low	Low	Moderate
September	0.314	0.343	0.257	0.257	0.371	0.486	0.714	0.371	0.371	0.429
	Low	Low	Low	Low	Moderate	Moderate	High	Moderate	Moderate	Moderate

). Compared with June, the increase in the amount of solid waste in the sand in September led to a rise in the beaches’ urbanisation indices (Table S1a,b).

The environmental characteristics of the sandy beaches are shown in

Table 3. The environmental characteristics of the sandy beaches of Anmyeon Island in June and September.

Sandy Beach (Sampling Stations)	DO, mg/L	IL, %	Mean Grain, ∅	pH	Salinity, %	Soil Temperature, °C	Beach Type
Baekssajang (S1−S3)	5.8−8.2	0.4−0.7	2.8−3.2	7.1−7.7	10.1−27.9	19.1−25.6	Intermediate
Sambong (S4−S6)	4.2−8.3	0.2−0.9	2.6−3.2	7.2−7.5	11.9−28.7	19.1−25.8	Intermediate
Giji (S7−S9)	2.6−6.8	0.4−0.8	2.6−3.2	7.2−7.5	9.4−27.3	18.4−26.2	Intermediate
Anmyeon (S10−S12)	2−6.8	0.5−0.8	2.8−3.2	7.2−7.7	14.5−29.3	19.2−27.2	Intermediate
Batkkae (S13−S15)	2.4−5.8	0.2−1	2.8−3.1	7.2−7.6	20.9−25.9	18.6−28.9	/
Bangpo (S16−S18)	2.6−6	0.3−0.8	2.8−3.2	7.1−7.6	10.1−29.8	18.1−25.3	/
Kkotjji (S19−S21)	3.2−7.9	0.5−1.3	2.9−3.1	7.2−7.6	12.4−29.1	18.9−27	Intermediate
Saetppyol (S22−S24)	2.5−5.2	0.1−0.6	2.5−3	7.3−7.6	18.9−28.8	18−30.1	/
Jangsampo (S25−S27)	2.4−6.1	0.3−0.7	2.8−2.9	7.3−7.5	10.2−27.9	19.5−29.2	Intermediate
Baramarae (S28−S30)	1.6−4.8	0.2−0.9	2.2−2.9	7.2−7.3	15.7−28.7	18.2−27.5	/

Note: DO, dissolved oxygen; IL, ignition loss; /, lack of data.

. Based on a previous study [36], we confirmed that six beaches, excluding four beaches (i.e., Batkkae, Bangpo, Saetppyol, and Baramarae) due to lack of data, were considered intermediate beaches with a dimensionless fall velocity ranging from two to five (Table 3). In the results of the independent samples t-test, the pH values exhibited significant seasonal variation (p ≤ 0.05). The soil temperature values, in the Mann–Whitney U test, showed significant seasonal variation (p ≤ 0.05).

3.2. Macrobenthic Composition and Dominant Species

The macrobenthic species were identified for 84 taxa in the study. The most common taxa were Arthropoda (37 species, 44.05%), Annelida (22 species, 26.19%), Mollusca (20 species, 23.81%), Echinodermata (2 species, 2.38%), Hemichordata (1 species, 1.19%), Nematoda (1 species, 1.19%), and Nemertea (1 species, 1.19%) (Table S3). The average number of individuals and biomass for macrobenthos were 2725.6 ± 3530.99 ind./m² and 54.47 ± 81.77 g./m² in June, respectively. In September, the average number of individuals and biomass for macrobenthos were 1894.5 ± 1445.58 ind./m² and 45.99 ± 58.8 g./m², respectively (Table S3 and

Table 4. The dominant species and their dominance values for sandy beaches in June and September.

Season	Taxa	Species	Dominance Value
June	Amphipoda	Urothoe convexa	0.204
	Spionida	Pygospio sp.	0.092
	Scolecida	Armandia lanceolata	0.077
	Amphipoda	Haustorioides koreanus	0.053
	Amphipoda	Eohaustorius longiclactylus	0.038
September	Amphipoda	Urothoe convexa	0.190
	Spionida	Pygospio sp.	0.072
	Scolecida	Armandia lanceolata	0.065
	Amphipoda	Eohaustorius longiclactylus	0.054
	Spionida	Spiophanes sp.	0.043
	Decapoda	Scopimera globosa	0.025

). The number of species and their abundance and biomass at each station on the sandy beaches in June and September are shown in Figure 2. Five and six species were dominant in June and September, respectively. The species Urothoe convexa exhibited the highest dominance value in June and September (Table 4).

3.3. Benthic Indices

In the ecological group assignments for the AMBI, the EGI (disturbance-sensitive species) contained 35 species, the EGII (disturbance-indifferent species) contained 24 species, the EGIII (disturbance-tolerant species) comprised 12 species, the EGIV (second-order opportunistic species) comprised 3 species, and the EGV (first-order opportunistic species) included 1 species. A total of nine species were not assigned (NA). In the functional group assignment for the BPI, N1 (filter feeders or large carnivores) included 35 species, N2 (surface deposit feeders or small carnivores) included 19 species, N3 (subterranean deposit feeders) included 11 species, N4 included 2 species (opportunistic species), and NA included 16 species. Most macrobenthos were categorised into EGI and EGII (70.24%) in the AMBI ecological group. In contrast, the BPI classified them into N1 or N2 (64.29%) in the BPI functional group (Table S2).

The maximum AMBI was 4.078, the minimum was 0.013, and the average was 1.21 ± 0.88. The maximum BENTIX was 6, the minimum value was 2.22, and the average value was 4.69 ± 1.11. The maximum BPA was 0.296, the minimum was 0, and the average was 0.114 ± 0.09. The maximum BPI was 99.6, the minimum was 33.8, and the average was 69.35 ± 19.09. The maximum ISEP was 1.40, the minimum was 0.082, and the average was 0.423 ± 0.242. The maximum M-AMBI was 0.833, the minimum was 0.346, and the average was 0.582 ± 0.121. The maximum W-value was 0.441, the minimum was −0.305, and the average was 0.077 ± 0.178. Figure 3 shows the results regarding ecological quality and the values of the benthic indices for each station.

In June, 29 stations were classified as having an acceptable ecological quality per the AMBI, whereas 1 station was classified as unacceptable. A total of 24 stations were classified as being of acceptable ecological quality according to the BENTIX, whereas 6 stations were classified as unacceptable. Using the BPA, 24 stations were classified as being of acceptable ecological quality, and 6 stations were classified as unacceptable. Using the BPI, 24 stations were classified as being of acceptable ecological quality, and 6 stations were classified as unacceptable. For the ISEP, 13 stations were classified as being of acceptable ecological quality, and 17 stations were classified as unacceptable. According to the M-AMBI, 18 stations were classified as being of acceptable ecological quality, and 12 stations were classified as unacceptable. Using the W-value, 12 and 18 stations were classified as being of acceptable and unacceptable ecological quality, respectively (Figure 3 and Figure 4).

In September, 30 stations were classified as being of acceptable ecological quality using the AMBI, while 25 stations were classified as acceptable using the BENTIX, and 5 stations were classified as unacceptable. Using the BPA, 24 stations were classified as being of acceptable ecological quality, and 6 stations were classified as unacceptable. Using the BPI, 27 stations were classified as being of acceptable ecological quality, and 3 stations were classified as unacceptable. A total of 14 stations were classified as being of acceptable ecological quality according to the ISEP, whereas 16 stations were classified as unacceptable. Per the M-AMBI results, 22 stations were classified as being of acceptable quality and 8 stations were classified as unacceptable. Using the W-value, 9 and 21 stations were classified as being of acceptable and unacceptable ecological quality, respectively (Figure 3 and Figure 4).

3.4. Correlation and Kappa Analyses

According to the Spearman’s rank correlation coefficients, the AMBI was negatively correlated with the BENTIX, BPI, and M-AMBI. The AMBI was positively correlated with the BPA. The BENTIX was negatively correlated with the BPA, and the BENTIX was negatively correlated with the BPI and M-AMBI. The BPA was negatively correlated with the BPI. The ISEP was positively correlated with the M-AMBI and W-value. The M-AMBI was positively correlated with the W-value. The AMBI was negatively correlated with DO. The BENTIX was negatively correlated with soil temperature; the ISEP and W-value were positively correlated with IL and mean grain; and the M-AMBI was positively correlated with mean grain. However, the BPA and BPI showed no correlation with any of the considered environmental factors (Figure 5). The urbanisation index was positively correlated with the ISEP and W-value in June (

Table 5. The Spearman’s rank correlation coefficients between benthic indices and the urbanisation index of intermediate beaches in June and September.

Urbanisation Index	AMBI	BENTIX	BPA	BPI	ISEP	M-AMBI	W-Value
June	−0.102	0.099	0.315	−0.111	−0.569 *	−0.159	−0.560 *
September	−0.436	0.394	−0.353	0.245	0.358	0.426	0.369

Note: *, p ≤ 0.05.

).

In the results of the kappa analysis, the BENTIX and BPI values were the highest (0.741), and their level of agreement was very good. The AMBI and ISEP values were the lowest (−0.033), and their level of agreement was null. Despite a high match percentage of 78.3% between the AMBI and BPA, the kappa analysis revealed a null level of agreement (Table S5).

4. Discussion

4.1. Urbanisation Index

In September, the urbanisation index of most sandy beaches was higher than that of June (Table S1a,b). Our observations indicated that, after the tourism period of July and August, the quantity of solid waste found on beaches in September exceeded that of June. Although the urbanisation index was initially applied to sandy beaches in South America [33], it has been used on sandy beaches in Korea in this study for the first time. This index, which utilises seven indicators to assess beach urbanisation comprehensively, accurately evaluated the degree of urbanisation of Korean sandy beaches in this study. However, the index also has its limitations. Assessing its seven indicators relies heavily on the observer’s subjective judgment, making comparisons with other studies challenging. Furthermore, we observed artisanal invertebrate fisheries on all the sandy beaches; these types of fisheries significantly impact sandy beach ecosystems [5]. We recommend that, when employing the urbanisation index, it is essential to quantify the evaluation criteria of its seven indicators (for example, the evaluation of the ‘Frequency of Visitors’ should be based on the specific number of tourists on the sandy beaches), taking into complete account the human-activity-specific local characteristics of the sandy beaches (for example, in Korea, tourists are permitted to collect the invertebrates found on beaches).

4.2. Macrobenthic Community

In our study, on the sandy beaches of Anmyeon Island, the primary taxa identified were arthropods, which is consistent with the findings from previous research [44]. While a noticeable reduction in macrobenthic density from June to September is observed, suggesting changes in the macrobenthic community (Figure 2), the absence of an undisturbed sandy beach with the same morphodynamics for use as a control makes it challenging to conclusively attribute these changes solely to human activities. This underscores the importance of considering seasonal dynamics when assessing changes in these communities. In future research, incorporating an undisturbed sandy beach with the same morphodynamics as a control group will be crucial to distinguish between the effects of human activities and natural seasonal fluctuations on the macrobenthic community. Moreover, despite a slight decrease in biomass in September, there was a substantial reduction in the biomass of edible molluscs (like Dosinia (Phacosoma) japonica or Solen strictus).

4.3. Ecological Quality of the Sandy Beaches

Sandy beaches are often considered naturally stressed due to their specific environmental conditions and the natural stress caused by the dominant highly hydrodynamic regime [45]. Benthic indices are more likely to exhibit poor performance in naturally disturbed environments, potentially leading to overestimations or underestimations of ecological quality [46,47].

The results of the benthic index evaluations (except for the ISEP, M-AMBI, and W-value) indicate that the ecological quality of most of the sandy beaches was acceptable (80–100%) (Figure 4). The AMBI and BENTIX classify macrobenthos based on their tolerance to organic matter. However, on the sandy beaches of Anmyeon Island, the organic matter content was low, and most macrobenthos were assigned to EGI and EGII. This leads to an overestimation of the beaches’ ecological quality when using the AMBI and BENTIX.

While the BPA has been effectively utilised in various geographic regions to assess the impact of diverse human-induced pressures on the benthic environment [28], it is crucial to acknowledge that the BPA calculations consider amphipods and polychaetes. Therefore, the absence or excessive abundance of either of these groups can significantly influence the BPA value. In this study, the substantial presence of amphipods could result in overestimating ecological quality when utilising the BPA index.

The BPI, which was developed based on the ecological environment of the Yellow Sea and is widely employed in South Korea [48,49,50,51], also has its limitations. Notably, when there is a significant abundance of N1 and N2, the BPI tends to overestimate ecological quality [52]. In the study, most benthic organisms were classified as N1 or N2, indicating that BPI has led to the overestimation of ecological quality.

At station 29, in June and September, we observed the most contradictory ecological quality assessments. While the ecological quality was deemed acceptable based on five indices’ evaluations, the assessments from the ISEP and W-value indicated that its ecological quality was unacceptable. We discovered that, in June at station 29, the individual numbers of Armandia lanceolata and Urothoe convexa accounted for 51.1% of the total macrobethos, but their biomass was only 11.8%. Similarly, in September, the proportion of these two species in terms of individual numbers was 35.8%, but their biomass was only 6%. The ISEP and W-value are calculated based on the abundance and biomass of macrobenthos. When the abundance of a single species is excessively high, but its biomass is low, it can affect the accuracy of these indices in assessing ecological quality. This suggests that the ISEP and W-value underestimate the ecological quality of sandy beaches. Although the ISEP was developed based on the marine environment of the Yellow Sea, the threshold settings of the index were referenced from sublittoral zone data. Therefore, when using the ISEP to assess sandy beaches, it is necessary to recalibrate its threshold value [32]. The W-value is an index established based on contrasting larger animals with higher biomass and longer reproductive cycles (k-strategists) with smaller animals with lower biomass and shorter reproductive cycles (r-strategists) [53]. In this study, the accuracy of the W-value appears to have been affected due to the sampling being within recruitment periods. Furthermore, the effectiveness of the w-value in the intertidal zone necessitates further verification [54,55].

At station 7, in June and September, only the M-AMBI categorised the ecological quality as unacceptable. Establishing reference conditions is crucial for the computation of the M-AMBI [29]. Due to a lack of historical data, we found reference conditions for our study by adopting the approach used in others’ research, which involved selecting the highest observed richness and diversity values and increasing them by 15% [37,56]. However, these studies focused on setting reference conditions for subtidal zones or estuaries, which may have led to an inappropriate establishment of reference conditions in our study. This could result in the M-AMBI underestimating the ecological quality of sandy beaches.

4.4. Statistical Analysis

In the correlation analysis, it was observed that the ISEP and W-value exhibited correlations with the M-AMBI. In contrast, no significant correlation was found between the ISEP, W-value, and the other indices. This discrepancy can be attributed to the fact that the ISEP and W-value are computed based on the abundance and biomass of benthic organisms, which diverge from the underlying principles of the other indices. Moreover, the noteworthy correlation between the ISEP and W-value and their very good level of agreement in the kappa analysis suggests that the ISEP could serve as a viable alternative index to the W-value.

The characteristics of the sediment and the morphodynamics of the sandy beach significantly influence the structure of the benthic community [57]. The ISEP and W-value responded to the content of organic matter (IL%) and mean grain size in the correlation analysis. This indicates that the ISEP and W-value better reflect the changes in benthic communities than other indices.

In June, the ISEP and W-value correlated with the urbanisation index at intermediate beaches, but no benthic indices correlated with the urbanisation index in September. The ISEP and W-value are calculated based on the abundance and biomass of macrobenthos. This leads to a situation where the values of these indices are easily influenced by juveniles when surveys are conducted during recruitment. For example, in this study, Kkotjji Beach, with a higher urbanisation index value, saw a significant presence of Armandia lanceolate and Pygospio sp. juveniles. This likely led to the observed correlation between the ISEP and W-value and the urbanisation index. In September, with the extensive mortality of juveniles, the ISEP and W-value were not correlated with the beach’s urbanisation index. Overall, the ISEP and W-value correlations with the urbanisation index in June are likely coincidental. It is necessary for subsequent studies to determine whether the significant occurrence of juveniles of Armandia lanceolate and Pygospio sp. is associated with the urbanisation of beaches. Although there was no correlation between the benthic indices and the urbanisation index, the effectiveness of the benthic indices in assessing the sandy beach ecosystem is still being determined. As previously mentioned, the urbanisation index only reflects a portion of the impact of human activities on beaches.

4.5. Shortcomings of the Study

Due to the inherent instability, complexity, and natural pressure of sandy beaches, distinguishing between the changes caused by biological processes and those induced by human activities in these environments is challenging [58]. However, previous studies have demonstrated that macrobenthic species are effective indicators for assessing the impact of human activities on beach environments. For example, in Brazil, clams can serve as indicators to evaluate the urbanisation of beaches [59]; in Ghana, ghost crabs are effective indicators for assessing the impact of beach sand mining [60]; and in Italy, Talitrus saltator can serve as an indicator of the effects of human trampling on beaches [61]. However, these studies necessitate long-term investigations and are difficult to generalise due to their broad focus on specific species.

Benthic indices, which utilise macrobenthic communities, have been valuable in evaluating the effects of human activities on sandy beaches. However, their effectiveness has sometimes fallen short of expectations. Recent research highlights the significant role that the morphodynamics of sandy beaches play in determining the success of these benthic indices [22]. In our study, although we have considered the impact of beach morphodynamics on macrobenthic communities, none of the seven benthic indices responded to the urbanisation index at intermediate beaches. We speculate that benthic communities may possess a certain resilience, allowing them to maintain stability in response to human activities. To confirm this speculation, establishing control groups (undisturbed sandy beaches with similar morphodynamics) in future research is necessary, as this would enable us to precisely assess the specific impacts of human activities on sandy beaches.

5. Conclusions

In our study, the AMBI, BENTIX, BPA, and BPI overestimated the ecological quality of the sandy beaches of Anmyeon Island, Korea. At the same time, the ISEP, M-AMBI, and W-value underestimated the ecological quality. Although we considered the morphodynamics of these beaches, none of the seven indices responded to the urbanisation index. It will be necessary to include control groups in subsequent studies to further verify the validity of benthic indices on Korean beaches. In addition, our research provides vital information for the conservation and management of sandy beach ecosystems on Anmyeon Island, and it also serves as a reference for subsequent sandy beach ecological quality studies in Korea.