*4.1. Genotype* × *Environment Interactions*

Understanding genotype by environment (G × E) interactions is a necessary step for identifying and selecting poplar clones used for phytoremediation and associated phytotechnologies [68]. Poplar phenotypes are a function of their genotype, environment, and genotypic response to specific site conditions [49]. Phyto-recurrent selection has been used to choose superior poplar genotypes in the Midwestern United States [32,68]. Using both generalist and specialist genotypes enhances ecosystem services provided by phytoremediation applications. Deploying generalists with low G × E interactions and robust productivity across the region may be beneficial for cost and operational efficiencies [47,56], while specialists with high G × E interactions may maximize productivity, phytoremediation potential, and overall benefits of ecosystem services [26,28,49]. In the current study, sixteen phytoremediation buffer systems (i.e., multi-environmental trials (MET)) were established to evaluate trends in G × E interactions and identify generalist and specialist poplars in order to reduce runoff and clean groundwater (Figure 1).

In this study there were significant main (buffer, clone, and year) and interaction effects on tree health and growth parameters. In particular, interactions involving the buffer main effect were major factors governing clonal productivity. Buffer effects reflect tree responses to combined edaphic and local climatic conditions, and influence clonal performance traits such as: stem biomass production [69,70]; foliar fungal microbiomes [71]; leaf characteristics [72]; and diameter, height, and wood volume production [49,73]. In the current study, the broad spectrum of MET buffers with varying soil and climate conditions led to a wide range in clonal performance related to changes in both genotypic magnitude and ranks over time. More favorable climatic conditions (i.e., warmer, more precipitation; Table 1) and adequate soils for poplar cultivation (e.g., suitable texture, water-air properties, pH; Table 2) likely led to greater volume of most clones at Menomonee Falls, Slinger, Manitowoc and Caledonia. On the other hand, lower performance of clones at Marquette and Ontonagon can be attributed to less favorable climatic and soil conditions such as lower precipitation, temperature, and soil pH (i.e., more acidity). Similar results were

obtained by Hansen et al. [74] who showed that soil water availability played a key role in the productivity of woody biomass plantations. In the current study, sites with irrigation or shallow water tables exhibited the greatest wood volume, which was further corroborated through poplar biomass productivity modeling in the Midwestern United States [59,60]. Overall, G × E interactions resulted in mean annual increment (MAI) of four-year-old trees from the current 2017 buffer group ranging from 1.1 to 7.8 Mg ha−<sup>1</sup> yr−<sup>1</sup> (mean = 3.2 Mg ha−<sup>1</sup> yr−<sup>1</sup> ), which agreed with results for similarly-aged (i.e., 3 to 5 years) poplars in the region, whose MAI ranged from 0.6 Mg ha−<sup>1</sup> yr−<sup>1</sup> to 7.1 Mg ha−<sup>1</sup> yr−<sup>1</sup> [49,75].

Multi-environmental trials are key tools for defining gains achieved through identifying genotype characteristics, their stability, and relevance of their interaction with varying environmental conditions (i.e., G × E interactions) [76]. Although G × E interactions can have a significant impact on the precision of breeding value estimates, often resulting in decreased genetic gain [53], matching superior species and clones to particular site and growing conditions has been critical in maximizing the productivity of SRWC plantations [70]. Tree age is an important factor shown to govern G × E interactions for poplars. Although Riemenschneider et al. [48] found significant G × E interactions first occurring at three years after planting, Semerci et al. [52] recorded significant G × E interactions for growth and phenology traits in one-year-old poplar clones grown on sites with different water availability in Turkey. Similarly, in the present study, all tested traits exhibited G × E interactions after the first year in all three buffer groups (e.g., for one- to four-year-old trees). In contrast to the results presented here, greenhouse phyto-recurrent selection experiments with soils from the six phyto buffers of the 2017 buffer group [Bellevue (West), Caledonia (East), Menomonee Falls (East), Menomonee Falls (West), Slinger, and Whitelaw] showed a lack of G × E interactions for root-shoot ratio and growth performance index of many poplar clones tested at the current MET. Nevertheless, there were significant G × E interactions for tree health [58]. Such differences may be attributed to variability in environmental conditions between the greenhouse and field buffers and/or the length of the experiment (i.e., months versus years).

Regardless, such results have indicated that G × E interactions in poplar clones vary during the life cycle of the trees. Zalesny and Headlee [25] found significant G × E interactions for biomass and carbon production in both 10- and 20-year-old poplar plantations, despite negligible genotypic effects on both traits for 20-year-old trees. The presence/absence of G × E interactions within clones during the production cycle also can be expressed by variability in growth patterns across clones. Netzer et al. [77] recorded that some clones had greater biomass productivity in the second half of the stand rotation, while Ghezehei et al. [78] recorded both lack and presence of significant differences in clonal productivity of four- and eight-year-old poplars.

#### *4.2. Generalist and Specialist Response Groups*

Phenotypic responses determine comparative genotypic performance, resulting in some clones growing well and providing higher levels of ecosystem services across a broad range of soil, climate, and/or contaminant conditions (i.e., generalists). On the other hand, specialists optimize their growth and physiological processes when subjected to specific site conditions [43,44]. We identified both generalist ('DM114', 'NC14106', '99038022', '99059016') and specialist ('7300502', 'DN5', 'DN34', 'DN177', 'NM2', 'NM5', 'NM6') clones, along with others that exhibited volume consistent with both response groups across buffers and years ('9732-11', '9732-24', '9732-31', '9732-36', 'DN2') (i.e., those that shifted from generalists to specialists as trees aged; see below) (Table S13). Classification of these clones has important practical implications in reducing uncertainties associated with fielddeployment of these genotypes for multiple applications, including phytoremediation.

Changes in both magnitude and ranks across buffer × clone × year combinations defined the G × E interactions of the current study [43,44]. Different classifications were found for volume versus mean annual increment (MAI), which supports the need for longterm monitoring throughout plantation development [68]. One explanation for differences

between these traits may be related to the age at which the trees were measured. As noted previously, in the Southeastern United States, clonal rankings in poplar wood volume production changed with increasing stand age [70,79]. In this study, all volume estimates were from one- to three-year-old trees, while MAI was determined for trees after their fourth growing season, the start of the mid-rotation growth stage for poplars used for phytoremediation [80]. Similar changes were also apparent when evaluating measurement years within individual buffer groups. That is, oftentimes clonal rankings dramatically changed as trees aged (e.g., '7300502' had the greatest volume at Slinger during the establishment year only to have the least volume at this buffer after two and three growing seasons). A second explanation for differences in classifications between volume and MAI may be related to individual clones expressing higher levels of genetic variation and phenotypic plasticity as they responded to highly variable and changing soil conditions both within and across growing seasons at the phyto buffers. Guet et al. [72] reported that *P. deltoides* and *P. nigra* (which were the most common species used as parents in the current study) exhibited high levels of such genetic variation and plasticity, allowing them to better adapt to site-level spatial and temporal heterogeneity. These responses could lead to a greater propensity for specialist growth performance, and may explain why some NRRI clones ('9732-11', '9732-24', '9732-31', '9732-36') shifted from generalists to specialists in the current study (Table S13). In contrast, Nelson et al. [46] identified most of these clones as being geographically robust across both latitudinal and longitudinal gradients in North America. One of these clones, '9732-36', had very consistent volume production in our 2017 and 2018 buffer groups (as well as MAI in our 2017 buffer group), yet trended towards specialist responses at phyto buffers established in 2019. This may have been due to a negligible relationship between phenotypic plasticity and G × E interactions. As in our interpretation (and its definition), Des Marais et al. [81] linked G × E interactions to changes in clonal ranking and growth performance (i.e., variance-changing interaction) [56].

This concept of variance-changing interaction also supports differences in classifying generalist and specialist clones of the current study within individual measurement years associated with buffer × clone interactions for specific buffer groups. Specifically, individual response group designations for buffer × clone × year combinations (from Tables S6, S9 and S12) may differ somewhat from final classifications listed in Table S13, given the need to assess stability and magnitude of ranks within years and over time (i.e., classifying the clones holistically). With the exception of the NRRI clones, most other genotypes from the *P. deltoides* × *P. nigra* 'DN' genomic group (with 'DN2' being the only exception), as well as the *P. deltoides* 'D' clone '7300502' and all clones of the *P. deltoides* × *P. maximowiczii* 'DM' and *P. nigra* × *P. maximowiczii* 'NM" genomic groups exhibited consistent classifications across buffer groups. Nevertheless, of particular interest was that individual clones within genomic groups (or breeding groups, for the 'DN' hybrids) performed similarly, indicating that selection of genomic groups may be effective for early phyto-recurrent selection cycles (i.e., when choosing base populations for testing). Such genomic group trends have been reported for the same or related genotypes used in other phytoremediation applications [35,82]. Overall, the preponderance of specialist clones in the current study supports the need for phyto-recurrent selection in order to match genotypes to sites for small-scale applications with location-specific requirements (i.e., see variation in ranks for 'NM5' from Tables S6, S9 and S12), as well as the parallel need for continued testing of new genetic material, such as the NRRI clones, to select robust genotypes with minimal G × E interactions that can be used for large-scale, commercial applications at a justifiable cost while providing a multitude of ecosystem services across the rural to urban continuum [46].

**Supplementary Materials:** The following are available online at https://www.mdpi.com/article/10 .3390/f12040430/s1; **Table S1**: Dates of planting for each phytoremediation buffer system (i.e., phyto buffer); **Table S2**: Probability values from analyses of variance for health and mean annual increment (MAI); **Table S3**: Probability values from analyses of variance for height, diameter, and volume; **Table S4**: Height for the buffer × clone × year interaction (2017 buffer group); **Table S5**: Diameter

for the buffer × clone × year interaction (2017 buffer group); **Table S6**: Volume clone rank for the buffer × clone × year interaction (2017 buffer group); **Table S7**: Height for the buffer × clone × year interaction (2018 buffer group); **Table S8**: Diameter for the buffer × clone × year interaction (2018 buffer group); **Table S9**: Volume clone rank for the buffer × clone × year interaction (2018 buffer group); **Table S10**: Height for the buffer × clone × year interaction (2019 buffer group); **Table S11**: Diameter for the buffer × clone × year interaction (2019 buffer group); **Table S12**: Volume clone rank for the buffer × clone × year interaction (2019 buffer group); **Table S13**: Final classification of clones into generalist and specialist response groups; **Figure S1**: Health for the buffer × clone interaction measured in 2018 (2017 buffer group); **Figure S2**: Health for the buffer × clone interaction measured in 2019 (2017 buffer group); **Figure S3**: Health for the buffer × clone interaction measured in 2019 (2018 buffer group).

**Author Contributions:** Conceptualization, E.O.B., L.B., J.G.B., B.S.D., R.A.H., C.-H.L., M.P., R.S., A.H.W., and R.S.Z.J.; methodology, E.O.B., J.G.B., B.S.D., R.A.H., C.-H.L., B.G.M., A.P., E.R.R., R.S., R.A.V., A.H.W., and R.S.Z.J.; validation, A.P., E.R.R., and R.S.Z.J.; formal analysis, R.S.Z.J.; investigation, E.O.B., J.G.B., B.S.D., R.A.H., C.-H.L., B.G.M., N.D.N., A.P., E.R.R., R.A.V., A.H.W., and R.S.Z.J.; resources, L.B., J.G.B., R.A.H., C.-H.L., B.G.M., N.D.N., M.P., R.S., and R.S.Z.J.; data curation, B.S.D., E.R.R., R.A.V., A.H.W., and R.S.Z.J.; writing—original draft preparation, A.P., E.R.R., and R.S.Z.J.; writing—review and editing, E.O.B., L.B., J.G.B., B.S.D., R.A.H., C.-H.L., B.G.M., N.D.N., M.P., A.P., E.R.R., R.S., R.A.V., A.H.W., and R.S.Z.J.; visualization, E.R.R., and R.S.Z.J.; supervision, J.G.B., B.S.D., C.-H.L., N.D.N., E.R.R., A.H.W., and R.S.Z.J.; project administration, J.G.B., C.-H.L., and R.S.Z.J.; funding acquisition, A.H.W. and R.S.Z.J. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work was funded by the Great Lakes Restoration Initiative (GLRI; Template #738 Landfill Runoff Reduction).

**Acknowledgments:** The findings and conclusions in this publication are those of the authors and should not be construed to represent any official USDA or United States Government determination or policy. This work was funded by the Great Lakes Restoration Initiative (GLRI; Template #738 Landfill Runoff Reduction). We are grateful to S. Johnson, C.H. Perry, and N. Vrevich of the USDA Forest Service for GLRI support. We thank B. Dean (formerly), K. Corcoran, and M. Mujaddedi of the International Programs Office of the USDA Forest Service for administrative support of scientific exchanges between A. Pilipovi´c and R. Zalesny, as well as their administration and support of International Forestry Fellows: D. Karlsson (Sweden), P. Muñoz Gomez de la Serna (Spain), and A. Peqini (Albania). We thank the following site managers for access to their field sites: J. Forney and B. Pliska (Waste Management, Inc.); K. Dorow, D. Koski, and K. McDaniel (City of Manitowoc, Wisconsin); D. Henderson (AECOM Technical Services, Inc.); B. Austin and J. Wales (Marquette County Solid Waste Management Authority); D. Pyle (Delta County Solid Waste Management Authority); and S. Coron (Great American Disposal), as well as T. Beggs (Wisconsin Department of Natural Resources) for regulatory assistance and B. Sexton (Sand County Environmental) for insight, knowledge, and guidance on sustainable phytoremediation systems. We appreciate laboratory, greenhouse, and field technical assistance from: B.A. Birr, T. Cook, S. Eddy, F. Erdmann, C. Espinosa, D. Karlsson, R. Lange, P. Manley, M. Mueller, C. Munch, P. Muñoz Gomez de la Serna, D. Nguyen, A. Peqini, M. Ramsay, E. Schmidt, J. Schutts, and M. Wagler. We are grateful to J.M. Suvada for developing Figure 1, as well as R. Klevickas and P. Bloese for providing cutting material. We thank W.L. Headlee for statistical mentoring and reviewing earlier versions of this manuscript.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

