Explaining naturalization and invasiveness: new insights from historical ornamental plant catalogs

Abstract We identified plant attributes associated with naturalization and invasiveness using century‐old ornamental plant catalogs from Québec (Canada). We tested the hypothesis that naturalization is determined by fewer factors than invasiveness, as the latter also requires dispersal, which introduces additional complexity. The approach we used took into account not only plant attributes as explanatory factors, but also propagule pressure, while accounting for phylogenetic relationships among species. Museum collections were used, in combination with scientific journal databases, to assess invasiveness. Particular attention was given to species that never escaped from gardens and thus represent cases of “failed” invasions. Naturalization in cold‐temperate environments is determined by fewer factors than invasion, but only if phylogenetic links between species are taken into account, highlighting the importance of phylogenetic tools for analyzing species pools not resulting from a random selection of taxa. Hardiness is the main factor explaining naturalization in Québec. Invasion requires dispersal, as shown by three significant variables associated with the spread of diaspores in the invasiveness model (seed weight, hydrochory, number of seed dispersal modes). Plants that are not cold‐hardy are likely to disappear from the market or nature, but the disappearance phenomenon is more complex, involving also seed dispersal abilities and propagule pressure. Factors contributing to naturalization or invasiveness may differ greatly between regions. Differences are due in part to the plant traits used in the models and the methodology. However, this study, conducted in a cold‐temperate region, sheds new light on what is likely a context (climatic)‐dependant phenomenon.

negative economic or environmental impacts; see Richardson, Pyšek, & Carlton, 2011 for definitions). However, the proportions (introduced vs. naturalized vs. invasive or weedy) greatly differ between regions (Richardson & Pyšek, 2012) and species groups (Pemberton & Liu, 2009). For instance, in Australia and New Zealand, the proportion of introduced species that naturalized varied from 9% in some families to 76% in others (Diez et al., 2009). In Britain, 68% of the species sold in nurseries from 1885 to 1985 escaped from cultivation (Dehnen-Schmutz, Touza, Perrings, & Williamson, 2007). In Ireland, 48% of the exotic plants found in nature after 1970 have well-established populations, and 19% are truly invasive (Milbau & Stout, 2008). Of the 1112 exotic species introduced (accidentally or deliberately) in the continental part of the United States and classified as invasive, 36% are considered noxious weeds (Lehan, Murphy, Thorburn, & Bradley, 2013). At the other end of the spectrum, only 10% of the 887 exotic species naturalized in Québec are weeds (Lavoie, Guay, & Joerin, 2014). In Hawaii, 5% of the 7866 ornamental species cultivated between 1840 and 1999 naturalized, and <1% became weeds (Schmidt & Drake, 2011).
These statistics indicate that predicting how many and which species will naturalize and eventually become invasive or weedy is an extremely difficult and context-dependent task. Consequently, there is an urgent need to better understand the interactions between plant attributes and the processes which facilitate naturalization and invasiveness, to reduce uncertainties associated with predictions. This information will help plant biologists to develop efficient tools that can be used by environmental managers to prevent detrimental invasions.
The horticultural industry is a major player in the world plant market, with sales of about USD 109 billion in 2011 (Gyan Research and Analytics 2012). This industry is largely responsible for the introduction of exotic species in new regions or continents (Mack & Erneberg, 2002;Reichard & White, 2001). For instance, of the 671 invasive plants deliberately introduced in the continental United States, 426 (64%) were imported for ornamental purposes (Lehan et al., 2013). A large proportion of these species were introduced in the 19th century and in the first half of the 20th century (Lavoie, Saint-Louis, Guay, Groeneveld, & Villeneuve, 2012;Mack, 1991), but the emergence of new horticultural trading partners from tropical regions, the Middle East, and Eastern Europe could be responsible for a new wave of plant invasions, underscoring the need for efficient risk assessment tools (Bradley et al., 2012).
Nursery catalogs can be extremely useful for identifying the characteristics of plants likely to naturalize or to become invaders (Dehnen-Schmutz et al., 2007;Pemberton & Liu, 2009). They offer an excellent record of plants sold (although not necessarily bought by customers), and by comparing a list of catalog species with a list of naturalized species, those that escaped from gardens (successful naturalizations) can be easily distinguished from those that did not ("failed" invasions).
Old (>100 years) catalogs are especially relevant for building models explaining naturalization, since the species sold for more than a century and that are still not found in nature are unlikely to naturalize in the future, at least under the present-day climate.
The main objective of this study was to identify plant attributes associated with naturalization and invasiveness using century-old nursery catalogs. This is not the first study of this kind (although there are only a few: Dehnen-Schmutz et al., 2007;Pemberton & Liu, 2009;Skou, Pauleit, & Kollmann, 2012), and attempts to link invasiveness with plant attributes are multiplying (for reviews and debates on their relevance, see Pyšek & Richardson, 2007;van Kleunen, Weber, & Fischer, 2010;van Kleunen, Dawson, & Dostal, 2011;Thompson & Davis, 2011;Leffler, James, Monaco, & Sheley, 2014). However, we propose a new approach that takes into account not only plant attributes as explanatory factors, but also propagule pressure, while accounting for the nonindependence of the species analyzed due to their phylogenetic relationships. Museum collections were used, in combination with scientific journal databases, to assess invasiveness.
We paid a particular attention to the species that never escaped from gardens and were thus potential cases of "failed" invasions. We tested the hypothesis of Richardson and Pyšek (2012) that naturalization is determined by fewer factors than invasion, as the latter also requires dispersal, which introduces additional complexity.

| Taxon selection
This study was conducted using the ten nursery catalogs that were published in the province of Québec (Canada) in the 19th century, from 1817 to 1894, and that were still available from library archives (see  for the complete list). The list of taxa sold in each catalog was first extracted. There were significant changes in nomenclature (in Latin, English, and/or French) over the last 200 years.
Only taxa, including species, subspecies, varieties, and hybrids, for which the identification was certain, were retained. The taxonomic nomenclature was standardized using the Canadian Biodiversity Information Facility (2015) or Tropicos (Missouri Botanical Garden 2015) for taxa not listed in the former database.
Plants unable to escape from cultivation, that is, indoor taxa from tropical or equatorial regions (often listed as "greenhouse plants"), and taxa sold exclusively for human food production (fruits, vegetables) and with no ornamental value, were eliminated. The taxa were identified using various sources, such as ornamental plant guides, nursery catalogs, and agricultural or horticultural Web sites. The remaining taxa were mostly ornamental plants, but several were also probably sold for other purposes (e.g., medicinal plants).
We then identified the taxa no longer sold (in 2015) in Québec. (iii) the search engine tool of the Association québécoise des producteurs en pépinière du Québec indicating which nurseries in the province produce a particular ornamental taxa, and (iv) the updated list of all trees available in Québec nurseries (Dumont, 2014(Dumont, , 2015. Taxa not found in at least one of the different plant lists were checked by two professional horticulturists cumulating 55 years of experience for detecting other taxa that were available to customers in 2015.
Preliminary tests using only the number of specimens as dependant variable in multiple regression models had a poor performance for explaining invasiveness. This performance was likely related to the under-representation of species that spread mainly during the last 30-40 years, a period with a very low specimen collection effort in Québec . We nevertheless estimated that the number of specimens was a reliable source of data, as long as it was combined with another indicator of invasiveness, the scientific research effort. This effort, estimated using the number of published scientific papers, provides an indirect measurement of invasiveness: The more invasive the species, the more it attracts the attention of scientists, and the more papers focussing on this species are likely to be published (Lavoie et al., 2014). The scientific research effort was estimated using the Web of Science ™ database (Thomson Reuters 2013; last query: 10 December 2013) with the name of the taxa (in Latin) and the keyword "invasive" or "invasion" in the "title" or "topic" research fields, to extract the associated papers. Each paper was checked for relevance. Only studies clearly related to the taxa of interest and conducted in northeastern North America (the area covered by the flora of Gleason & Cronquist, 1991), that is, in a region roughly similar to Québec from a climatic and vegetation point of view, were retained.
Data on the number of specimens and the number of papers were first cubic-root transformed to normalize their distribution and then standardized on a 0-1 scale to give equal weight to the variables before analysis. On these two sets of variables, a k-mean clustering algorithm (iterated 100 times) was run in R software (R Development Core Team 2013) to partition the naturalized taxa group into k = 2 subgroups (invasive or noninvasive).

| Plant attributes and propagule pressure
A database of plant attributes was generated for the taxa listed in the catalogs. Retained attributes were those readily accessible from online databases or the scientific literature and available for all taxa (Table 1). For perennials, the hardiness zone variables were derived using the methodology developed by Lavoie et al. (2013), essentially based on the overlap between the geographic distribution of the taxa in the native and exotic ranges and hardiness zone maps. For annuals, hardiness zones are less relevant to horticulturists. However, several ornamental plant guides and Web sites provide information on the lowest temperature a taxon can tolerate, which was used to estimate the coldest hardiness zone; warmer zones were assumed to be tolerated by the taxon. A similar approach has been successfully used in the past to compare the invasion probability of annuals and perennials from historical catalogs (Dehnen-Schmutz et al., 2007).
Two variables were used as indicators of propagule pressure (sensu Lockwood, Cassey, & Blackburn, 2005), that is, the number of catalogs in which the taxon was listed and the number of years elapsed since its first mention in a catalog. We hypothesized that a taxon available from more nurseries and sold for a longer period of time would be more widely planted, thus producing more propagules with the potential to escape from gardens and to contribute to naturalization and/or invasiveness (Dehnen-Schmutz et al., 2007;Pemberton & Liu, 2009;Pyšek, Křivánek, & Jarošík, 2009;Skou et al., 2012).

| Statistical models and phylogenies
Three logistic regression models (Hosmer & Lemeshow, 2000) were constructed for this study. The naturalization (or not) of a taxon included in at least one of the ornamental plant catalogs was the dependent variable of the first model (the naturalization model). Whether a taxon from the catalogs became invasive (or not) was the dependent variable of the second model (the invasiveness model). That a taxon from the catalogs was neither sold nor naturalized in 2015 (or still sold and/or naturalized) was the dependent variable of the third model (the disappearance model, in that the taxon was no longer found in Québec in nurseries or in nature, albeit potentially still present in gardens).
The remaining variables (plant attributes, propagule pressure) were integrated in a first series of models as independent (or explanatory) variables. Prior to including plant attributes in the analyses, categorical data were coded into binary dummy variables, continuous data were log-transformed to normalize their distribution, and discrete and continuous data were standardized by subtracting the mean of each variable and dividing by two times its standard deviation to facilitate comparisons with the dummy variables (see Gelman, 2008 for details).
Linearly dependent variables and variables that showed high collinearity (VIF > 3) in the full models were removed. A forward stepwise model selection was then performed to construct the logistic regression models and finally select the best models based on the Akaike information criterion (AIC; Burnham & Anderson, 2002). All models were run in R software (R Development Core Team 2013).
Plants found in catalogs are by no means a random selection of species. To verify whether the models elaborated in this study were taxonomically or phylogenetically biased, logistic regressions correcting for phylogenetic correlations in the residuals of the models (Ives & Garland, 2010) were also performed (hereafter named phylogenetic logistic regressions). Phylogenies for the taxa were obtained from the online tool PhyLomatic version 3 (Webb & Donoghue, 2005), which is based on the APG III classification system (Angiosperm Phylogeny Group 2009). Node ages were calibrated with data from Wikström, Savolainen, and Chase (2001), and branch lengths were adjusted using the bLadj tool in PhyLocom (Webb, Ackerly, & Kembel, 2008). The R package phylolm (Ho & Ane, 2014) was then used to run a second series of logistic regression models with the whole phylogeny included as a covariance structure. This approach assumes that the residual variation follows a homogeneous model of evolution across the branches of the phylogenetic tree, and a violation of this assumption could lead to unacceptable type I error rates and/or reduced statistical power (Mazel et al., 2016). For each model, this assumption was tested by looking for rate shifts in the residuals of a standard logistic regression along the phylogeny using the auteur approach (Eastman, Alfaro, Joyce, Hipp, & Harmon, 2011) from the geiger package in R, as recommended by Mazel et al. (2016). Rate shifts were detected in all models, but the vast majority (naturalization), or almost all (invasiveness) if not all (disappearance), of these shifts occurred on branches sustaining either one species or one genera. Consequently, we concluded that the rate shifts did not have phylogenetic structure in the residuals and decided to perform the regressions with the unmodified phylogeny. As for the nonphylogenetic models, a forward stepwise selection procedure based on the AIC was used to select the best models. These models, with and without phylogeny, were analyzed side by side-to-see how incorporating phylogeny affected the significance of species attributes directly. McFadden's pseudo R 2 , correcting for the number of parameters included in the model (R 2 adj ), was estimated for each model.
T A B L E 1 Plant attributes that were used to explain the naturalization, invasiveness, and disappearance of taxa listed in nursery catalogs published in the 19th century in Québec (Canada).

| DISCUSSION
The models constructed with the historical nursery catalogs published in Québec show that naturalization in cold-temperate environments is determined by fewer factors than invasion. However, this conclusion was reached only when phylogenetic relationships were taken into account, highlighting the importance of phylogenetic tools for analyzing species pools not resulting from a random selection of taxa. This is especially true for plant catalogs, given the strong preference of horticulturists for certain families and genera with high ornamental value (e.g., Iris, Rosa, Primula, and Lilium).
Hardiness is the main factor explaining naturalization in Québec; plants tolerating a wider range of hardiness zones are also less likely to naturalize, but regardless of the number of zones, if a plant is not cold-hardy, its establishment and survival chances are low. In a coldtemperate region such as Québec, cold hardiness as an explanatory variable is unsurprising, but that hardiness is the only significant attribute for naturalization is especially revealing. Cold hardiness is not a plant trait by itself: It is an indicator of a combination of morphological and physiological traits allowing plants to survive cold temperatures, and especially frost (United States Department of Agriculture 2015c).
In Québec, being cold-frost resistant is necessary for naturalization and for the transition from naturalization to invasiveness, but other attributes not included in our models probably help the establishment and expansion of populations over large areas, such as a long flowering time, a large specific leaf area, and the presence of adequate pollinators (Bufford & Daehler, 2014;Gallagher, Randall, & Leishman, 2014). The importance of cold tolerance is highlighted by the analysis of ornamental plants that never naturalized and are no longer sold in Québec, which, as a group, are much less cold-hardy than the other plants. Nurseries and horticulturists of the 19th century T A B L E 2 Attributes (mean or median values, proportion of taxa) characterizing plant taxa listed in nursery catalogs published in the 19th century in Québec (Canada), according to whether the taxa were still sold (or not) in the province in 2015, and according to their status (naturalized or not). South America (% of taxa) 5.5 7.7 1.7 13.9 T A B L E 3 Standard logistic regression and phylogenetic logistic regression models explaining the naturalization, invasiveness, and disappearance (plants neither sold nor naturalized in 2015) of plant taxa listed in nursery catalogs published in the 19th century in Québec (Canada). Only significant variables are shown. T A B L E 3 ( Continued) America (Guo, Qian, Ricklefs, & Xi, 2006;van Kleunen et al., 2015;Rejmánek, 2014;Stohlgren et al., 2011), so it is not surprising to see European and Asian origin as significant explanatory variables in the invasiveness model. Plants from Asia (temperate) and/or Europe form 60% of the catalog taxa pool, but 86% of the invasive species identified in this study.
Plants that are not sufficiently cold-hardy are likely to "disappear" from the market or nature in Québec, but the disappearance phenomenon is more complex, involving low or short-distance seed dispersal abilities (e.g., autochory) and low propagule pressure (not widely sold).
To our knowledge, this study is the first to analyze the disappearance phenomenon from a large pool of ornamental plants. It is noteworthy that propagule pressure-as estimated from plant sale data-is often identified as a an important, if not the most important, determinant of naturalization for ornamental or cultivated plants (Dehnen-Schmutz et al., 2007;Moodley et al., 2014;Pemberton & Liu, 2009;Skou et al., 2012). However, in a cold region such as Québec, not being sold (low propagule pressure) is at least partially dependant on a lack of cold hardiness.
Factors contributing to naturalization or invasiveness may differ greatly between regions, as indicated by similar studies conducted in Australia, Central Europe, Hawaii, or Ireland (Gallagher et al., 2014;Milbau & Stout, 2008;Moodley, Geerts, Richardson, & Wilson, 2013;Schmidt & Drake, 2011). Differences rely in part on plant traits used in the models-some are almost always used (e.g., maximum height, seed mass), others rarely (specific leaf area), some included phylogeny, others not, etc. However, this study, conducted in a cold-temperate region, sheds new light on what is likely a context (climatic)-dependant phenomenon.
Using the naturalization model constructed in this study as a tool to predict the naturalization of a newly introduced plant would be risky, because only about a quarter of the variation was explained.
The invasiveness model-the most important from an environmental management perspective-is of limited use for the industry, as nurseries in Québec do not typically sell outdoor plants that are not, for instance, cold-hardy. On the other hand, it highlights the challenge this industry will face in an ever warming world: Hardiness zones are likely to shift northward over the next decades (Bradley et al., 2012), and several species currently sold could soon transition from casual to naturalized to invasive, causing additional pressure on native plant diversity. Regularly updating the cold hardiness zone maps would help rapidly flag new potential invaders, and banning the sale of invasive and weedy species in Québec-a list of such species has recently been compiled (Lavoie et al., 2014)-could be part of a solution. Unfortunately, there is actually no political will in the province to tackle this problem.