Early succession on slag compared to urban soil: A slower recovery

Slag, waste from the steel-making process, contains large amounts of calcium, magnesium, iron and other heavy metals. Because of its composition, high pH and low water retention ability, slag is considered inhospitable to plants. Nevertheless, the spontaneously generated plant communities on slag are surprisingly diverse, but the assembly and structure of such communities are poorly studied. Previous studies suggest reduced rates of succession due to low growth rate and slow accumulation of topsoil. To investigate whether slag communities display similar patterns, we used two former industrial sites on the South Side of Chicago, IL, both with high pH (8–9.2) sand content (80%) and calcium concentration (> 9000 ppm). We removed all vegetation from both slag and non-slag plots to test whether recovery differed over one growing season (4 months). To directly assess plant growth, selected focal species were planted on both sites and harvested. We show that recovery from removal differed at slag and non-slag sites: the recruitment process on slag, measured by percent vegetative cover and number of species in plots, was significantly slower at 6–8 weeks of the manipulation and beyond, suggesting a potential stage-dependent effect of slag on plant growth. Certain slag plots recorded less cover than non-slag plots by >30% at maximum difference. Functional trait analysis found that graminoid and early successional species preferentially colonized slag. Overall, slag plots recovered more slowly from disturbance, suggesting a slow succession process that would hinder natural recovery. However, slag also has the potential to serve as plant refugia, hosting flora of analogous habitats native to the area: one of our industrial sites hosts nearly 80% native species with two species of highest Floristic Quality Index (10). Restoration efforts should be informed by the slow process of natural recovery, while post-industrial sites in urban areas serve as potential native plant refugia.


Introduction
Human activities have drastically modified natural landscapes, creating many uniquely anthropogenic systems. One such example is the urban ecosystem, which encompasses the PLOS  that integrate into the soil, increasing organic material content [11][12][13]. Therefore, the ability of the primary successional community to provide organic material for soil is important for both the survival of later successional species and the recovery of an ecosystem from disturbance (facilitation model [14]). Without soil formation, many natural habitats display extremely slow or arrested primary succession. For instance, alvar habitats, comprised of thin topsoil on limestone bedrock, have high pH, shallow topsoil, and low nutrient availability similar to slag; due to the low growth rate of early successional plants, the community is "arrested" at a primary successional stage and does not succeed to later stages [15,16]. Similarly, reclamation efforts of forests on old mine sites are incomplete due to the failure of late successional and specialist species to colonize [17,18]. A slag site can be considered a habitat after disturbance, with most plants killed and soil covered by dumping. As disturbance ends (i.e. no further dumping) and erosion breaks the contiguous slag surface, primary succession takes place and the present community assembles. Previous studies have estimated successional trajectories of other industrial sites, predicting recovery to the original climax community in 75-150 years [19,20]. However, it is possible that plant communities on slag, though they have been relatively undisturbed for long periods of time, may have arrested primary succession. Because slag sites are usually characterized by a thin layer of topsoil with high pH, they might be comparable to alvars in displaying "arrested" soil and plant community development. The environmental stress on slag may lower growth rates of early successional plants, resulting in less accumulated plant litter and decreased soil formation [21]. This hypothesis is supported by the fact that some sites, including Big Marsh and Van Vlissingen Parks, have been undergoing succession from 1977 and 1927, respectively [5], but could not support survival of later-successional species, such as most trees and shrubs, except with an experimental supplement of compost [6]. Nevertheless, more evidence is needed to determine whether such slag plant communities are arrested in early successional stages by the lack of nutrients.
We tested the following hypotheses regarding the different succession process on slag: that 1) slag sites differed from non-slag sites by soil composition and community structure, 2) the recovery rate from disturbance, measured in percent cover, growth rate and number of recolonized species, differed for slag and non-slag sites, and 3) plants on slag differed in their functional traits that reflect adaptation to early successional habitats. Over the course of one growing season (June to October 2018), a combination of plant surveys, soil testing and controlled experiments on two types of sites, Slag and non-slag or urban soils (hereafter Reference), were used to test the three hypotheses. Overall, we expected to see a slower recovery of plants and a non-random assembly of early successional communities on slag.

Site selection and characterization
The study was conducted at two locales: Big Marsh Park (BM) and Van Vlissingen Park (VV), both Chicago Park District properties located on Chicago's Southeast Side (Fig 1; permit #1766). Neither site contains any remnant natural habitat but is part of a larger, highly modified wetland complex. The locales were chosen so that they both contained Slag (S) and nonslag Reference (R) sites, where R sites had urban soil deposited on top of fill of mixed origins [5], and were in close proximity to each other (within 1.9 km), to minimize variance in local climate. Slag at Big Marsh (BM-S) was deposited between 1965-1977 while slag at Van Vlissingen (VV-S) was deposited between 1902-1927 [5]. The reference site at BM (BM-R) has been seeded and managed for invasive species by Chicago Park District (See S1 File for details), and VV-R is surrounded by sparsely spaced cottonwood (Populus deltoides). All sites are in full sun, except for a few plots at VV-R which encounter partial shade for no more than 2 hours per day. Both slag sites contain slight depressions which allow standing water to accumulate after heavy rainfall. The experimental site at VV-S is surrounded by Phragmites-dominated shallow slag-bottomed wetlands.
The study spanned 4 months, from Jun 5, 2018 (planting focal species) to Oct 6, 2018 (final harvest). Weather data of study sites were drawn from the National Weather Service [23]. Highest temperature ranged from 17.8 (Jun 22) to 36.1˚C (Aug 4), and low temperature ranged from 8.9 (Jun 6) to 24.4˚C (Jul 1). Highest precipitation was 6.00 cm (Aug 7). Typically, weather at the study sites was highly variable temporally and geographically; during the study period, the weather was characterized by consecutive hot and sunny days, continuous rainfall over prolonged periods and occasional local thunderstorms.
Plant surveys were conducted on Jul 9, 2018 at VV-S and Jul 31, 2018 at BM-S using the censusing procedure from the Northwest Indiana Restoration Monitoring Inventory (NIRMI). At each site, a 50 m × 20 m survey plot was set up, and a comprehensive list of  species and cover data was recorded according to NIRMI procedure [24]. At BM-S, parts of the experimental plots were included in the survey. Species lists produced at BM-R and VV-R were in close proximity of experimental plots (� 10m), providing a regional species list for identification of species potentially germinating in experimental plots. Floristic Quality Assessment (FQA, [25]) was also performed on species lists using the Universal FQA Calculator [26] to evaluate conservation values at each site. FQA is used to estimate habitat quality based on the coefficient of conservation for each species present; the coefficient is determined by whether the species is disturbance-adapted (low score) or associated with undisturbed natural area (high score; [25]).
Samples of soil at each site were drawn from 10 different spots < 0.5 m from experimental plots, with an approximate depth of 10 cm. Samples were also taken from the commercial topsoil used in germination plots (see Experimental Setup: Experiment 1: Germination). Two distinctive types of soil were observed at VV-S, and were therefore sampled separately and designated VV-S-1 (block A-D) and VV-S-2 (block E). Samples from each site were combined and~400 mL soil was drawn from each mixture for laboratory analysis. A total of 6 samples (one per each site for BM-S, BM-R, VV-R; two for VV-S; one from commercial topsoil) were analyzed by A&L Great Lakes Laboratories (Fort Wayne, IN) using the environmental variables in Table 1.

Experimental setup
The basic unit of replication was the block, which consisted of five types of experimental plots to evaluate two general aspects of plant performance: three examined growth of established focal species and two tested germination success (Fig 2A). Each site hosted five blocks, named A-E. The order of plots within block was randomized. Two blocks on each Reference site contained only three focal species plots due to spatial constraints. Overall, each Slag and Reference site hosted 25 plots and 21 plots, respectively. Experiment 1: Germination. To evaluate which species are involved in recolonization after disturbance, all above-ground vegetation was removed from experimental plots (1.0 × 1.0 m, with additional buffer zone of 0.2 m on each edge), followed by either tilling and removal of 5-10 mm of extant soil and plant material (hereafter removal) or removal plus covering by 5-mm thick layer of commercial topsoil (New Plant Life All Purpose Topsoil, Markman Peat Corp., Le Claire, IA; hereafter topsoil). Removal plots simulated the recruitment of plants into newly disturbed habitat with no current residents; topsoil plots characterized the "background dispersal rate" because the effect of growth from seed bank or below-ground vegetative structures remaining in the soil was reduced.
To measure the establishment and growth rate of plants in germination plots, monitoring photos were taken weekly during the course of 16 weeks. An aluminum quadrat (0.25 × 0.25 m) was placed at one corner of the 1.0 × 1.0 m plot area, and a photo was taken for each corner. Each photo was then imported to ImageJ [27]. The percentage of green vegetation within quadrat (hereafter cover) was extracted; cover of a plot was then determined by the arithmetic average of cover values from all four photos. To calculate species diversity within germination plots, species presence/absence data for experimental plots was recorded from weekly monitoring photos selected from three dates, Aug 22, Sep 5 and Sep 26 (day 77, 91 and 112). A species was considered present if observed in any of the three records. To discover possible time trends in diversity patterns, the number of species present in each plot over the whole experimental period was also counted and recorded from monitoring photos. Experiment 2: Focal species. To directly measure the growth rate of plants on Slag and Reference sites, three native focal species with high conservation value were manually planted in experimental plots. For sideoats grama (Bouteloua curtipendula, BC) and showy goldenrod (Solidago speciosa, SS), each 1.0 × 1.0 m plot (with buffer zone of 0.2 m on each edge) hosted 16 plants in 4 × 4 pattern (Fig 2B); due to a sourcing constraint, only five common milkweed (Asclepias syriaca, AS) were planted in a 0.5 × 0.5 m plot (with buffer zone of 0.1 m on each edge; Fig 2C). Seedings less than a month old and were sourced from Cardno Native Plant Nursery (Walkerton, IN). All focal plants were watered daily during the first two weeks and at least once a week afterwards. No other manipulations (intensive weeding, application of fertilizer or pesticide) or physical protection were done.

Biomass harvest
The aboveground biomass was quantified twice for focal species plots and once at the end of season for germination plots. For focal species plots, an initial harvest was done at week 2 (Jun 20, 2018) where two plants (one for AS plots) were harvested for each focal species plot, and a final harvest was done at week 17 (Oct 6, 2018), when all plants in focal species plots were harvested except those harvested art week 2 and re-sprouted. Therefore, a maximum of 14 (4 for AS plots) plants were harvested from each plot, with the exact number depending on the survivorship. For germination plots, the central 0.25 × 0.25 m portion was harvested so that an equal portion of each corner that contributed to the cover measurement was included. After harvest, all biomass was left at room condition for no more than 12 hours before drying at 80˚C for at least 96 hours. The biomass in germination plots was measured as total biomass of all species in the plot; each focal species plant was measured individually.

Data analysis
All data analyses and plotting were done in R version 3.5.1 [28]. See Supplementary Methods in S1 File for a full list of packages used.
Soil, cover and biomass. Soil data were scaled, centered, and analyzed by principal component analysis (PCA). The resulting sets of vectors characterizing each soil sample were then extracted as independent variables in a linear regression of cover and biomass.
We tested whether site, locale or treatment of germination plots affected cover and biomass by one-way ANOVA. The relationship between cover and species number at each timepoint in removal plots was investigated by linear regression. Because the experiment involved a block design, linear mixed effects (LME) models were also fitted for block effect as a source of random error. Biomass estimates were divided into three groups: initial harvest, final harvest of focal species, and germination plots. All groups were analyzed using Welch's t-test and Wilcoxon's test when applicable, with biomass on Slag lower than on Reference as the alternative hypothesis. Data from final harvest of germination plots were further analyzed with two-way ANOVA of site and treatment effects.
We tested whether plant cover and biomass were related to soil properties using linear regression, with soil measurements extracted from PCA as independent variables. To compress the time series data of cover, the maximum was taken for each plot. Missing samples due to mortality were designated 0 in biomass analyses. Linear regression analyses used a Bonferroni correction; with 16 soil measurements, the adjusted p value was 0.05/16 = 0.003125.
Species diversity and community composition. Species presence-absence data from germination plots and on both Slag sites was compiled from three days of monitoring photos and the plant survey, respectively (see Experiment 1: Germination and S1 File). For both datasets, species richness (α diversity), plot dissimilarity (Whittaker's β diversity; [29]) and native status were analyzed. Non-metric multi-dimensional scaling (NMDS) analysis was performed on presence-absence data in germination plots to further evaluate the dissimilarity among sites and locales using the R package vegan [30]. Species number from each plot over the experimental period was also compiled from monitoring photos and subsequently tested for any site effects using ANOVA.
We also tested whether functional traits differed in communities that colonized removal plots on Slag and Reference, selecting traits that are important to plant production and reproductive strategies ( [13,31,32]; Table 2). Functional traits and native status for species present in germination plots were obtained from Hilty [33], the TRY database [34] and Grime et al. [35]. Canonical (constrained) correspondence analysis (CCA; [36]) was used to detect the potential correspondence between environmental variables and species distribution in germination plots. Model selection was based on maximum likelihood method using the R function ordistep in the package vegan [30] to obtain environmental variables that best explain species composition. To evaluate the association between functional traits and environmental variables, an RLQ model [37] was constructed with the package ade4 ( [38]; see S1 File for details).

Plant survey of slag sites
Slag sites at BM and VV hosted 66 and 44 species, respectively. Most species were herbaceous and relatively short in stature, forming a sparse cover. Dominant species included rough false pennyroyal (Hedeoma hispida), Dichanthelium spp. and rosette-forming forbs including fleabanes (Erigeron spp.), goldenrods (Solidago spp.) and plantains (Plantago spp.). Invasive European buckthorn (Rhamnus cathartica) and glossy buckthorn (R. frangula) are examples of shrubs on slag. Although seedlings of several trees were observed sporadically on slag, the only established trees are cottonwood (Populus deltoides), mulberry (Morus alba) and sumac (Rhus spp.; not in surveyed area); trees display visibly stunted stature compared to those growing on non-slag soil. Many invasive wetland species, such as common reed (Phragmites australis), cattail (Typha × glauca), and purple loosestrife (Lythrum salicaria) have either established monoculture or contribute largely to the plant community in wet depressions on slag.
FQA showed that VV-S is a higher quality habitat by hosting more native and high-conservation value species than BM-S. At BM-S, percentage of native species was 55.6%; the adjusted Floristic Quality Index (FQI) for the community was 19.4, with Solidago rugosa having the highest FQI of 6. At VV-S, percentage of native species was 79.4%; the adjusted FQI for the community was 40.1, with Carex crawei and Eleocharis elliptica both having the highest FQI of 10. See S2 File for the full species list and associated FQI.
We did not conduct surveys on both Reference sites because our focus was on slag communities. Reference sites have been seeded and/or undergone intense management by Chicago Park District. We have obtained lists of planted species on BM-R and a basic plant survey of VV-R (together in S10 File) conducted by Chicago Park District. Table 3 shows the main results of soil analysis. pH and Ca content of both Slag sites are considerably higher than those of both Reference sites (8-9.2 comparing to 7.2-7.9). Sand contents are also higher on Slag sites (72-80% compared to 50%), indicating lower water retention rate. BM-S had higher N and P content than all other sites and an organic matter content comparable to other Reference sites, in contrast to the expectation that slag sites are lower in nutrients. Heavy metals showed no consistent patterns between Slag and Reference sites, also counter to our expectation. PCA of soil measurements showed significant differences between Slag, Reference and commercial topsoil (CTRL; Fig A in S1 File). The first two principle components explained 42.99% and 29.57% of the total variance. Specifically, Slag samples were characterized by high Ca, Mn, K, Zn, Cr, sand content and higher CEC and pH ( Table 4).

Slag effects on recovery
Cover. Overall, Slag plots showed lower plant cover than Reference plots at both locales. ANOVA on cover showed a significant site effect between Slag and Reference with both treatments at both locales ( Table 5), but the effect of site was not distinctive until later in the experiment (day 49 to 63), as shown by time series plots of cover between sites of the same locale and treatment (a "time threshold" ; Fig 3A-3D; Table C in S1 File). Removal plots on both slag sites showed slower increase of cover with increasing species number, indicated by smaller slope of regression lines (Fig 4; Table E in S1 File). It is worth noting that there was a significant locale effect on cover across sites and treatments, except between removal plots on Reference sites ( Table 5), suggesting that the two locales cannot be treated as replicates. Furthermore, using linear mixed-effect (LME) models, both removal and topsoil plots at BM showed a significant block effect (σ block effect /σ overall = 2.422/7.308 and 2.552/8.163, respectively), indicating a high heterogeneity within these sites. See S5 File for the raw cover data. Species number. Removal plots on slag hosted fewer species at both locales (ANOVA on site effect, F = 26.18 and 47.82 for BM and VV, respectively; both p < 0.01), although the difference between Slag and Reference did not diverge until day 63 at BM, a pattern similar to the "time threshold" reported for the site effect of cover (Fig 3E and 3G, Table C in S1 File; see S6 File for the raw cover data.).
Biomass. Harvested aboveground biomass on Slag sites was lower than Reference at both BM and VV for both germination and focal species plots. The difference was significant for germination plots using Welch's t-test on log-transformed data (Fig 5D-5E; Table F in S1 File). Two-way ANOVA on site (Slag or Reference) and treatment (removal or topsoil) showed a consistent effect of slag (F = 21.335 and 33.53, respectively; both p < 0.001) but not treatment (between removal and topsoil; see Table G in S1 File).
During initial harvest, focal species plots showed little effect of site or locale: while the biomass of BC on Slag was higher than on Reference (both Welch's and Wilcoxon's tests, p < 0.05); the significance was not robust using log-transformed biomass data (Table F in S1 File). At the end of the growing season, the biomass on Slag sites was consistently lower than on Reference soils for BC and SS using both Welch's and Wilcoxon's tests (p < 0.01; Fig 5B  and 5C). Notably, the high mortality (43/70, or 61.4%) of SS on BM Slag resulted in a median biomass of 0. The strong effect of slag was still detectable using only surviving SS biomass. The high mortality of AS across all sites and locales (54/80, or 67.5%) led to little difference in final biomass and therefore was not informative (Fig 5A). No mortality was observed for BC at any site. BC at both BM and VV and SS at VV showed a significant block effect (σ block effect /σ overall = 0.7642/1.929, 1.382/3.538 and 0.38/1.607, respectively).
Correlation with soil measurements. Slag soil variables mostly correlated negatively to plant growth, measured in both cover and biomass when cover and biomass were fitted to each soil variable with linear regression (LR) models using both untransformed and log-transformed soil measurements. The transformation of soil variables yielded no qualitative differences on LR results. With only a few exceptions, most variables in category I (Slag; see Table 4) negatively correlated with plant growth, measured in either cover or biomass, of both experiments; most variables in category II (Reference) and III (Topsoil) positively correlated with plant growth (Table J in S1 File).

Structure of recolonized communities
In the original experimental design, topsoil plots were set up to minimize germination of the local seed bank and to characterize the background dispersal rate. However, the commercial topsoil was not sterile, biasing species diversity and community composition estimates. Analyses on community structure (see below) excluded data from topsoil plots (see Supplementary Methods in S1 File). Table 5 Early succession on slag β Diversity. β diversity (as Whittaker's β) showed that removal plots on the same site had high similarity (Fig B in S1 File), in contrast to Slag or Reference plots across locales. Removal plots differed in species composition across sites (Fig 6A). Permutational Multivariate ANOVA (PERMANOVA) results suggested a significant difference between four sites (R 2 = 0.726, p = 0.001) but not among blocks within each site (see Table H in S1 File for full statistics). Fig D in S1 File shows the distribution of species and sites in NMDS space. Most species resided in clusters, indicating their constrained distribution at some sites; cottonwood (Populus deltoides, PODE) and false pennyroyal (Trichostema brachiatum, TRBR) were present in plots across multiple sites and therefore were distant to all existing site clusters.

. ANOVA results of site (between slag and non-slag sites) and locale effects (between BM and VV) on cover in germination plots
Functional traits. Results from CCA showed significant clustering for plant species (Fig  6B). Similar to the results from NMDS, each cluster of plots had a distinctive collection of species. Model selection using all environmental variables identified soil K, Mn and N as the three most significant variables determining species distribution. Limiting environmental variables to category I (Slag) and category II + III (Reference; see Table 4) yielded different results: for Slag, a collection of soil Ca, pH and Mn explained most of the variation in species distribution; for Reference, the explanatory environmental variables changed to organic matter, N and Mg.
Fourth-corner analysis showed some signal of environmental filtering on functional traits (Fig 6C). Graminoid species (Funct.G) and species with widespread seeds (Regen.W) associated most strongly with organic matter, P, Mg, Fe, As (negative) and K (positive). Late-summer flowering species (after July; Pheno.L) were positively associated with soil N, Zn and Cr; the reverse was true for summer-flowering species (before July; Pheno.S). Specific leaf area (SLA) was negatively associated with soil pH and Cr and positively associated with K. Although some individual associations were significant, the overall association was not significant as determined by the S RLQ statistic (p > 0.50). Seed bank longevity and lateral spread Percent cover and the number of species in germination plots between Slag (dashed line) and Reference (solid line) sites at each locale are plotted. x axes show time from start of experiment in days. Shaded areas denote one standard deviation. A dip in cover at Big Marsh can be seen at day 70, shown by arrows. Notation above each graph shows the site-effect ANOVA result using data from that specific day. Significance levels: p > 0.05 (ns), p < 0.05 ( � ), 0.01 < p < 0.05 ( �� ), p < 0.001 ( ��� ). See Tables C-D in S1 File for complete statistics.
https://doi.org/10.1371/journal.pone.0224214.g003   Table F in S1 File for complete statistics. All photos were taken by the authors. https://doi.org/10.1371/journal.pone.0224214.g005 Early succession on slag were excluded from this analysis because of not enough data. See S4 File for complete functional traits data.
Native status. Although Reference removal plots had higher species richness than Slag removal plots (mean: 8.33 and 4.90, respectively), the latter had a higher proportion of native species (Welch's t-test, p = 0.0357), a difference especially apparent between VV-R and VV-S removal plots (Fig 7A). Topsoil plots generally had a higher number of nonnative species, likely due to seed contamination, and the difference between Slag and Reference was not significant (Welch's two-sided t-test, p = 0.3449; Fig 7B). Removal plots on both Reference sites had a higher number of native species than topsoil plots on slag sites (p < 0.001). Additional investigation showed that topsoil was a source of nonnative species (see S1 File).

Discussion
Slag sites, among many other fragmented urban habitats with unique environmental conditions, offer insight into the development of novel communities through environmental filtering and adaptation to the urban environment [1]. Community composition and growth of self-assembled vegetation on contaminated sites informs conservation methods such as phytoremediation [3,6,7,20]. Through manipulative experiments and functional trait analysis, we characterized the plant community assemblage on slag from two perspectives: 1) how environmental variables affect plant growth and subsequently succession, and 2) how environmental variables affect plant community assembly. Although the heterogeneity between our locales and blocks within each slag site contributed significantly to cover and biomass results, several consistent patterns were observed. Plant growth and community recovery, measured in percent cover, species number and biomass, were significantly hindered by slag, and community composition on slag did not represent communities on non-slag soil of close proximity. Further studies of community processes on slag would benefit greatly from longer study period, a larger number of soil samples and in situ measurements of functional traits.

Slag effect on plant growth
Plant growth, measured by both percent cover and biomass, was negatively affected by slag. Overall, Slag and Reference sites differed for all measurements of plant growth, both cover and biomass. (Figs 3 and 5). Moreover, recovery of Slag removal plots was slower in terms of both Early succession on slag cover and number of recolonized species (Fig 4). Combining these results with previous studies that discovered arrested primary succession on natural habitats with similar soil profiles [15,16], it is reasonable to infer that with such low growth rate, slag vegetation might also experience arrested, or at least delayed primary succession. Linear models showed that environmental variables characterizing slag (category I in Table 4) consistently correlated negatively with plant growth. Nevertheless, the small environmental sample size resulted in low R 2 and p values (Table J in S1 File); the patterns could also arise from other sources, including positive association between As and plant growth might be due to a high As content in nutrient-rich commercial topsoil.
In contrast, germination and early growth was not affected by slag. A "time threshold," before which cover and species number did not show a significant site effect, was observed at day 49-63 (Fig 3; Table C-D in S1 File). Previous literature has equivocal evidence on the inhibitory effect of heavy metal ions on seed germination and seedling growth [39][40][41][42][43], though these studies were conducted in a lab setting with one or a few species. Environmental conditions in the field could be more complex, and populations on slag might be more tolerant to heavy metal contamination because of previous exposure.
Additionally, a significant dip in BM-S cover indicates a heat wave on Aug 4 (arrows in Fig  3A and 3B). which also caused a more than 50% mortality of SS on that site (Fig 5C). Under such heat, low water retention rate caused by high sand content at BM-S might have resulted in greater desiccation [12].

Slag effect on community structure
The concept of environmental filtering describes community assembly shaped by environmental conditions of habitats [44,45]. Theories on environmental filtering predict that environmental factors shape community composition by selecting for species that confer higher fitness, resulting in clustering [46][47][48][49]. Studies have shown that this clustering effect is more visible for the functional traits a species possesses than species identity [47], and at larger rather than smaller spatial scales [48][49][50]. Thus, if slag imposes environmental filtering, species composition or functional traits, or both, should be clustered.
Although species composition of removal plots was more similar within each site, they were different between Slag sites (Fig 6A and Fig B in S1 File). This observation is partly explained by the CCA result (Fig 6B), which showed significantly different soil compositions between the two Slag sites. This effect due to locale difference potentially selected for different assemblages in removal plots. The clustering was also observed in functional traits (Fig 6C). Graminoid species (grasses, sedges) and species with high SLA aggregated with high K content, which characterizes BM-S in CCA. Graminoids are iconic early-successional species [11,12], and a high SLA implies high investment in primary production, often a characteristic of fast growing, early successional species [13]. Therefore, BM-S hosted more early-successional species than VV-S. Late-summer flowering species aggregated with high N, Zn and Cr content, all of which characterize VV-S. Late flowering is associated with either competitive or ruderal plant species [51] and has been shown to associate with higher relative reproductive success [52]. Generally, plants colonizing removal plots on Slag were early successional species with high growth rates, consistent with the prediction. Compared to Slag plots, Reference plot species are not characterized by early-successional traits, likely due to regeneration from the existing, later-successional seed bank.
However, it is worth noting that the functional trait data analyzed were obtained from a database, not in situ. Therefore, failure to detect strong evidence of environmental filtering might be due to changes of plant functional traits on slag due to plasticity [53,54]. For instance, all plants in experimental plots showed slower growth on slag, suggesting traits associated with production such as SLA and LDMC might not be accurately captured by measurements from an online database. Trait measurements through time would provide more information on community assembly, resilience and the importance of functional trait plasticity on slag.

Implications for conservation
The succession status of slag sites could be important to conservation and management practices. In previous phytoremediation efforts, low survival undermined the ability of plants to take up contaminants [6]. In contrast, studies have shown that industrial ecosystems such as sand-gravel pits, alkaline waste, and limestone quarry floors have been colonized by species suitable for these habitats and show a normal succession process [3,16,19,55]; these authors state that these industrial sites are able to provide habitat for natural successional processes without remediations such as seeding or topsoil capping. Indeed, Smith et al. [19] predicted that a slag dump from 1918 would have accumulated enough organic material to support a pine forest in 75 years, while Ř ehounková and Prach [20] predicted that 25 years of natural succession would restore gravel-sand pits aged 1-75 years back to grassland, woodland or wetland. However, these conclusions are derived from systems that contain enough organic material [19], where the recovery of community is to the target climax community [20], or when successful establishment of several species can start succession [55]. None of them explicitly addressed the slow growth rate of recolonized species and its potential effect on soil formation argued by Stark et al. [15].
Although many species successfully colonized slag, it remains unknown whether the community would continue along a successional trajectory. Empirical evidence from both germination and the focal species experiments suggests that species on slag developed more slowly than on non-slag soils; slag plots were colonized by fewer species and accumulated less biomass. Therefore, the succession process on slag is expected to be slow. Furthermore, the goal of restoration for slag sites is hard to define. Embedded in the urban matrix, the "natural community" of slag sites is not readily identifiable. In our region of study, the pre-development communities were likely tallgrass prairie, wetlands, or woodlands. Rehabilitation of slag to its "original state" could be difficult, costly, and incur risks due to the introduction of weedy species in commercial topsoil, as we demonstrate here (Fig C in S1 File).
Nevertheless, urban and industrial habitats might serve as plant refugia, providing an alternative motivation for restoring slag habitats. Because of its similarities to naturally occurring habitats such as alvar or dolomite prairie, slag sites have the potential of hosting rare plants that could only inhabit such environments [3]. For instance, Tomlinson et al. [16] found significant overlap between the flora of alvar habitat and artificial quarry floors through natural colonization. Furthermore, species that are adapted to natural habitats with high metal content, such as serpentine soil or mining sites, could establish in industrial ecosystems [3,56,57]. Inhabitable by many competitors, habitats such as slag could provide a refuge for those rare species. A relatively high proportion of native flora and many species with a very high conservation index were identified at both BM-S and VV-S; more strikingly, the latter has a high adjusted FQI, suggesting habitat quality comparable to natural habitats (40.1 compared to 30-40 locally [58][59][60]; see, however, [61] for a discussion and potential drawbacks of using FQA). Depressions on VV-S hold water during the rainy season, supporting many native fen and marsh species such as Carex crawei and Eleocharis elliptica; the native orchid Spiranthes cernua occurs at both VV-S and BM-S; native plants with high conservation values planted on VV-S, Bouteloua curtipendula and Solidago speciosa, both displayed 100% survival. These results suggest that VV-S has potential as an urban refuge for rare plants, especially native species on analogous habitats such as wetlands and dolomite prairie. Thus, habitat reconstructions that mimic a natural analog of slag sites, such as alvar or dolomite prairie [3,16] might have the greatest success. Those habitats formed by limestone outcrops are native to the Great Lakes region, hosting many species that are adapted to their environments [62][63][64]. Given the similarity of environmental conditions and many overlapping species such as Carex spp., Eupatorium spp. and Panicum virgatum, Verbena spp., Bouteloua curtipendula, VV-S could be reclaimed with a plant assemblage that mimics these natural habitats, while retaining its unique species of high conservation values. Future studies may introduce native plants specialized on similar ecosystems to slag and evaluate the feasibility of slag as native plant refugia.

Conclusions
Although slag has traditionally been viewed as a contaminated wasteland, it is a unique urbanindustrial ecosystem that has the potential to host a unique flora. We showed that unfavorable environmental conditions significantly lowered the growth and recovery rates of slag communities in terms of percentage cover, biomass and number of recolonizing species. The composition of slag communities generally corresponds to early successional communities, but the low growth rate may significantly reduce the rate of succession. Unlike many other industrial systems such as sand-gravel pits, natural recovery of slag sites might not be feasible in the short term, requiring active restoration efforts. While topsoil capping might be the most effective method of increasing organic matter, it risks introducing nonnative species, further lowering habitat quality. Although challenging for restoration and remediation, slag could potentially host rare plants from natural, analogous habitats, including flora found at native Midwest habitats such as dolomite prairie and alvar. Therefore, in addition to "radical" remediation such as topsoil capping, burning, and intensive weeding of undesirable species, further efforts to maintain resilient urban ecosystems should also consider retaining slag sites as potential refuges for native plants, which in turn will contribute to the conservation of regional biodiversity.
Supporting information S1 File. Supplemental information: Methods and results. Detailed history of management for both Reference sites and methods for data analysis, including the list of R packages used. All supplemental results ( Table A- S4 File. Functional traits of plant species observed in the survey used for data analysis. Species abbreviations follow those used in S3 File. If one species has multiple measurements of one single trait, measurements were either averaged (numerical) or entered together (categorical). Numerical data with multiple measurements were averaged. Missing data were left blank and substituted with NA in actual analyses. For all analyses, "Type", "Growth Form" and "Native?" were taken out; in addition, "Seedbank Longevity" and "Lateral Spread" were taken out when performing RLQ analysis and fourth-corner method because of low availability of data. Data were unavailable for the two unidentified Dicanthelium species and Cyperus bipartitus.
(CSV) S5 File. Cover data. Cover in all removal and control plots is measured from weekly monitoring photos. Time is in days from start of the experiment.