Recovery stages over long-term monitoring of the intertidal vegetation in the ‘Abra de Bilbao’ area and on the adjacent coast (N. Spain)

Abstract

Few studies have documented the recovery of phytobenthic marine assemblages following the reduction of long-standing sewage discharges. This paper reports on the main changes in intertidal macroalgal communities after the gradual application from 1984 to 2006 of a sewerage plan for the metropolitan area of Bilbao. Sampling sites along a pollution gradient were surveyed eight times during the course of the above period. Improvements in water quality were followed by noticeable changes in species composition and vegetation structure. Species richness significantly increased throughout the study area, while algal cover only increased at the most degraded sites. Pollution removal promoted the development of morphologically more complex species. Intertidal vegetation at the degraded sites became progressively more similar to that at the reference site. Five recovery stages discriminated by different species (SIMPER routine) were characterized from ordination (MDS) analyses: (i) extremely degraded–Gelidium pusillum is the most abundant species which is accompanied by Bachelotia antillarum at the low intertidal level (0.75 m); (ii) heavily degraded–Gelidium pusillum remains dominant and accompanied by Caulacanthus ustulatus at the high intertidal level (1.4 m); (iii) moderately degraded–Corallina elongata becomes dominant, C. ustulatus remains abundant at the high level; (iv) slightly degraded–C. elongata remains dominant in both tidal levels, Chondracanthus acicularis and Lithophyllum incrustans are abundant at the high level, whereas the latter, Pterosiphonia complanata and Stypocaulon scoparium become abundant at the low level; (v) reference stage–Lithophyllum incrustans and Laurencia obtusa are abundant together with C. elongata at the high level, whereas Stypocaulon scoparium dominates the low level, with Bifurcaria bifurcata, Jania rubens and Cystoseira tamariscifolia as abundant species. Thus, this study reveals that phytobenthic communities are useful indicators of water quality and provide real data that contribute to the assessment of the ecological status of rocky open shores on the Basque coast.

Keywords:

• benthic monitoring
• community parameters
• ecological quality assessment
• long-term changes
• macroalgae
• marine bioindicators
• pollution abatement
• seaweed

Introduction

Pollution results in a number of environmental stresses that can profoundly influence the distribution of species (Fairweather, 1990). Domestic sewage, industrial effluents and diffuse pollution (e.g. agriculture drainage) alter natural environmental conditions by increasing the availability of nutrients and organic matter, water turbidity, siltation rates and toxic chemical concentrations. Nutrient enrichment affects marine ecosystems by promoting a shift from dominance of slow-growing perennial algae to dominance of fast-growing ephemeral algae (Cloern, 2001), while an increase in organic matter implies that filter-feeders could become more competitive than macrophytes (Kautsky et al., 1992). Turbidity results in a generally lower light regime that limits macroalgal vegetation development in the sublittoral zone (Kautsky et al., 1986; Eriksson et al., 1998). Sediment smothering reduces the availability of oxygen, nutrients and firm substrata for settlement and recruitment (Devinny & Volse, 1978; Eriksson & Johansson, 2003) and especially limits the establishment of species with sexual reproduction (Eriksson et al., 1998). The toxicity of sewage has been documented in laboratory tests (Kevekordes, 2001), but chronic low level effects in the field are difficult to detect and may result in long-term changes which may be hard to separate from natural spatial and temporal variability (Walker & Kendrick, 1998).

In general, studies on phytobenthic communities subject to pollution have documented a reduction in species richness (Brown et al., 1990; Munda, 1993; Díez et al., 1999), a decline in large perennial algae (Eriksson et al., 2002; Thibaut et al., 2005), and an increase in turf-forming algae (Gorgula & Connell, 2004) and ephemeral species (Johansson et al., 1998; Eriksson et al., 2002). It might be expected that all these responses in community structure would reverse as pollution decreases, but the recovery processes of benthic assemblages are poorly documented (Hardy et al., 1993; Bonk et al., 1996; Gorostiaga & Díez, 1996; Archambault et al., 2001; Soltan et al., 2001).

The Nervión River, which flows into the study area, has been used for the disposal of high loads of domestic and industrial wastewater since the nineteenth century (Belzunce et al., 2004). Industrial effluents derived from heavy metal processing acid, fertilizer and chemical manufacture and coke furnaces, as well as domestic wastewater from a metropolitan area with about one million inhabitants greatly deteriorated the quality of estuarine waters (Azkona et al., 1984) and the adjacent coastal waters. All these activities led to a severe water oxygen deficiency and the virtual elimination of all heterotrophic life within the estuarine sediments (Saiz-Salinas & González-Oreja, 2000). This scenario started to change in the 1980s. Some of the most polluting industries and mines closed down and, more importantly, in 1983 the local water authority, Consorcio de Aguas Bilbao Bizkaia, started up a sewerage scheme (Plan Integral de Saneamiento) whose main objectives were the recovery of a good environmental status and the restoration of good aesthetic and sanitary conditions along the estuary. In this respect, a large proportion of the raw wastewater discharge was intercepted and treated from 1990 onwards, resulting in a progressive and significant reduction in the pollution load reaching the estuary: 51.8% in biochemical oxygen demand, 70.9% in ammonia nitrogen and 81.9% in faecal coliforms (data for 2003) (García-Barcina et al., 2006).

After the sewerage scheme started, the local water authority implemented several monitoring programmes to assess the extent of the biological recovery of the estuarine and marine coastal environments. Long-term monitoring programmes have proved to be useful, efficient tools for detecting the effects of pollution and assessing recovery processes (Jan et al., 1994). The ability of benthic communities to reveal spatial and temporal changes in environmental conditions has made them a target of the EU Water Framework Directive (WFD) 2000/60/EC. The monitoring programmes carried out in the Nervión estuary and in the study area have already revealed changes in benthic communities and conceptual recovery models have been proposed (Gorostiaga & Díez, 1996; Pagola-Carte & Saiz-Salinas, 2001; Gorostiaga et al., 2004; Borja et al., 2006)

The aims of this paper are (i) to study the structural changes in phytobenthic intertidal communities over the last two decades; (ii) to define different stages in the recovery process; and (iii) to provide real data to help assess the ecological status of the semi-exposed and exposed rocky shores of the Basque coast using macroalgae.

Materials and methods

Study area

The study area extends over approximately 20 km in the north of the Iberian Peninsula (from 43°21′N; 3°7′W to 43°24′N; 2°58′W), on the south-eastern corner of the Bay of Biscay (Fig. 1). Half of the shoreline corresponds to the Abra de Bilbao, a semi-enclosed bay with a river discharge in its innermost part. The river plume is usually carried eastwards by currents influenced by the prevailing north-westerly winds. The transport process dilutes polluted waters and establishes a pollution gradient. It is mostly rocky and exposure to wave action ranges from semi-exposed places to exposed ones. Tides are semidiurnal with a spring tidal range of 4.5 m. Mean water surface temperature off the Basque coast ranges between 12°C in February and 22°C in August (Valencia et al., 2004).

Fig. 1. Location of the study area and sites.

Sample collection

Starting in 1984, eight sampling surveys were carried out (in 1984, 1992, 1996, 1998, 2000, 2002, 2004 and 2006) at four sites in summer. Ordered from the most polluted to the least polluted, the sites were (Fig. 1): Arrigunaga (Site 1), Azkorri (Site 2), Meñakoz (Site 3) and Matxilando (Site 4). Since 1996 a control site (Kobaron: Site 5) not subjected to the plume of pollution was established. In the low intertidal zone two sampling levels (0.75 and 1.4 m, respectively) were established at each site. Replicates were taken following a systematic sampling design, using 50 × 50 cm quadrats. In all cases, the replicates (n = 5) were on comparable flat or slightly sloped hard substratum surfaces. A non-destructive sampling strategy was implemented which consisted of visually assessed estimates of algal cover in % at specific level following the scale: + (<1%), 1 (1–5%), 2 (5–25%), 3 (25–50%), 4 (50–75%) and 5 (75–100%) (Braun-Blanquet, 1951). Mean cover of species was calculated for each replicate using the categorical mean of each range.

Data analysis

Mean algal cover and diversity were studied for each intertidal level. Species richness (S) and Shannon diversity (H′, log_e) were calculated by means of applying the DIVERSE routine of the PRIMER software package (Clarke & Gorley, 2006). For those indices based on ‘N’ (number of individuals), this value was replaced in the formulae by the cover of the species (Magurran, 1988). Differences in algal cover and diversity measures between sites and sampling surveys were expected. Two-way analyses of variance (ANOVA) were applied to test hypotheses on spatial and temporal patterns using Statview 4.5 (Abacus Concepts Inc., Berkeley, USA, 1996). Prior to the analyses, the homogeneity of variances was checked by F-tests. Data were transformed appropriately when necessary to stabilize variances. In several comparisons variance remained slightly heterogeneous even after transformation and analysis was run at a = 0.001 for significance test. Student Newman-Keuls (SNK) tests were applied for a posteriori comparisons of means.

For each intertidal level, temporal and spatial patterns in phytobenthic assemblages were explored by applying multivariate techniques (Clarke & Warwick, 2001). Non-parametric multidimensional scaling (MDS) ordinations were done using PRIMER software package on the basis of Bray–Curtis dissimilarity matrices calculated from square-root-transformed data. The Bray–Curtis coefficient (Bray & Curtis, 1957) is very widely employed by ecologists and environmental scientists. One of its properties is that the exclusion or inclusion of taxa not present in either sample does not affect the resemblance between two samples (Clarke et al., 2006). As there were too many observations to visualize the multivariate patterns, 38 centroids (Site x time) were examined. To calculate centroids in the multivariate space defined by the Bray–Curtis measures, principal coordinates were calculated from the Bray–Curtis dissimilarity matrix among all pairs of observations using the computer program PCO.exe (Anderson, 2003). Centroids, as arithmetic averages, were calculated using these principal coordinates. The Euclidean distance between each pair of centroids was then calculated and used as the input distance matrix for the nMDS (Terlizzi et al., 2005). The indicator species that contribute most to the multivariate patterns detected were determined by means of the SIMPER routine (Clarke & Warwick, 2001). This routine computes average dissimilarities between all pairs of inter-group samples, followed by computations of the contributions of each taxon, giving standard deviations for the values.

Likewise, temporal and spatial patterns of algae distribution according to their morphological complexity were explored by ordination (MDS) analyses. For each replicate, intertidal level and sampling site species cover data were aggregated into 11 algal functional groups following a modification of the scheme proposed by Steneck & Dethier (1994). Erected species were classified into nine groups: (i) microalgae (mainly cyanophytes), (ii) filamentous algae (uniseriate), (iii) polysiphonated algae, (iv) foliose non-corticated algae (one or two layers), (v) foliose corticated algae, (vi) terete slightly corticated macrophytes (cortex with two layers), (vii) terete corticated macrophytes, (viii) leathery macrophytes and (ix) articulated calcareous algae. Two groups of crustose algae were considered: calcareous and non-calcareous species. Mean cover values were square root transformed and Bray–Curtis dissimilarity was applied to this data set. The SIMPER routine was applied to identify the discriminating morphological types. To measure how closely related the multivariate patterns resulting from species data were to those obtained from functional groups, the RELATE routine of the PRIMER software package was applied. This routine is based on Spearman correlation coefficients (ρ) between the corresponding elements of two triangular matrices of rank dissimilarities.

Results

Changes in algal cover

High level of the low intertidal zone: 1.4 m

The total amount of substratum covered by algae ranged from 52.7 to 176.5%, with an overall average of 118.6% (SE = 2.4, n = 190). Spatial and temporal differences in mean total percentage cover were tested by a two-way ANOVA. Since the control site was only sampled from 1996, the only period considered for the analysis was 1996–2006. Differences between sites depended on the sampling year as indicated by the significant interaction of time x site (Table 1). Temporal variability was higher in the polluted sites than in the control site (Fig. 2). According to the SNK a posteriori comparisons of means, mean algal cover was significantly higher in the control site in most of the cases (exceptions were found for Site 1 in 2004, Site 2 in 2000, Site 3 from 1998 to 2002, and Site 4 in 2006). Likewise, algal cover was higher at the end of the study in all sites, but only the most polluted site (Site 1) showed an increasing trend.

Fig. 2. Species richness S, Shannon diversity H′ and algal cover (in %) through time (from 1984 to 2006) of assemblages from 1.4 m intertidal level at each sampling site. No data were collected at Site 5 in 1984 and 1992. Bars show standard error (SE) (n = 5).

Low level of the low intertidal zone: 0.75 m

The total amount of substratum covered by algae ranged from 31.5 to 185%, with an overall average of 137.3%. Differences between sites depended on the sampling year as indicated by the significant interaction of time x site (Table 1). With respect to the control site, only the most polluted sites, Site 1 (from 1984 to 2004) and Site 2 (years 1984 and 1992), showed a significantly lower mean algal cover. Likewise, Site 1 and Site 2 exhibited the highest temporal variability and significant higher values at the end of the study (Fig. 3).

Fig. 3. Species richness S, Shannon diversity H′ and algal cover (in %) through time (from 1984 to 2006) of assemblages from 0.75 m intertidal level at each sampling site. No data were collected at Site 5 in 1984 and 1992. Bars show standard error (SE) (n = 5).

Changes in diversity

During the study, a total of 128 taxa were recorded: 1 Cyanophyta, 87 Rhodophyta, 23 Phaeophyceae and 17 Chlorophyta. The total number of species recorded per site showed a clear increase in all cases: Site 1 (8 in 1984 vs 29 in 2006), Site 2 (20 vs 44), Site 3 (25 vs 42), Site 4 (29 vs 50) and Site 5 (43 in 1996 vs 51 in 2006).

Table 1. Summary of results from the ANOVAs comparing the species richness S, Shannon diversity H′ and algal cover (in %) at different sites and times for each of the two intertidal levels studied.

	1.4 m Intertidal level				0.75 m Intertidal level
Source of variation	df	MS	F	p	df	MS	F	p
Algal cover
Sites	4	9817.93	29.77	<0.0001	4	12710.34	86.41	<0.0001
Time	5	3259.08	9.88	<0.0001	5	2087.21	14.19	<0.0001
Sites × time	20	1908.65	5.79	<0.0001	20	850.63	5.78	<0.0001
Residual	120	329.85			120	147.09
Species richness S
Sites	4	402.69	122.40		4	234.12	61.88	<0.0001
Time	5	76.66	23.30	<0.0001	5	68.99	18.24	<0.0001
Sites × time	20	23.87	7.26	<0.0001	20	32.40	8.56	<0.0001
Residual	120	3.29		<0.0001	120	3.78
Shannon diversity H′
Sites	4	0.86	22.30	<0.0001	4	0.11	2.75	0.0316
Time	5	0.36	9.32	<0.0001	5	0.18	4.39	0.0010
Sites x time	20	0.19	5.01	<0.0001	20	0.20	4.93	<0.0001
Residual	120	0.04			120	0.04

High level of the low intertidal zone: 1.4 m

Differences in species richness between sites was the main source of variation (F = 122.40, p < 0.0001) but it depended on the sampling year as indicated by the significant interaction of time x site (Table 1). There was greater variability through time at the polluted sites than at the reference site (Fig. 2). The number of species remained smaller at the polluted sites than at the control site but increased with time. With respect to Shannon diversity, differences between sites varied with time (Table 1). There were no differences between the control site and Sites 3 and 4 in most of the cases (exceptions were found for Site 3 in 1984 and 1992, and for Site 4 in 1992 and 2002). On the other hand, Sites 1, 2 and 3 exhibited a significant increase in diversity with respect to the first sampling survey.

Low level of the low intertidal zone: 0.75 m

For species richness, significant interactions of time x site were also identified at this level (Table 1). An increase in the number of species was detected in all polluted sites (Fig. 3). At the end of the study (2004–2006), Sites 3 and 4 even reached values comparable to or higher than those of the control site. Shannon diversity results showed that the effect of the interaction of time x site was higher than the differences between sites and between sampling surveys (Table 1). According to the SNK a posteriori comparisons of means, only Site 1 in 1984 and 1992, Site 2 in 1992, Site 3 in 1984 and 1992, and Site 4 in 1992 exhibited significant lower values than those of the control site. Likewise, only Site 1 showed a significant increase in diversity with respect to the first sampling survey.

Spatial and temporal trends in community composition

High level of the low intertidal zone: 1.4 m

The resulting diagram from the MDS ordination (Fig. 4A) shows the spatial and temporal relationships between all sites from 1984 to 2006. Each site is represented by a symbol. From right to left a pollution gradient is detected. The most degraded sites (Sites 1 and 2) are on the right, whereas Sites 3 and 4 are located between Sites 1 and 2 and the control site (Site 5). The displacement of each site with respect to its initial position reflects the changes in the community composition over time (from 1984 to 2006). A net displacement of all the sites towards the control site becomes evident. Sites 1 and 2, which were heavily polluted at the beginning of the study, exhibit the greatest displacements. Points from the control site display the lowest dispersion (i.e. the least changes). It should be noted that since 2002, assemblages from Sites 3 and 4 are progressing in a different direction in relation to the control site. Five levels of degradation (extreme, heavy, moderate, slight and control) were differentiated. Discrimination between moderately and slightly degraded communities was subjective, since the MDS plot shows a continuum. Species that contribute most to the discrimination of these groups were identified by means of the SIMPER routine (Table 2). The degradation stages are characterized as follows: (i) extremely degraded–characterized by Gelidium pusillum, (ii) heavily degraded–Caulacanthus ustulatus becomes abundant together with G. pusillum, (iii) moderately degraded–Corallina elongata is the most abundant and C. ustulatus remains with high cover, (iv) slightly degraded–C. elongata is dominant, (v) the reference stage–Lithophyllum incrustans and Laurencia obtusa are abundant together with C. elongata.

Fig. 4. Non-metric multidimensional scaling ordinations plots based on species abundances (percentage cover) showing the separation of assemblages according to sites and time of sampling. Lines reflect the displacement of each site with respect to its initial position (from 1984 to 2006). (A) Assemblages from 1.4 m intertidal level. (B) Assemblages from 0.75 m intertidal level.

Table 2. Species percentage contributions to the similarities within each of the five degradation levels identified from the MDS ordination pattern (1.4 m intertidal level) according to the routine SIMPER.

Degradation stage	Extreme	Heavy	Moderate	Slight	Control
Gelidium pusillum	86.34	30.22	5.77
Caulacanthus ustulatus	5.97	15.52	22.82	6.85	0.98
Ulva rigida		9.1	12.26	11.69	2.76
Ceramium ciliatum		6.9	2.45	0.95
Chondria coerulescens		5.29		2.13
Corallina elongata		5.22	42.82	45.73	33.57
Cladophora lehmanniana		5.04		1.47	3.13
Gastroclonium reflexum		4.69		1.35
Gelidium pulchellum		4.58	1.4	3.32
Spirulina sp.		2.13
Pterosiphonia pennata		2.05
Ceramium flaccidum			1.87
Ralfsia verrucosa			1.49	4.16	2.14
Lithophyllum incrustans				4.49	12.91
Chondracanthus acicularis				4.3
Falkenbergia rufolanosa				1.65	5.11
Codium decorticatum				1.24
Petrocelis cruenta				0.84
Laurencia obtusa					13.67
Stypocaulon scoparium					4.41
Boergeseniella thuyoides					2.94
Colpomenia peregrina					2.69
Jania rubens					1.95
Cladostephus spongiosus					1.67
Ophidocladus simpliciusculus					1.65
Sphacelaria cirrosa					1.02
Note: The 26 listed taxa have a cumulative contribution over 90% of the total similarity of each stage. The most contributing taxa (contribution >20%) are in bold type.

Low level of the low intertidal zone: 0.75 m

A pollution gradient from right to left is detected in the MDS ordination plot (Fig. 4B). The most degraded site (Site 1) is on the right and the control site (Site 5) is on the left. The net displacement of the polluted sites is towards the control site. Communities from Sites 1, 2 and 3 exhibit great changes. By contrast, points from Site 4 display a low dispersion, similar to that of the control site. Five levels of degradation were differentiated from this analysis and the species that mostly contributed to the discrimination of these groups were identified using the SIMPER routine (Table 3): (i) extremely degraded: characterized by Gelidium pusillum and Bachelotia antillarum, (ii) heavily degraded–dominated by Gelidium pusillum, (iii) moderately degraded–dominated by the calcareous Corallina elongata, (iv) slightly degraded–Corallina elongata remains dominant and other species such as Lithophyllum incrustans, Pterosiphonia complanata and Stypocaulon scoparium become abundant, and (v) the reference stage–characterized by the dominance of Stypocaulon scoparium and high abundances of Bifurcaria bifurcata, Jania rubens and Cystoseira tamariscifolia.

Table 3. Species percentage contributions to the similarities within each of the five degradation levels identified from the MDS ordination pattern (0.75 m intertidal level) according to the routine SIMPER.

Degradation stage	Extreme	Heavy	Moderate	Slight	Control
Bachelotia antillarum	34.06
Polysiphonia atlantica	6.77
Ulva rigida	5.58	6.32	6.24	4.22	2.30
Gelidium pusillum	43.61	32.96
Ceramium ciliatum		7.99	5.00
Chondria coerulescens		8.61	7.75	5.44
Gastroclonium reflexum		6.36
Corallina elongata		6.00	41.36	29.73	6.52
Ceramium flaccidum		6.69	1.87
Ceramium botryocarpum		4.71	1.06
Cladophora lehmanniana		3.99	1.94	1.16	2.23
Gelidium attenuatum		4.24	6.59
Caulacanthus ustulatus		3.74	3.15
Mesophyllum lichenoides			6.39
Lithophyllum incrustans			3.70	11.16	3.75
Falkenbergia rufolanosa			2.40	6.09	2.35
Codium decorticatum			1.89	4.82
Chondracanthus acicularis			1.04
Acrosorium venulosum				2.75
Pterosiphonia pennata				1.62
Stypocaulon scoparium				8.86	28.22
Pterosiphonia complanata				10.08	1.52
Jania rubens				1.53	8.61
Plocamium cartilagineum				2.24	1.47
Halurus equisetifolius				1.14	2.93
Bifurcaria bifurcata					15.07
Cystoseira tamariscifolia					7.69
Ophidocladus simpliciusculus					2.83
Pterosiphonia ardreana					3.66
Cladostephus spongiosus					2.50
Note: The 30 listed taxa have a cumulative contribution over 90% of the total similarity of each stage. The most contributing taxa (contribution >20%) are in bold type.

Functional group distribution

Figure 5 shows the diagrams resulting from MDS ordinations based on the functional groups data sets. The RELATE routine indicates that the spatial–temporal patterns based on these functional groups are closely related (1.4 m intertidal level: ρ = 0.83, p = 0.1%; 0.75 m intertidal level: ρ = 0.86, p = 0.1%) to those obtained from the full species data sets. The broad underlying trend is similar in both cases, thus functional group abundances are different in heavily polluted, moderately polluted and pristine conditions. However, a more detailed scale such as the four degradation levels identified from the full species MDS ordination plots (Fig. 4) cannot be differentiated. This point is more noticeable in the MDS of the high level (Fig. 5A) where the displacements of Sites 3 and 4 with respect to their initial positions in time are not towards the control site.

Fig. 5. Functional group MDS ordination plots on the basis of the Bray–Curtis dissimilarity measure of centroids of each site at each time of sampling. Lines reflect the displacement of each site with respect to its initial position (from 1984 to 2006). (A) Assemblages from 1.4 m intertidal level. (B) Assemblages from 0.75 m intertidal level.

Discussion

For decades the study area received huge amounts of non-treated wastewaters, which caused a dramatic deterioration of the benthic communities. Since the early 1980s, the local water authority started rehabilitation work to remedy this ecologically unsustainable situation. A large proportion of the raw wastewater discharges began to be treated in 1990, leading to an improvement in water quality. The overall conclusion of this study is that following pollution abatement there was a partial recovery of intertidal phytobenthic assemblages. Species richness significantly increased throughout the study area, though algal cover only increased at the most degraded sites. Pollution removal promoted the development of morphologically more complex species. Intertidal vegetation at the degraded sites has become progressively more similar to that of the reference site. Based on multivariate patterns we have characterized five recovery stages.

Algal cover

The results of our study show that only the most degraded assemblages (Sites 1 and 2 at the beginning of the study) experienced a significant increase in algal cover, revealing that this structural community parameter is not relevant in distinguishing between moderately degraded and unaltered vegetation. In this context, several papers (Belsher, 1974; Littler & Murray, 1975; Brown et al., 1990; Cormaci & Furnari, 1991) report less algal cover in the vicinity of sewage effluent outfalls, whereas Soltan et al. (2001) found no significant differences in algal cover in the vicinity of the Marseille sewage outfall 8 years after the setting up of a wastewater treatment plant. These results suggest that a significant reduction in algal cover takes place when a threshold of pollution intensity is exceeded. Likewise, the degree of water motion, depth (light penetration), salinity and the nature of the pollution discharged seem to play major roles in algal cover response. Shallow enclosed or semi-enclosed water bodies exhibit worldwide an increased algal growth due to eutrophication (Schramm & Nienhuis, 1996); in these environments, pollution abatement could give rise to a reduction in algal cover (Schernewski & Neumann, 2002). However, eutrophication results in a generally lower light regime and in increased organic sedimentation, factors that seem to entail algal cover decline below the first few metres in depth (Eriksson et al., 2002). Similarly, Krause-Jensen et al. (2007) found that total macroalgal cover responds to differences in water quality and salinity between Danish coastal areas, increasing in clear waters and with salinity. Thus, to correlate changes in water quality with algal cover it would be necessary to develop a model that considers the intensity and nature of the pollution involved, along with other factors causing variations such as wave exposure, sedimentation rates, salinity and light penetration.

Diversity

Species richness dramatically decreases under heavy pollution. It also drops in moderately degraded communities, although this parameter cannot discriminate here between slightly degraded assemblages and unaltered ones. The pollution abatement that took place in the study area between 1984 and 2006 was accompanied by a significant increase in species richness. However, Shannon diversity only showed a significant increase at the most degraded site. The loss of both species richness and diversity in the vicinity of sewage outfalls has been reported on many occasions (Borowitzka, 1972; Littler & Murray, 1975; May, 1985; Tewari & Joshi, 1988; Munda, 1993; Arévalo et al., 2007). Similarly, several studies (Hardy et al., 1993; Bonk et al., 1996; Gorostiaga & Díez, 1996; Archambault et al., 2001; Soltan et al., 2001) point out that the number of species in marine assemblages increases following water quality improvements.

Recovery stages

Five levels of degradation in the macroalgal assemblages from both intertidal levels were characterized on the basis of ordination (MDS) analyses, which were also identified as different recovery stages. This conclusion is based on the variation over time in communities, since vegetation at the most degraded sites became progressively more similar to that of the reference site. Nevertheless, some exceptions are found in the MDS plot of the 1.4 m tidal level. Assemblages from Sites 3 and 4 are progressing since 2002 in a different direction in relation to the control site. This result does not imply that water quality at Sites 3 and 4 became deteriorated during the period 2002–2006, since assemblages from their lowest level experienced a recovery. Rather, it reveals the weakness of the sampling design that has an unreplicated control site. It cannot be expected for all assemblages of the polluted sites to be similar to those of a single control site if pollution stops.

With this limitation in mind, the recovery stages are characterized as follows:

(i) extremely degraded–Gelidium pusillum is the most abundant species which is accompanied by Bachelotia antillarum at the low intertidal level (0.75 m); (ii) heavily degraded–Gelidium pusillum remains dominant and accompanied by Caulacanthus ustulatus at the high intertidal level (1.4 m); (iii) moderately degraded–Corallina elongata becomes dominant, C. ustulatus remains abundant at the high level; (iv) slightly degraded–C. elongata remains dominant in both tidal levels, Chondracanthus acicularis and Lithophyllum incrustans are the most abundant accompanying species at the high level, whereas the latter, Pterosiphonia complanata and Stypocaulon scoparium are abundant at the low level; (v) reference stage–Lithophyllum incrustans and Laurencia obtusa are abundant together with C. elongata at the high level, whereas Stypocaulon scoparium dominates the low level, with Bifurcaria bifurcata, Jania rubens and Cystoseira tamariscifolia as abundant species.

Filamentous or sheet-like chlorophytes (Fairweather, 1990; Bellgrove et al., 1997; Archambault et al., 2001) have been reported in degraded intertidal habitats close to domestic outfalls. These species are considered as ruderals (opportunistic) that occupy highly disturbed areas (Arévalo et al., 2007). Filamentous algae such as Bachelotia antillarum and the slightly corticated Gelidium pusillum are dominant in the most degraded stage of vegetation in this study. As pollution decreases Caulacanthus ustulatus and Corallina elongata come into the vegetation. The identification of G. pusillum and C. elongata as bioindicator species at two different pollution levels coincides with the findings of several other studies. Although calcareous red algae are considered as domestic pollution-tolerant species (Bellan & Bellan-Santini, 1972; North et al., 1972; Kindig & Littler, 1980; Arévalo et al., 2007), several authors highlight their partial replacement by Gelidium pusillum in the surroundings of the outfalls (Littler & Murray, 1975; May, 1985; Brown et al., 1990). It seems likely that these differences in Corallina spp. and Gelidium pusillum abundances are related to the intensity and nature of pollution, with the latter being more tolerant to anthropogenic disturbance. Finally, our results suggest that, in the study area, the first sign of degradation of natural communities is the loss of large perennial macrophytes (Cystoseira tamariscifolia, Stypocaulon scoparium, Bifurcaria bifurcata). The decline of large perennial algae under the effect of pollution has been extensively reported, e.g. for several species of kelps (Burrows & Pybus, 1971; Borowitzka, 1972; Littler & Murray, 1975), Cystoseira (Golubic, 1970; Munda, 1980, 1993; Sfriso, 1987; Rodríguez-Prieto & Polo, 1996; Cormaci & Furnari, 1999; Benedetti-Cecchi et al., 2001; Thibaut et al., 2005; Pinedo et al., 2007), Fucus (Russell, 1974; Vogt & Schramm, 1991; Bonk et al., 1992; Kautsky, H. et al., 1992; Eriksson et al., 2002), Sargassum (Hirose, 1978; Munda, 1993; Thibaut et al., 2005), Ascophyllum (Rueness, 1973; Bonk & Lein, 1978), Hormosira (Brown et al., 1990) or Gelidium sesquipedale (Gorostiaga & Díez, 1996; Díez et al., 2003). These species are good competitors in habitats with low stress and disturbance, but are progressively replaced by stress-tolerant species when pollutants are introduced into the environment.

Despite the noticeable changes in communities in the study area, complete restoration has not been achieved at any site. Although some authors (Smith et al., 1981) report a rapid response of benthic assemblages following the cessation of sewage discharge, in our opinion recovery processes of benthic communities seem to be much longer than their degradation processes. Once turf-forming algae are established, they appear to inhibit invasion of canopy species (Airoldi, 1998). Likewise, for those species with large zygotes and low dispersion (e.g. Cystoseira spp.), the restoration of their populations would not be possible if the remaining natural populations were far away (Bellgrove et al., 1997; Thibaut et al., 2005). In this study, the main difference between the reference vegetation and the slightly degraded communities was the absence of Bifurcaria bifurcata and Cystoseira tamariscifolia in the latter. However, whereas B. bifurcata was found in the surroundings of Sites 3 and 4, C. tamariscifolia did not appear after intense visual surveys, in spite of being a very frequent macrophyte in unpolluted stretches of the Basque coast. In the same way, Soltan et al. (2001) point out that 8 years after the setting up of a treatment plant, the Cystoseira amantacea community was not still restored, even in the site farthest from the outfall. As they suggest, the slow recovery process of Cystoseira populations could be related to the low range of egg dispersal and the distance of parent individuals. Several studies have documented the decline and even disappearance of Cystoseira species in areas influenced by pollution (Cormaci & Furnari, 1999; Thibaut et al., 2005), leading scientists to propose species of this genus as biological indicators of environment quality (Panayotidis et al., 1999; Pinedo et al., 2007).

Functional/morphological groups

Several authors (Littler & Littler, 1980; Steneck & Dethier, 1994) have proposed that functional form/group models could be used as a tool to predict changes in algal community composition resulting from disturbance. The review of the relevant literature by Padilla & Allen (2000) indicates that there are enough studies not supporting the functional-form group approach to call its generality into question. These authors also point out that there are questions such as nutrient uptake efficiency where these models would be useful. In this way, Eriksson et al. (2002) suggest that small, thin, filamentous species were favoured at the expense of large, complex algae on the Swedish Skagerrak coast as result of an eutrophication process. There is also considerable documentation (Munda, 1993; Rodríguez-Prieto & Polo, 1996; Arévalo et al., 2007; Pinedo et al., 2007) and experimental support (Benedetti-Cecchi et al., 2001; Gorgula & Connell, 2004) for the loss of large canopy-forming algae in favour of turf-forming algae as a consequence of anthropogenic disturbances. Likewise, in our work, results indicate that pollution abatement promotes the development of more complex algae. Thus, under heavily altered conditions, filamentous and slightly corticated species are the main components of vegetation. As pollution decreases, they are replaced by articulated calcareous algae, and finally, coarsely and leathery species come into the communities. Recently, several authors (Orfanidis et al., 2001; Panayotidis et al., 2004; Wells et al., 2007) have proposed the use of a functional-form group system for environmental quality assessment. However, we agree with Arévalo et al. (2007), who point out that two species belonging to exactly the same functional group can display completely different responses to pollution. They show Corallina elongata and Jania rubens as an example, which is also a valid case for our study. They propose that for environmental quality assessment purposes, methods based on abundances of indicator species seem to outperform functional-form group indices. We agree with this suggestion. Abundances of different functional groups differ from heavily polluted, moderately polluted and pristine conditions, but a more detailed scale in environmental assessment can be only achieved on the basis of species data.

From an environmental management viewpoint, this study could help to classify the ecological status (EcoQ) of the open shores on the Basque coast in relation to the European Union's WFD 2000/60/EC. Macroalgae comprise the dominant group of organisms thriving on rocky shores, and one of the key biological elements to be considered in the determination of the EcoQ of any given coastal water in the framework of the WFD (Pinedo et al., 2007). Our study reveals five environmental quality levels that under WFD nomenclature they would correspond to high, good, moderate, poor, and bad ecological status.

Acknowledgements

This work was supported by the local water supply and sewerage authority ‘Consorcio de Aguas Bilbao Bizkaia’. We thank two anonymous reviewers for helpful comments and criticism on the manuscript, and we are also grateful to our zoologist colleagues M. Bustamante and J. Tajadura.