The estuary of Bilbao (43.35° N, 3.03° W), once considered among the most polluted estuaries in Europe, initiated in the late 1980s a water quality restoration program that allowed the re-colonization of the zooplankton community. The surrounding neritic species successfully colonized and became autochthonous but with a polyhaline - euhaline habitat preference (Uriarte & Villate 2004, 2005; Albaina et al., 2009; Aravena et al., 2009). More in detail, the inner part of the system changed from nearly no plankton inhabitants to dominance of the opportunistic Acartia clausi Giesbrecht 1889, and, within a few years, to a diverse “Acartia species complex” fitting in salinity driven species-specific niches (Albaina et al., 2009; Aravena et al., 2009). However, since 2003 the invasive species Acartia tonsa Dana 1849 (first cited in 2001) displaced the rest of congeneric species becoming predominant in inner waters of the Bilbao estuary (Aravena et al., 2009). The combination of empty niches, suitable environmental conditions, and availability of proper vectors (high vessel traffic with massive ballast water exchanges, Butrón et al., 2011) made the estuary of Bilbao a perfect recipient for Non Indigenous Species (NIS) (Paavola et al., 2005). Today, this estuary is among the main ports in the North-East Atlantic (NEA) with connections to > 800 ports worldwide. Shortly after its dominance in the Bilbao estuary, the nearby estuary of Urdaibai (43.41° N, 2.69° W), where A. bifilosa Giesbrecht 1881 had been dominating since 1980s (e.g. Uriarte & Villate 2005), was also colonized by A. tonsa (Aravena 2009). And since its first appearance in 2003, this species has showed a huge increase in number, becoming a codominant species nowadays (data not shown).

The introduction of NIS to estuaries and coastal embayment is of great concern and regular monitoring of these particularly sensitive ecosystems for NIS and its ecological impact is paramount for a timely answer and any further mitigation/adaptation policy design. In the present study, we assessed zooplankton communities in the inner part of the estuaries of Bilbao and Urdaibai during the summer of 2013. Our objectives were to 1) examine the genetic information of the known NIS A. tonsa in both estuaries and gain insights into the invasion source by comparing mitochondrial DNA (mtDNA) sequences obtained from the native and introduced ranges (e.g. Downie 2002; Cristescu 2015 for a review) and 2) discover any novel invasive species to these Basque estuaries.

The zooplankton sampling was carried out monthly, from July to October, in consecutive days during neap tides in waters with a salinity of ~30. Samples were obtained from horizontal tows (below halocline) of a 200 μm mesh plankton net (mouth diameter: 0.25 m) equipped with a flowmeter (Hydrobios). Net samples were transferred to the laboratory where subsamples were preserved either in 4 % borax buffered formalin seawater solution for taxonomic and abundance analysis or in absolute Ethanol for genetic analyses. The qualitative and quantitative analysis of zooplankton was carried out under a stereoscopic microscope (Olympus SZX12). Identification was made to species or genus level following mainly the ICES Identification Leaflets for Plankton (http://goo.gl/qfKg9V). In each sample, a minimum of 100 individuals of the most abundant taxa were counted before finishing sub-sampling; if not possible, the whole sample was identified. An Olympus IX70 microscope was used for Pseudodiaptomus marinus Sato 1913 identification confirmation. This species was identified following the keys in Grindley & Grice (1969) and Brylinski et al. (2012). Five to ten individuals per estuary (A. tonsa and P. marinus) were sorted for DNA sequencing of mitochondrial cytochrome c oxidase subunit I gene (MT-CO1).

Finally, as commercial ships’ ballast water exchange is generally considered as the main vector for NIS introduction to the marine environment, Bilbao port ship traffic data from the 2001–2013 period was obtained from the port authorities. The number of ships entering the Bilbao port and their previous port of origin were considered.

In the summer of 2013 the copepod community of two Basque estuaries (namely the estuary of Bilbao and the estuary of Urdaibai) was monitored to discover any non-indigenous species and any changes of known NIS Acartia tonsa, with a goal to gain insights into NIS invasion patterns following a DNA barcoding approach (e.g. Downie 2002). Apart from the expected presence of abundant A. tonsa in both estuaries, we surprisingly recorded Pseudodiaptomus marinus for the first time in the estuary of Bilbao as well as the entire Iberian Peninsula waters. The MT-CO1 gene was selected because it is widely used as the standard DNA barcoding gene with a high species-discriminating power and a vast amount of publicly available sequences.

Both morphological characters and barcode sequences confirmed the presence of P. marinus (Fig. 1a), a species considered to be native to the Northwestern Pacific Ocean and that had been recently cited in NEA waters (see Brylinski et al., 2012). Interestingly, although only four individuals were sequenced, three distinct haplotypes were recorded (Fig. 1b). In comparison, we found one A. tonsa haplotype within the Bilbao (n = 10) and two in Urdaibai (n = 8) (Fig. 2b). Extremely low mtDNA haplotypic diversities are common in invasive species due to a strong founder effect (e.g. Holland 2000; Chandler et al., 2008). However, the degree of reduction depends on several factors such as the number (and scale) of introduction events and the diversity of the invasion origin (Cristescu 2015). However, the lack of P. marinus MT-CO1 records in the public databases prevents any discussion. Concerning this species, the first record in European waters was in 2007 in the north Adriatic Sea (De Olazabal & Tirelli 2011) and since then P. marinus has been cited in other locations along the Mediterranean Sea (Sabia et al., 2014 and references therein). Regarding North-East Atlantic (NEA) waters, it has been only reported in the southern North Sea in 2010 and in 2011, in two independent surveys (Brylinski et al., 2012; Jha et al., 2013). Taken into account an unpublished record at the Gironde estuary (Brylinski et al., 2012), to our knowledge, present work represents the fourth citation for this NIS in NEA waters and the first one in the Iberian Peninsula waters. However, the 243 ind. m^-3 abundance for P. marinus in July 2013 in the estuary of Bilbao is unprecedented as previous NEA records were always below 4 ind. m^-3. Moreover, the presence of females carrying eggs (5 out of 37 individuals sorted in July) suggested an established P. marinus population in the estuary of Bilbao. Finally, further analysis of archived samples (2005 to 2012 period; monthly surveys in waters with a salinity of ~30) allowed us to date the first recorded presence of this species back to October 2010 (10/30/2010, 28 individuals with no ovigerous female).

On the other hand, A. tonsa, native to the American and Indo-Pacific area, is a well-known NIS in European waters since the early 1900s (Brylinski 1981), with ships’ ballast water exchange proposed as the most likely introduction vector (e.g. Paavola et al., 2005). In summer 2013, A. tonsa was consistently more abundant (~2000 ind. m^-3, considering only adults; Table 2) in the estuary of Bilbao, where it was the only species of the genus Acartia, than in the estuary of Urdaibai, where this species was at the same abundance levels of A. bifilosa. The rich amount of MT-CO1 sequence data in GenBank for A. tonsa allowed us to compare the herein obtained MT-CO1 sequences against a large amount of representative sequences from both the species’ native and invaded area (Table 1, Fig. 2). Chen & Hare (2011), covering 20 estuarine systems along the whole USA Atlantic coast, reported three sympatric genetic lineages for A. tonsa (namely X, F and S) for both a mitochondrial and a nuclear DNA loci, respectively, MT-CO1 and ITS (Internal Transcribed Spacer). These three clades are also shown in Fig. 2a. Every European record, including a total of eight different haplotypes, clustered together in one clade (“A” in Fig. 2a) that corresponded to the “X” lineage in Chen & Hare (2011). The fact that, to date, only one of the three A. tonsa native lineages has been reported in European waters deserves further discussion. While an insufficient sequencing effort in European waters compared with North American ones could represent a putative explanation (Table 1), it is also likely that an undetermined physiological pre-condition such as a broader tolerance range in lineage “X” has allowed the colonisation of European waters by this lineage. In this sense, Chen & Hare (2011) reported a higher gene flow corresponding to the “X” lineage in North American estuaries including preference for a higher and larger range of, respectively, salinities and latitudes than the other two native lineages. This could represent a higher resilience capacity and an obvious advantage when facing longer displacements such as trans-oceanic ones and, more importantly, for successful settlement in the invaded systems. Environmental stress associated with ballast water environment and/or a better adaptation to European water conditions (either for adults or resting eggs) in those individuals could explain the reported pattern. In this sense, the lower haplotype diversity (1 haplotype out of 10 individuals) in Bilbao may be related to the selective effect of pollution, whereas higher haplotype diversity (2 out of 8) to the relatively pristine condition of Urdaibai, an UNESCO Biosphere reserve since 1984, with mainly recreational use and a relatively low urbanization pressure. However without further experiments, such as common-garden ones, this questions will remain unsolved.

Interestingly, the two haplotypes recorded in the Bilbao and Urdaibai estuaries (H-1 and H-2 in Fig. 2b) had been previously reported at both the native and invaded European locations. While H-1, the dominant haplotype in Europe (including the ten specimens of the estuary of Bilbao along with six out of eight from the estuary of Urdaibai; Fig. 2b), was also the predominant haplotype in the native area’s lineage “X”, the H-2 had been previously cited at three estuaries in the Eastern USA coast (Gang Chen, personnel communication) and at an undetermined Italian location (GenBank accession number HE647798). More interestingly, the haplotype network for European A. tonsa MT-CO1 records, showing closely related haplotypes in Basque waters but a well-diversified pattern in North Sea area, suggests the latter as putative center of origin (Fig. 2b). This fits also the chronology of A. tonsa’s records in NEA waters, with earliest records centered in the North Sea area (Figure 1 in Brylinski 1981). Although non-European origins of A. tonsa invaded to Basque estuaries cannot be ruled out, both the previously cited presence of H-1 and H-2 within European waters and the haplotype network shape points to a secondary introduction from an European source as the most likely process. Moreover, this is further supported by the Bilbao port maritime traffic patterns when considering the putative source of NIS. In this sense, during the 2001–2013 period, from an average of 3400 commercial ships per year that entered the port of Bilbao, only 3 and 2 % of the ships arrived from a North American and Asian (Indian and Pacific) port of origin, respectively, whereas 53 % arrived from Atlantic European ports.

Paavola et al. (2005) have cited that low species richness environments are more prone to NIS successful establishment and this is supported by present data showing the current preponderance of NIS in the Bilbao estuary following pollution abatement. Shortly after A. tonsa, P. marinus also invaded the system and both became settled in a few years’ time (Aravena et al., 2009, Uriarte et al., 2016). Both species had been recorded earlier in nearby NEA systems and, corresponding to a broad range of tolerance to environmental stressors such as salinity and pollution, are capable of proliferating in harbours and other anthropized environments (see Brylinski 1981 and Brylinski et al., 2012 for, respectively, A. tonsa and P. marinus). While data are still scarce, A. tonsa and P. marinus previous records’ pattern along with present data suggest an ongoing colonization process in NEA waters. More surveillance of plankton environment, especially at potential NIS sources such as harbours or marinas is critical to anticipate further spreading of these species. Finally, the development of multi loci markers such as SSRs or SNPs for these NIS species is desirable as they have the potential to determine a more precise geographical origin of marine invasions (e.g. Reusch et al., 2010).