High genetic diversity, phenotypic plasticity, and invasive potential of a recently introduced calcareous sponge, fast spreading across the Atlanto-Mediterranean basin

Sponges are considered poor invaders, and no genetic studies on introduced sponges have been performed up to now. Paraleucilla magna is the first calcareous sponge introduced to the Mediterranean and Northeastern Atlantic. The study aimed at investigating the genetic makeup and connectivity of the introduced populations of P. magna and at exploring signs of local phenotypic adaptation, to gain insight on the species invasive potential. Ten populations along the species introduction range (Brazil, Açores, Madeira, and continental Europe) were genetically characterized by using nine microsatellite markers. Most populations were genetically structured as suggested by significant Dst and Fst values, significant differences among populations (AMOVA) and the presence of private alleles. The analyzed populations belonged to three genetically homogeneous groups (K) according to the Bayesian algorithm (structure software) and the UPGMA dendrogram. Genetic diversity within populations was higher than expected. Recurrent introductions of non-randomly selected individuals from the native sources may have contributed to the heterozygote deficit found in all populations by forming pedigree structures with mating among relatives. Moreover, the species biological cycle was monitored in a population established on native Mediterranean assemblages (41°40′27″N, 2°47′25″E) and compared with the species cycle in other habitats. Contrasting life spans, growth habits, and reproduction cycles, depending on the habitat conditions, were recorded. To summarize, high genetic diversity, phenotypic local adaptation, and high reproduction rates altogether allow predicting the fast proliferation of P. magna in newly colonized regions and point to its strong invasive potential.

Although most introduced marine species remain confined to particularly favorable habitats, such as harbors, marinas or sea farms (Zibrowius 1991;Robinson et al. 2005;Glasby et al. 2007;Tyrrell and Byers 2007;Dafforn et al. 2009;Bulleri and Chapman 2010) and only a small percentage (ca. 0.1 %) modifies the structure and function of native communities (Mack et al. 2000;Colautti and MacIsaac 2004;Blackburn et al. 2011;Thomsen et al. 1 3 123 Page 2 of 16 2011), they may cause important economic losses by fouling ships' hulls, clogging intake pipes, and competing for resources with cultured species and thus lowering culture yields (Ruiz et al. 1997;Pimentel et al. 2000).
The effects are even worse when introduced species become invasive, alter the biodiversity and structure of native assemblages (Ruiz et al. 1997;Naylor et al. 2000;Grosholz 2002;Coles and Bolick 2007;Molnar et al. 2008) and modify the evolutionary rates of native species (Mooney and Cleland 2001;Grosholz 2002;Shine 2012). Indeed, marine invasions are currently considered the second most important cause of diversity loss in the world oceans (Blakeslee et al. 2010).
The genetic traits of alien populations may provide invaluable information on their capacity of proliferation and resilience in their introduction area (Holland 2000;Grosberg and Cunningham 2001;Féral 2002;Turon et al. 2003;Rius et al. 2008;Geller et al. 2010;Pineda et al. 2011). The study of the genetic structure of introduced populations allows understanding dispersal patterns, whether natural or man-produced, identifying colonization events, and predicting the species invasive potential (e.g. Rius et al. 2012).
Phylogeography studies of marine organisms have been traditionally based on mitochondrial genes, in particular, the COI barcode marker (Avise 2009). However, the several partitions of COI used up to now for phylogeography studies seem to be extraordinarily conserved in the Phylum Porifera (Duran et al. 2004a;Wörheide et al. 2005;Uriz and Turon 2012; but see León-Pech et al. 2015) and thus are poorly informative for analyses of population genetics at ecological time scales. Conversely, hypervariable molecular markers such as microsatellites have proved suitable for studying recent historical events in demography and population genetics of sponges (see Uriz and Turon 2012 for a review; Chaves-Fonnegra et al. 2015) and may help in tracking the expansion patterns of recently introduced species.
Only few sponges introduced to foreign habitats have been reported up to now (Calcinai et al. 2004;Pérez et al. 2006;Longo et al. 2007;Van Soest et al. 2007; Avila and Carballo 2009;Henkel and Janussen 2011) and none of them have been studied genetically at large geographical scales. The calcareous sponge Paraleucilla magna Klautau, Monteiro, and Borojevic 2004 is the first known alien sponge established in native Mediterranean assemblages (Longo et al. , 2007 and emerges as an exceptional model for studying the invasive capacities of introduced, sessile marine species with poor natural dispersal capacity. It was first described from Rio de Janeiro, where it was considered an alien ) and has proliferated in the Atlantic and Mediterranean during the last years, with records in the central (Zammit et al. 2009) and western Mediterranean (Longo et al. 2007;Guardiola et al. 2012; this study), Adriatic (Cvitkovic et al. 2013), South of Portugal, and Madeira and Azores Archipelagos (this study). Adult individuals recorded on native Northwestern Mediterranean assemblages, disappeared in July after larval release (Guardiola et al. 2012, authors pers. obs.). However, presence along the whole year has been reported in South Italy  and Brazil (Padua et al. 2012;Lanna et al. 2014). Populations of P. magna appear to be genetically structured across time, and at very short distances (i.e., less than one hundred meters), which was attributed to the genetic drift associated with small population sizes .
Studies on population genetics of introduced species rarely have considered phenotypic and life history traits in parallel, although this type of information may help to foresee their success as invaders. Here, we analyzed the genetic features of ten populations of P. magna along the species introduction range and monitored the species life cycle and population dynamics on a native assemblage of Northwestern Mediterranean.
We did not attempt to establish the species origin since only introduced populations are currently known. To determine the species introduction pathway is also challenging because the species was simultaneously detected in several distant locations, which suggests parallel introduction events . Moreover, if several introductions have occurred in a short period of time, the phylogeographic signal might be weakened or even lost, which will difficult tracking the introduction pathway (Pineda et al. 2011). Rather, the study aimed at investigating the genetic makeup and connectivity of the introduced populations, exploring signs of local phenotypic adaptation in parallel to gain insight on the species invasive potential. Sagres (SGR)-Portugal-37°00′31″N, 8°56′35″W; and Rio de Janeiro-Brazil-(BRZ) 22°54′10″N, 43°12′27″W (Fig. 1). These locations were chosen to cover most of the species known range and include sites with biogeographic interest. All specimens were collected from rocky, vertical and sub-horizontal native assemblages at depths that ranged between 6 and 15 m, and dwelt on fleshy seaweeds, mussels, or directly on the rock.

Population genetics
Most population samples consisted of ca. 30 individuals, although one population (La Herradura), which was sampled exhaustively, only provided 14 individuals (Table 1). A fragment of 0.5 cm 3 was removed from the apical part of every individual to ensure not to take embryo-harboring tissue since the zone close to the oscula is oocyte free . A maximum of two alleles were obtained from each sample, confirming that we succeeded at obtaining parental tissue only. The samples were preserved in absolute ethanol and stored at −20 °C until DNA extraction.

DNA extraction, amplification and microsatellite genotyping
DNA was extracted using PureLink Extraction Kit (Real Pure) following the manufacturer's indications and the protocol described in Pascual et al. (1997). The nine polymorphic microsatellites were amplified using the following primers (cal_a-, cal_b-, cal_c-, cal_d-, cal_e-, cal_f-, cal_g-, cal_h-and cal_j forward and reverse) previously designed for the species (Agell et al. 2012). Forward primers were labeled with a fluorescent dye from Applied Biosystems: PET red (cal_a, cal_h), NED yellow (cal_b, cal_g), VIC green (cal_c, cal_e, cal_j) and FAM blue (cal_d, cal_f). The length of PCR products was estimated relative to the internal size standard GeneScan 500LIZ (Applied Biosystems) and determined using GeneMapper ® and PeakScanner ® software from Applied Biosystems. Three independent readers checked the results visually to avoid scoring errors.   Excoffier et al. 2005). The probability of recent effective population size reductions from allele data frequencies was tested under three mutation models of evolution: infinite allele mutation (IAM), two-phase mutation (TPM), and stepwise mutation model (SMM), using BOT-TLENECK software. The method is based on the assumption that in non-bottlenecked populations (close to mutation drift equilibrium), the value of expected heterozygosity at Hardy-Weinberg equilibrium (He) is equal to the expected heterozygosity at mutation drift equilibrium (Heq). The excess of He over Heq is the evidence of severe reduction in population effective size that is compatible with a bottleneck or founder event (Cornuet and Luikart 1996).
The genetic differentiation among populations was determined by means of the statistics F st (Weir and Cockerman 1984) with ARLEQUIN v.3.11 (1000 permutations, Excoffier et al. 2005) and D est (Jost 2008), calculated with RStudio (2012) using the package DEMEtics v.0.8-3 (Gerlach et al. 2010). The latter software estimates the p values (FDR corrected) and the confidence intervals according to Manly (1997) with a bootstrap method that distributes the alleles for a specific locus randomly when the populations are in Hardy-Weinberg equilibrium for this locus. If the populations are not in Hardy-Weinberg equilibrium, the genotypes instead of the alleles are distributed randomly among populations (Goudet et al. 1996). The p values indicate how different the allele distribution in the studied sample is from a randomly obtained one. A hierarchical analysis of molecular variance (AMOVA) was performed with ARLEQUIN (Excoffier et al. 2005). The sources of variation considered in the AMOVA were between regions (Atlantic and Mediterranean basins), among groups (Iberian Mediterranean localities plus Sagres, NE Atlantic localities plus La Herradura, and Brazil), among populations, among individuals within populations, and within individuals. The AMOVA was run twice, either clustering populations into basins according to their geographical location or according to the nature of the water masses bathing the populations (Bryden et al. 1989) (i.e., the Atlantic Sagres was included in the South Mediterranean group and the Mediterranean La Herradura formed part of the North Atlantic group).
The number of genetically homogeneous groups was inferred using a Bayesian algorithm with STRUCTURE software v.2.3.3. (Pritchard et al. 2000;Falush et al. 2003Falush et al. , 2007Hubisz et al. 2009) with the parameters: range of K = 1 to 4, MCMC repetitions = 50,000, burnin range = 100,000, runs = 20) and the ad hoc statistic ΔK, based on the rate of change in the log probability of data between successive K values (Evanno et al. 2005). The results from STRUCTURE were merged across the 20 runs for the best K with CLUMPP v.1.1.2b (Jakobsson and Rosenberg 2007), while DISTRUCT v.1.1 (Rosenberg 2004) was used to graphically represent the results.
The genetic distance Ds (Nei's standard distance 1972) between the populations was calculated with GENDIST, a program of the package PHYLIP v. 3.69 (Felsestein 1993). The matrix generated was then used in the program NEIGHBOR from PHYLIP to obtain the Unweighted Pair Group Method with Arithmetic Mean (UPGMA) tree.
Isolation by distance was tested with a Mantel test using the correlation between pairwise genetic differentiation (D est values) and the logarithm of the geographical distances between the populations (XLSTAT v. 2012 5.02 Microsoft Excel add-in, Addinsoft's Core Software).

Population dynamics
A population of P. magna was monitored from July 2007 to July 2008. The population was established on a semi-sciophilous seaweed assemblage, dominated by Halopteris scoparia, from 8 to 12 m of depth, at the Blanes litoral (NW Mediterranean, 41°40′27″N, 2°47′25″E). Sampling periodicity was monthly, except on two occasions when the sea conditions prevented diving (i.e., no sampling was conducted in February 2008, and the sampling periodicity was 6 weeks instead of 4 from January to March). For the population dynamics' study, randomly selected, 20 × 20 cm surfaces (three replicates) were scraped to reveal the rocky substrate monthly from July 2007 to April 2008. The samples were taken to the laboratory and examined under a stereomicroscope. Calcareous sponges from 1.5 mm in size on, purportedly belonging to P. magna, were isolated and their identity confirmed by amplifying three to five species-specific microsatellites . The height of those individuals unambiguously ascribed to P. magna was measured. In May-June, when the 20 × 20 cm 2 quadrats proved inappropriate to represent adequately the size classes present in the population because of insufficient sponge density, the monitoring was conducted in situ by randomly placing three 1 m 2 quadrats on the vertical rocky wall where the population was established. Individuals of P. magna within the quadrats were counted, and their maximum height measured. Then, counts were scaled to the 20 × 20 cm 2 area used in the previous months to make data comparable.
From July to September (2007), no individual in the scrapped quadrats could be ascribed to P. magna recruits. These months without recruits are not included in the graphic representation of the population dynamics, which started in October 2007 (Fig. 6).
Measuring individual height across the whole sponge cycle as an estimator of growth was a compromise solution. Height was a good proxy of increase in size for recruits and juveniles, which mainly grew up with a cylindrical shape. Although this measure was less appropriated for measuring larger irregular individuals, it also captured some growth, as irregular individuals also increased in high. Moreover, the study focused on uncovering the sponge growth during the first months after recruitment and on assessing how long the individuals took to reach their maximum size in the population, and these two goals were acceptably achieved with the measure used.

Life cycle
The sponge life cycle was analyzed from December 2007 (when some individuals reached 16 mm in height) to July 2008. Five, ≥16 mm high, randomly selected specimens were collected monthly for detecting the presence of reproductive elements by means of histological sections and ultrastructure observations. From April 2008 (when oocytes were first detected in the sponge tissues), sampling was conducted every 15 days to capture more accurately gamete formation, embryo hatching, and larval release.
Sponge samples for light microscopy were preserved in 10 % formaldehyde in seawater for 5 h and then changed to 4 % formaldehyde in seawater. Once fixed, the samples were decalcified in 5 % EDTA for 2 h and subsequently rinsed in distilled water, dehydrated in an ethanol series (increased concentrations from 70 to 100 %, 30 min each), cleared in toluene/ethanol 1:1 (30 min) and then in pure toluene (15 min) and embedded in paraffin before being cut in 5-µm thick sections using a Autocut 2040 Reichert-Jung microtome. After deparaffining with xylene, the sections were stained with hematoxylin and Schiff reagent. Photographs were taken with a digital camera mounted on a Zeiss Axioplan II microscope.
Samples were also processed for transmission electronic microscopy (TEM), to confirm the presence and the nature of the reproductive elements detected through light microscope. Samples of ca. 2 mm 3 in size were fixed in 1 % OsO4 and 2 % glutaraldehyde in 0.45 M sodium acetate buffer (pH 6.4) with 10 % sucrose (1:3) (Leys and Reiswig 1998) for 12 h at 4 °C. Samples were repeatedly rinsed in the same buffer, dehydrated in an ethanol series and embedded in Spurr's resin. Ultrathin sections were stained with uranyl acetate and lead citrate, and observed with a TEM (JEOL 1010), implemented with a Bioscan system (Gatan) for image digitalization (Microscopy Unit of the Scientific and Technical Services of the Universitat de Barcelona).
The inbreeding coefficient values were positive and significant (p < 0.01) for all populations (Table 1). The exact tests for Hardy-Weinberg equilibrium confirmed these results showing significant deviations for each population (Table 1). Analysis with Micro-Checker indicated the presence of null alleles for all loci in all populations (Table 2a), which might be the cause of the deviation from Hardy-Weinberg equilibrium detected. The percentage of failed amplifications is reported in Table 2b. No evidence of scoring errors due to stuttering or large allele drop out resulted from the Micro-Checker analysis. No identical multilocus genotypes were found in any of the populations analyzed. Hence, asexual reproduction is rare or inexistent in the introduced populations.
Bottleneck analysis showed that from three to seven populations, depending on the mutation model considered, differed significantly from mutation drift equilibrium (Table 3) and thus underwent a recent founder effect. However, only Blanes, Cabo de Gata and Sagres showed bottleneck under the TPM model, while Estartit, Flores, and Brazil did it under the SMM model. The populations from Port Lligat and Madeira did not show bottleneck independently of the mutation model considered (Table 3).

Genetic differentiation and structure
The genetic differentiation values (F st and D est ) were significant for all pairwise comparisons (D est ) or for all the comparisons but three (F st ), after the FDR correction (Table 4).
This indicates that the genetic diversity was spatially structured in most of the introduced populations studied. The values of D est ranged from moderate to high in most cases being the highest between Brazil and the remaining populations but Port Lligat, and between Madeira and the Mediterranean populations. The F st values were low in general being the highest between Madeira and most of the Mediterranean populations.
The hierarchical AMOVA results (Table 5) revealed that most of the genetic diversity was due to variation within individuals (55.34 %) and among individuals within populations (34.52 %). However, differentiation among populations, which explained 6.51 % of the variance, was also significant (p < 0.001). No significant differences were found between the Mediterranean and Atlantic basins or among population clusters formed according to their geographical proximity (data not shown). However, between basins and among groups differences turned out to be significant and explained 1.18 and 2.5 % of the variance, respectively (p < 0.05 in both cases), when the groups were established according to the nature of the water masses bathing the coasts of the sampling sites (i.e., placing Sagres within the Southwestern Mediterranean group, and La Herradura in the North Atlantic group) ( Table 5).
The Bayesian clustering method (STRUCTURE) indicated the highest likelihood for the model with three genetically homogeneous groups of individuals (k) according to the highest increment of k between the models with two and three groups (Fig. 2a). The graphical representation of the Structure's coefficient memberships for the populations indicated that these three groups corresponded to the populations from the Mediterranean basin, the North Atlantic,   (Fig. 2b). The Atlantic population of Sagres (continental Portugal) shared more alleles with the Mediterranean populations than with the other Atlantic populations (Fig. 2b).

Genetic distances
The UPGMA clustering method based on genetic mean distances between populations, recovered three main clades corresponding to: (A) Mediterranean populations plus that of Sagres (South Portugal) (B) the Brazil population plus that of Port Lligat, and (C) North Atlantic insular populations (Flores and Madeira) (Fig. 3). The Flores and Madeira populations appear as the sister clade of the remaining populations. The Mantel test indicated a low but significant correlation between population genetic differentiation and geographical distance (Fig. S1). When the Port Lligat population was removed from the analysis, the correlation value increased as well as the level of significance (Mantel test, r = 0.44, p < 0.01).

Population dynamics and reproductive cycle
The population monitored dwelt on a seaweed assemblage dominated by Halopteris scoparia, placed on a vertical rocky wall exposed to water flow, at 5-10 m depth. Individuals settled directly on the rock or on the seaweeds, and the latter were then subjected to continuous movement. Individuals <2 cm in height mainly formed single tubes (Fig. 4a-c), while larger individuals tended to be roughly massive with short osculifer processes (Fig. 4d-h).
Early recruits, 1.5 mm ± 0.8 mm in height (mean ± SE), were first recorded in October (2007). Recruitment pursued for 3 months as recruits of the 1.5-2 mm size class were observed until December (2007) (Fig. 5). The individuals grew to higher size classes with time and some representatives of the largest size class (32.1-65 mm in height) were already recorded in March 2008. Individuals reached the maximum size 8 months after the first recruits were observed. From March to June, the frequencies of the individuals belonging to the three largest size-classed increased progressively (Fig. 4) but growth stopped in April 2008 for members of the largest size class, concomitantly with the beginning of the reproductive cycle (Fig. 5).
Oocytes (15-20 µm in size) were recorded in the sponge tissue during ca. 4.5 months (from April to July 2008), and one 15 µm in diameter spermatic cyst was observed in April (Fig. 6). Oocyte formation was synchronous within individuals but asynchronous between individuals. From May to July, embryos and larvae at different developmental stages were present in the population sharing the sponge tissue with oocytes (Fig. 6). Ready-to-release larvae were observed within the sponge tissue in June-July (Fig. 6). Sponge size was highly variable within the population, which indicates asynchronous larval release. Nonetheless, only one three-month long recruitment event was recorded (Fig. 5), suggesting one reproductive cycle per year. In July 2008, most larvae had already been released and adult individuals started to collapse.
Larval characteristics were the typical of Calcarea. Cleavage occurred inside a follicle "placental membrane" and gave rise to a stomoblastula (Figs. 6, 7). Micromeres produced nucleolate ciliated cells with the cilia inwards ( Fig. 7a, b) and continued to cleave forming a single layered epithelium (Fig. 7d, e) before the inversion of layers formed a hollow amphiblastula with the cilia outwards (Fig. 8). Larvae remained within the "placental membrane" until release (Figs. 6, 8a, b). A couple of granulose cells, purportedly of maternal origin, were observed in the larval central zone together with numerous bacteria (Figs. 6,  8b). Four symmetrically placed cells (cells in croix) were present between the larval ciliated cells (Fig. 8b). The cells in croix were packed with multi-vesicular dark bodies (Fig. 8c). These and the ciliated cells enclosed yolk reserves and abundant collagen-containing vesicles. Some of them open to the extracellular space between adjacent cells, purportedly corresponded to spicule formation centers (Fig. 8a, d).
Some relevant differences in the species biological cycle were evident between the monitored population and those from other localities previously studied (Table 6). Populations forming part of fouling assemblages subjected to eutrophic conditions reproduced continuously along the year and remain in the same place for years (Table 6), while those established on native assemblages under noneutrophic conditions were clearly seasonal and disappeared after reproduction (Table 6; Fig. 6).

Genetic features
Although several sponge species had been studied using microsatellites at geographical scales (Duran et al. 2004b;Blanquer and Uriz 2010;Dailianis et al. 2011;Chaves-Fonnegra et al. 2015), none of them belonged in the Calcarea class or was allochthonous in its known distribution range. Calcareous sponges have particularly short-living, short dispersal larvae  and thus expectations would point to genetic differentiation among geographically distant populations where only natural dispersal was acting. However, man-mediated dispersal modifies natural gene flow and the introduced populations of P. magna might be connected.
The microsatellite markers designed for P. magna (Agell et al. 2012) proved informative to establish the genetic structure and relationships of its introduced populations, which were located from tens to thousands of kilometers apart. These markers were more polymorphic (mean Na per locus = 18.5) in the populations analyzed than microsatellites previously studied for native demosponges (mean Na = 3.9, Na = 5.9, Na = 15.1, for populations of Scopalina lophyropoda, Crambe crambe, and Spongia officinalis, respectively), covering similar geographical ranges (Duran et al. 2004b;Blanquer and Uriz 2010;Dailianis et al. 2011). Heterozygote deficit and, thus, inbreeding was shown in all the populations studied. The presence of null alleles for all of the loci in all populations may have contributed to the deviation from Hardy-Weinberg equilibrium observed, but the homozygote excess may also be the result of mating among relatives and the formation of pedigree structures subsequently to recurrent introductions of related (non-randomly selected) individuals from the native sources.
The genetic diversity of these populations, with a mean number of alleles per population (all loci included) ranging from 45 to 90, was higher than generally expected for a recent introduction . Recently introduced species are classically reported to suffer from reduced within-population genetic variation due to founder effects and population bottleneck (Grosberg 1987;Carlon 1999). However, empirical studies have challenged the classical genetic scenario of invasive species by reporting higher than expected genetic variation in introduced populations, due to multiple introductions from native-range sources (Kolbe et al. 2007), which could be the case for P. magna and may explain the high genetic within-population diversity and the mutation drift equilibrium found in most populations under the TPM and SMM models.
Moreover, most populations were strongly structured as proved by significant D est and F st values, significant differences among populations (AMOVA), and the presence of private alleles, which may be caused by both the genetic drift and recurrent introductions from genetically differentiated populations within the native range. Nevertheless, the Mantel test showed significant isolation by distance, which became more significant when the Port Lligat population, genetically close to the distant population of Brazil, was removed. This isolation by distance pattern suggests that either the introductions followed a stepping stone model, which is not likely, or that the species distribution was not only man-mediated but also produced by some natural dispersal.
The AMOVA analysis showed significant genetic differences between the Mediterranean and Atlantic basins only when they were defined as a function of the water masses (either Atlantic or Mediterranean) bathing the populations: Sagres (South of Portugal) is under the influence of Mediterranean water masses entering the Atlantic trough the Gibraltar Strait, while La Herradura (Mediterranean) is bathed by large, tidally induced pulses of Atlantic water entering the Mediterranean (Bryden et al. 1989;Estrada 1996). This suggests some gene exchange between populations at both sides of the Strait of Gibraltar favored by the water circulation in the area.
The analyzed populations belonged to three genetically homogeneous groups according to the Bayesian algorithm of Structure and the UPGMA dendrogram, which corresponded to the North Atlantic (including La Herradura), Mediterranean (including Sagres), and Port Lligat Brazil populations. The Sagres population shared alleles with those from the Mediterranean, and the same occurred between La Herradura and the North Atlantic populations.

Phenotypic adaptation traits
Several phenotypic characteristics of the species vary across the colonized habitats. Shape ranges from long tubes in sheltered eutrophic habitats Longo et al. 2012;Cvitkovic et al. 2013) to massive clumps with a series of oscula on the top in exposed clean habitats (current study). Reproduction is almost continuous along the year in eutrophic areas of Brazil (Padua et al. 2012;Lanna et al. 2014) and South Italy (Pierri et al. 2010;Longo et al. 2012), while just one reproduction event per year was recorded in clean environments of Northwestern Mediterranean (current study).
Investment in reproduction also varied as a function of the habitat conditions. High temperatures seem to enhance the species reproductive effort in South Italy where reproduction extended along the whole year and the maximum gamete production occurred in August-September with  . However, food supply has also been proposed to influence the species reproduction in Brazil, where water temperature undergoes minor changes between summer and winter and P. magna reproduced in summer (South Hemisphere) concomitantly with the runoff of nutrients into the seawater produced during the raining season (Lanna et al. 2007). The species reproduces in the Northwestern Mediterranean in winter-spring (this study), when temperature range between 14 °C (March) and 21 °C (June) (Blanquer et al. 2008) and primary production is high due to river runoff (Estrada 1996). Moreover, recruitment occurred in autumn (this study) alongside with the second peak of primary production in the area (Estrada 1996), after thermocline breakage and the arrival of bottom nutrients to the photic zone. On the other hand, TEM studies showed that, in both eutrophic  and noneutrophic (this study) habitats, larvae are well equipped with reserves and collagen Figs. 6, 7), which confer them strength and might facilitate fast spicule formation at settlement.
Population span also seems to vary according to food supply in the introduction habitats. The adult populations completely disappeared after larval release in August in native assemblages of Northwestern Mediterranean (this study) when most filter-feeders experience food depletion (Coma et al. 2000), while overlapping of generations gave rise to stable populations in eutrophic habitats such as the Mar Piccolo de Taranto-South Italy- , which is heavily exploited for mussel farming (Caroppo et al. 2012).
Likewise, the species growth rates also seem to vary according to the trophic conditions of the colonized habitat. Early, 1.5-2 mm high recruits employ six months to reach the maximum size (32.1-65 cm high size class) in well-preserved assemblages of Northwestern Mediterranean (this study), while 5 mm 3 in volume (ca. 2 × 1.75 × 1.5 mm in size) specimens (approximated estimate) spend only Fig. 6 Life cycle of Paraleucilla magna at the Blanes (Northwestern Mediterranean) population. Oocytes were present during 4 months. Spermatogenesis was extremely short since our monthly sampling only allowed us to detect one spermatic cyst. Embryos and larvae were recorded from May to July. Larvae were abundant m specimens higher than 3 cm in height three months to reach the highest size class of 100 ml (ca. 10.5 cm × 3.5 cm × 3 cm) in the eutrophic Mar Piccolo de Taranto ).

General conclusions
The massive proliferation of P. magna around shellfish farms, in particular Mytilus cultures (Longo et al. 2007; current study) points to shellfish fouling as the most likely introduction pathway (Longo et al. 2007), although introduction by vessels cannot be totally discarded. Epibiont removal, previously to placing mussels into depuration tanks authors' obs.), may have caused the arrival of brooding individuals/fragments of P. magna to native assemblages, where larval release and settlement could have occurred. Moreover, given the species fragility, rafting of brooding fragments from established populations may also be produced by natural causes such as storms.
The species origin remains undetermined, as it has been identified as an alien in its currently known distribution range. The genus Paraleucilla was unknown out of the Indo-Pacific (http://www.marinespecies.org/porifera/) before P. magna was recorded in the Mediterranean ) and Atlantic , which points to its Indo-Pacific origin. High temperatures seem to favor the species reproductive effort , suggesting a warm ocean as its native source. Moreover, there are several Indo-Pacific countries with an intensive production of Mytilus and Crassostrea, which export around the globe (Chew 1990) and have been responsible for the distribution of many unwanted species (Haupt et al. 2009) including some of the few introduced sponges known up to now (Henkel and Janussen 2011).
Adaptation to disturbance, wide environmental tolerance, rapid growth, and high reproductive capacity are characteristics of foulers (Murray et al. 2012) shared by P. magna (Lanna et al. 2007;Longo et al. 2007Longo et al. , 2012Zammit et al. 2009; current study). Thus, the species might be an opportunistic fouler also in its native range. A fouler nature and the scarce taxonomical studies of fouling sponges in the Indo-Pacific might explain why the species has not been recorded in its purported native habitat yet.
Although very small populations of P. magna in native communities were unstable across time because of the genetic drift, lack of natural gene flow, and stochastic recruitment failure , the genetic traits of the introduced populations indicate that the species is not in risk of extinction, and constrain on adaptation does not seem to hinder the spread of the species in the Atlanto-Mediterranean region.
The genetic makeup and the high levels of phenotypic plasticity of P. magna altogether allow predicting expansion of its introduced populations (Brown et al. 2011) by adopting the best-suited phenotypes to local conditions. Moreover, long-living populations such as those of eutrophic habitats ) may operate as reservoir for dispersal and fast colonization of additional locations (Ridley and Ellstrand 2010).
To summarize, P. magna shows both, tolerance to contrasting local conditions and phenotypic plasticity, which are typical traits of invasive species (Matesanz et al. 2012). The species is expanding in the colonized regions, in particular under favorable trophic circumstances, where it may overgrow native seaweeds and filter-feeder invertebrates, and compete for trophic resources with shellfish cultures, lowering culture yields and hampering culture manipulation and maintenance. The future spread of P. magna in its introduced range is likely to be fueled by phenotypic adaptation to local conditions and high fitness across diverse habitats, rather than by local evolution of adapted genotypes.
Additional introductions can be expected in the near future reinforcing population stability and expansion since no adequate policies can be implemented until the introduction vector is unambiguously identified.