Population BiologyOpen Access icon OPENOpen Access license

Spatial and Temporal Analysis of Populations of the Sudden Oak Death Pathogen in Oregon Forests

    Affiliations
    Authors and Affiliations
    • Z. N. Kamvar
    • M. M. Larsen
    • A. M. Kanaskie
    • E. M. Hansen
    • N. J. Grünwald
    1. First and fourth authors: Department of Botany and Plant Pathology, Oregon State University, Corvallis; second and fifth authors: Horticultural Crop Research Unit, United States Department of Agriculture–Agricultural Research Service, Corvallis, OR; third author: Oregon Department of Forestry, Salem; and fifth author: Department of Botany and Plant Pathology and Center for Genome Biology and Biocomputing, Oregon State University, Corvallis.

    Abstract

    Sudden oak death caused by the oomycete Phytophthora ramorum was first discovered in California toward the end of the 20th century and subsequently emerged on tanoak forests in Oregon before its first detection in 2001 by aerial surveys. The Oregon Department of Forestry has since monitored the epidemic and sampled symptomatic tanoak trees from 2001 to the present. Populations sampled over this period were genotyped using microsatellites and studied to infer the population genetic history. To date, only the NA1 clonal lineage is established in this region, although three lineages exist on the North American west coast. The original introduction into the Joe Hall area eventually spread to several regions: mostly north but also east and southwest. A new introduction into Hunter Creek appears to correspond to a second introduction not clustering with the early introduction. Our data are best explained by both introductions originating from nursery populations in California or Oregon and resulting from two distinct introduction events. Continued vigilance and eradication of nursery populations of P. ramorum are important to avoid further emergence and potential introduction of other clonal lineages.

    Sudden oak death (SOD) emerged as a severe epidemic disease on coast live oak (Quercus agrifolia) and tanoak (Notholithocarpus densiflorus) in California in the mid 1990s and reemerged shortly thereafter on tanoak in Oregon in the early 2000s (Everhart et al. 2014; Grünwald et al. 2008a; Hansen et al. 2008; Rizzo et al. 2005). SOD is caused by Phytophthora ramorum Werres, De Cock & Man in’t Veld and is considered to be one of the top two oomycete pathogens based on its scientific and economic importance (Kamoun et al. 2015; Werres et al. 2001). The Oregon epidemic was first detected during aerial surveys in 2001 on tanoak but likely derived from initial introductions in the late 1990s. The Oregon Department of Forestry has since monitored the epidemic and sampled symptomatic tanoak trees since 2001 (Hansen et al. 2008). Strains sampled from infected sites in forest or nursery environments have been genotyped in several labs using a range of microsatellite loci (Grünwald et al. 2009; Ivors et al. 2006; Prospero et al. 2004, 2007, 2009).

    P. ramorum has emerged repeatedly around the world as four distinct clonal lineages found in North America (lineages NA1, NA2, and EU1) and Europe (EU1 and EU2) (Grünwald et al. 2012; Ivors et al. 2006; Van Poucke et al. 2012). The lineages have been named by the continent on which they first appeared (i.e., North America = NA and Europe = EU) and are numbered in order of discovery (Grünwald et al. 2009). The NA1 clonal lineage was first discovered in California causing SOD on tanoak and coast live oak and is the one currently found in Curry County, OR (Mascheretti et al. 2008). The EU1 and NA2 populations were discovered later in nursery environments and are currently only found in California, Oregon, Washington, and British Columbia while the NA1 clone has been shipped with nursery plants from the west to the southern and southwestern United States (Goss et al. 2009, 2011; Grünwald et al. 2012; Ivors et al. 2006; Mascheretti et al. 2008; Prospero et al. 2009). The EU1 clonal lineage is the one first discovered in Europe but, in 2007, the new EU2 lineage emerged in Northern Ireland and has since migrated to western Scotland (Van Poucke et al. 2012; Werres et al. 2001). EU1 was first introduced to Europe and eventually migrated to the Pacific Northwest of North America (Goss et al. 2011).

    P. ramorum populations sampled in Oregon forests to date belong exclusively to the NA1 clonal lineage (Hansen et al. 2008; Prospero et al. 2007). Given that NA2 or EU1 clones have been found in California, Oregon, Washington, or British Columbia in association with nursery plant movements, introduction of NA2 or EU1 from nursery environments to Curry County forests is a plausible scenario (Goss et al. 2009, 2011; Grünwald et al. 2012; Prospero et al. 2007, 2009). Thus, our present work monitors populations and potential emergence of novel lineages in Oregon forests.

    Our main objectives here are to describe the spatial and temporal pattern of the populations and clonal dynamic of the SOD pathogen in Curry County in southwestern Oregon from 2001 to the present. Specifically, we asked whether (i) novel lineages have been introduced into the forests in Curry County, (ii) multiple introductions occurred, and (iii) introduction might have come from nursery populations. We sampled infected tanoak trees between 2001 and 2014 and characterized populations using microsatellite analysis.

    MATERIALS AND METHODS

    Location.

    The SOD-infested areas are located in the Siskiyou Mountains of Curry County in southwestern Oregon near the town of Brookings (42.0575° N, 124.2864° W) on the coast (Fig. 1) (Prospero et al. 2007). The Siskiyou mountains form part of the Klamath Mountain range (Franklin and Dyrness 1988). The vegetation in SOD-infested areas is a mosaic of different vegetation types, including mixed-evergreen, redwood (Sequoia sempervirens), and Douglas-fir (Pseudotsuga menziesii) forests, with tanoak as the dominant SOD host.

    Fig. 1.

    Fig. 1. Spatial distribution of the sudden oak death epidemic and multilocus genotypes of Phytophthora ramorum in Curry County, OR. Crosses mark tanoak trees found positive for P. ramorum during aerial surveys. In total, 70 multilocus genotypes have been identified between 2001 and 2014 and are marked as shown in the legend. Abbreviations for regions shown in the map: JHallCr = Joe Hall Creek; NFChetHigh = North Fork Chetco; Coast = Coastal Region; HunterCr = Hunter Creek, Cape Sebastian; and PistolRSF = Pistol River South Fork. Inset shows the placement of Curry County (dot) in southwestern Oregon.

    Download as PowerPoint

    Sampling.

    Commencing in 2001, two to four aerial surveys per year were conducted over the tanoak range by the Oregon Department of Forestry and the United States Department of Agriculture Forest Service in Curry County. The survey detects recently killed tanoak trees based on the reddish-brown color of foliage (Hansen et al. 2008). All trees identified by aerial surveys were ground checked and geographically referenced using a hand-held GPS instrument (Garmin GPS 12XL or 60CX; Garmin International, Olathe, KS). Bark or foliage samples were collected for determination of Phytophthora ramorum presence by culturing in the field and laboratory. Host plants within the area of the delimitation survey, generally 300 feet, were also inspected and sampled if they were symptomatic. Maps of distribution were prepared using ArcView GIS version 3.3 and ArcMap version 10.2 (Environmental Systems Research Institute, Redlands, CA).

    Isolation, identification, and DNA extraction.

    Isolations were made from symptomatic plant tissue onto selective CARP agar (Difco corn meal agar, natamycin at 10 ppm, NA-ampicillin at 200 ppm, and rifampicin at 10 ppm) (Prospero et al. 2007). Candidate Phytophthora cultures were transferred onto corn meal agar with β-sitosterol at 30 ppm. P. ramorum identification was confirmed by microscopic inspection for presence of characteristic chlamydospores and deciduous sporangia (Werres et al. 2001). Genomic DNA was extracted using either the FastDNA SPIN kit (116540600; MP Biomedicals, LLC) (Goss et al. 2009) or the cetyltrimethyl ammonium bromide-chloroform-isopropanol method (Winton and Hansen 2001). An overview of strains collected by year and region (following regions as shown in Figure 1) is provided in Table 1.

    TABLE 1. Summary of Phytophthora ramorum isolates sampled in Oregon forests and multilocus genotypes (MLG) observed across regions and years

    Genotyping, data validation, and harmonization.

    Five microsatellite loci were utilized in this analysis: PrMS6, Pr9C3, PrMS39, PrMS45, and PrMS43 (Grünwald et al. 2008b, 2009; Prospero et al. 2004, 2007). Genotyping (see specific protocols in the Supplementary Text) of P. ramorum strains collected from 2001 to 2012 occurred over several years and in several laboratories, with different protocols and sequencers. Consequently, a concerted effort was made to create a comprehensive dataset with identical allele calls. To detect errors, allele calls from all five genotyped loci were generated for a subsample of 40 isolates representing the most common multilocus genotypes (MLGs) from the culture collection, then compared with data from participating laboratories. Three of the five loci (PrMS6, Pr9C3, and PrMS39) had identical allele calls between laboratories for the subsampled isolates. The remaining two loci (PrMS45 and PrMS43) had allele calls that differed by a single base pair between laboratories. Therefore, data from PrMS45 and PrMS43 were corrected to allow consistent comparisons of allele calls. Given the varied nature of the genotyping data described above, genotyping of P. ramorum strains consists of two datasets: one including 5 loci (2001 to 2014) and the other (a newly developed, multiplexed method) including 14 loci (samples 2013 and 2014) (Table 2). Details on both genotyping methods can be found in the Supplementary Text and Supplementary Figure S1.

    TABLE 2. Newly multiplexed protocol for Phytophthora ramorum primer sequences of simple-sequence repeat (SSR) loci and final concentrations (conc.) used to determine multilocus genotypes for four clonal lineagesa

    Nursery populations.

    To determine whether forest populations cluster with different nursery populations from Oregon or California, we used previously published data from our work to determine relationships among nursery and Curry County forest populations (Goss et al. 2009, 2011; Grünwald et al. 2009; Prospero et al. 2007, 2009).

    Data analysis.

    All individuals genotyped for this effort belonged to the NA1 clonal lineage (Grünwald et al. 2009). Thus, all analyses presented here focused on describing the clonal dynamic using model-free approaches that avoid violation of population genetic theory. Samples were grouped into different MLGs defined by the unique combination of alleles across all observed loci from the consensus five simple-sequence repeat loci genotyped across all years. For identification purposes, unique MLGs were then assigned an arbitrary number from 1 to the total number of observed MLGs. Population genetic analysis was conducted using the computer and statistical language R (R Core team 2014) using various packages as well as R functions written specifically for this project (see github link below). Graphs and figures were created using the R packages ggplot2, ape, igraph, ggmap, and poppr (Csardi and Nepusz 2006; Kahle and Wickham 2013; Kamvar et al. 2014b; Paradis et al. 2004; Wickham 2009). Within-locus allelic diversity was analyzed across and within years and regions using the function locus_table from the R package poppr (Supplementary Table S1)(Kamvar et al. 2014b). To address the temporal and spatial aspects of the data, populations were analyzed both by year isolated and watershed region (Table 1; Fig. 1). Watershed regions were drawn with ArcMap (version 10.2; (Environmental Systems Research Institute). The regions represent drainages or portions of drainages in which infected trees were discovered as the disease progressed over time. In most cases, ridgelines dividing drainages formed the boundary of a region. These regions were saved as shapefiles and imported into R with rgdal (Bivand et al. 2014).

    Genotypic diversity was analyzed within and across years and populations with the Shannon-Wiener index (H) and the Stoddard and Taylor’s index (G), (Shannon and Weaver 1949; Stoddart and Taylor 1988). Both G and H measure genotypic diversity, combining richness and evenness. If all genotypes are equally abundant, then the value of G will be the number of MLGs and the value of H will be the natural log of the number of MLGs. Both G and H are used because they weigh more or less abundant MLGs more heavily, respectively (Grünwald et al. 2003). Evenness was calculated as E5, which is an estimator of evenness that utilizes both H and G that gives a ratio of the number of abundant genotypes to rare genotypes (Grünwald et al. 2003; Ludwig 1988; Pielou 1975). These were calculated with the R packages poppr and vegan (Kamvar et al. 2014b; Oksanen et al. 2013). Confidence intervals were calculated using the R package boot with 9,999 bootstrap resamplings (Canty and Ripley 2015). Richness, or the expected number of MLGs (eMLG), was calculated using rarefaction from the R packages poppr and vegan (Heck et al. 1975; Hurlbert 1971). Some statistics (analysis of molecular variance [AMOVA], genotypic diversity, index of association, allelic diversity, and Nei’s distance) were also performed on clone-censored data where each MLG was represented once per population hierarchy.

    Because the analysis of genotypic diversity, richness, and evenness is agnostic to specific alleles within MLGs, assessment of genetic relatedness between MLGs was performed using the function bruvo.dist using poppr, which calculates Bruvo’s genetic distance utilizing a stepwise mutation model for microsatellite loci (Bruvo et al. 2004; Kamvar et al. 2014b). Thus, this distance gives a more fine-scale picture of relationships between individuals than band-sharing models. These relationships were visualized with minimum spanning networks generated using the R packages igraph and poppr (Csardi and Nepusz 2006; Kamvar et al. 2014b).

    If the epidemic had a single origin, a correlation between genetic and geographic distance would be expected because populations acquire mutations over time and clonally diverge regardless of rates of spread. Divergence is affected by rates and distance of spread, where long-distance dispersal or low rates of mutation would lead to less divergence and, thus, lower correlations between genetic and geographic distances. This was tested by performing Mantel tests across all hierarchical levels in the data set utilizing the function mantel.randtest in the R package ade4 between Bruvo’s distance as described above and Euclidean distances between geographic coordinates (Dray and Dufour 2007; Mantel 1967). P values were calculated using 99,999 bootstrap replicates.

    Because eradication efforts destroy the immediate habitat in an infected area, one question that we wanted to address was whether or not genotypes were clustering to specific regions or evenly spread throughout Curry County (Prospero et al. 2007). This was tested using three methods: bootstrap analysis of Nei’s genetic distance, AMOVA, and discriminant analysis of principal components (DAPC) in the R packages poppr, ade4, and adegenet (Excoffier et al. 1992; Jombart et al. 2010; Kamvar et al. 2014b). The bootstrap analysis utilized 10,000 bootstrap replicates, treating loci as independent units with the function aboot in poppr, and was visualized as an unrooted neighbor-joining tree in figtree v. 1.4.2 (Fig. 2). AMOVA utilizes a distance matrix between genotypes for which hierarchical partitions are defined, and attempts to analyze the variation within samples, between samples, between subpopulations within populations, and, finally, between populations. In this case, we used the hierarchies of both samples within years within regions and samples within regions within years. DAPC is a multivariate, model-free approach to clustering based on prior population information (Jombart et al. 2010). This allowed us to analyze the population structure by assessing how well samples can be reassigned into previously defined populations. Both DAPC and AMOVA were run with and without Hunter Creek and Pistol River South Fork due to isolated genotypes and small sample size, respectively. For the DAPC analysis, these removed populations had their origins predicted from the DAPC object using the function predict.dapc in the R package adegenet (Jombart et al. 2010).

    Fig. 2.

    Fig. 2. Minimum spanning network based on Bruvo’s genetic distance for microsatellite markers for Phytophthora ramorum populations. Nodes (circles) represent individual multilocus genotypes. The 10 most abundant forest genotypes are labeled with their multilocus genotype (MLG) designation. Node colors/shades represent population membership proportional to the pie size. Node sizes are relatively scaled to log1.75n, where n is the number of samples in the nodes to avoid node overlap. Edges (lines) represent minimum genetic distance between individuals determined by Prim’s algorithm. Nodes that are more closely related will have darker and thicker edges whereas nodes more distantly related will have lighter and thinner edges or no edge at all. Reticulation was introduced by finding exact ties in genetic distance after Prim’s algorithm was run. Subgroups of >3 MLGs where all nodes are no more than one mutational step (d = 0.05) away from its neighbors are highlighted in various shades.

    Download as PowerPoint

    Because DAPC is sensitive to the number of principal components used in analysis, the function xvalDapc from the R package adegenet was used to select the correct number of principal components with 1,000 replicates using a training set of 90% of the data. The number of principal components was chosen based on the criteria that it had to produce the highest average percentage of successful reassignment and lowest root mean squared error (Jombart et al. 2010). Significant deviations from random population structure were tested in AMOVA utilizing the function randtest from the R package ade4 with 9,999 bootstrap replicates (Dray and Dufour 2007).

    All data and R scripts to reproduce the analyses shown here are deposited publicly on github (https://github.com/grunwaldlab/Sudden_Oak_Death_in_Oregon_Forests) and citable (Kamvar et al. 2014a).

    RESULTS

    Demographic pattern and genetic diversity.

    The epidemic has expanded over time from the initial focus in Joe Hall Creek northeast of Brookings, OR, mostly north (first to N Fork Chetco High) and northwest (Coast, Pistol River South Fork) but also east (Chetco Main and Winchuck) (Fig. 1; Table 1). To date, 70 MLGs have been found in forest populations (Table 1). MLG 22 is most abundant and the only MLG detected across the whole period (although it was not sampled in every year) (Fig. 3). MLG 59, the second most abundant MLG, was only detected in 2011 and has a high frequency due to the sampling design applied: all 2011 strains were sampled in one concentrated area in the northwestern sampling range, geographically distant from any other location (Fig. 1). Given that sampling strategies for some years were not comprehensive, samples from some years have to be interpreted with caution (e.g., 2005 to 2006 and 2010 to 2011). Samples from 2013 and 2014 are sampled from all regions and can be considered more representative.

    Fig. 3.

    Fig. 3. Rank distribution of multilocus genotypes (MLGs) of Phytophthora ramorum and recovery per year. The vertical axis denotes unique MLGs detected in the whole data set with decreasing abundance, as indicated by the bar plot on the right side. The horizontal axis indicates year of sampling. Each numbered circle represents the number of observations of each MLG, with lines connecting genotypes found in multiple years.

    Download as PowerPoint

    Allelic and genotype diversity within loci revealed that PrMS43 had, on average, the highest number of alleles (n = 18). All other loci had five or fewer alleles with a moderate to high amount of diversity (Supplementary Table S2). Nevertheless, the genotype accumulation curve showed a slight plateau, indicating that we have enough power in our data to describe a significant number of MLGs (Supplementary Fig. S2). Genotypic diversity (H = 2.98, G = 8.64), evenness (E5 = 0.41), and richness (eMLG = 7) were low, as expected for a clonal population slowly accumulating mutations over space and time (Supplementary Table S3). A pattern of increasing diversity across years (with number of MLGs not fewer than 10) was also observed. The minimum spanning network showed that MLGs 17, 22, and 28 clustered in the center of the network and had the highest number of connections to other genotypes in the forest populations (Fig. 2). Most genotypes were connected to their immediate neighbors by a genetic distance of 0.05 or the equivalent of one mutational step across five diploid loci.

    Spatial correlation.

    A Mantel test revealed significant correlations of genetic distance and geographic distance for most samples collected after 2003 (Table 3). When partitioned by year, this correlation appears to increase and become more pronounced with the progression of the epidemic. When partitioned by region, those that are closer to the origin of the epidemic (Joe Hall Creek and N Fork Chetco River) show significant correlation. When the overall mantel test was run without Hunter Creek, the correlation coefficient was reduced (0.175) but still significant (P = 0.0001).

    TABLE 3. Table of correlation coefficients generated across forest regions and years of Phytophthora ramorum isolates using the Mantel testa

    Population differentiation.

    Cluster analysis of populations with respect to year using Nei’s genetic distance showed no significant (>70%) bootstrap support for any clades but does show that these tend to cluster by region as opposed to year (Supplementary Fig. S3). AMOVA revealed significant population structure between regions on both clone-corrected (with respect to hierarchy) and uncorrected data sets (Table 4). Significant structure was found only between years within regions on the uncorrected data set. Both patterns were observed without Hunter Creek and Pistol River South Fork isolates. DAPC clustering showed that the first discriminant component separated Hunter Creek from all other regions and the second discriminant component shows a gradient from Joe Hall Creek to the coast (Fig. 4). This distinction was reflected in the percent correct posterior assignment of isolates to their original populations. Over the whole data set, there was an 81.5% assignment-success rate. Hunter Creek received 100% successful reassignment. Joe Hall Creek, Winchuck, and Coast all had >85% successful reassignment whereas Chetco Main, North Fork of the Chetco, and Pistol River South Fork all had <69% successful reassignment (Supplementary Fig. S4). The isolation of the Hunter Creek isolates in the DAPC analysis was found to be mainly driven by allele 493 at locus PrMS43 (Supplementary Fig. S5). The only other population to share this allele was Joe Hall Creek, where it was present in a total of four isolates; and only isolates found in the coastal region or North Fork Chetco contained the allele 489, which is one mutational step away in a stepwise mutation model of a tetranucleotide repeat locus. When DAPC was run without Hunter Creek and Pistol River South Fork data, percent successful reassignment for all regions did not change significantly. Prediction of sources for the Hunter Creek data revealed that over 98% of the genotypes were assigned to the Coast with a 99% probability.

    TABLE 4. Analysis of molecular variance table generated comparing Phytophthora ramorum isolates for two different hierarchies: year within region and region within yeara

    Fig. 4.

    Fig. 4. Scatterplot from discriminant analysis of principal components (DAPC) of the first two principal components discriminating Phytophthora ramorum populations by regions. Points represent individual observations. Lines and shapes represent population membership. Inertia ellipses represent an analog of a 67% confidence interval based on a bivariate normal distribution.

    Download as PowerPoint

    Clustering of forest with nursery populations.

    We used previously published data to determine whether nursery populations in California or Oregon could have been source populations for the Oregon forest epidemic (Goss et al. 2009, 2011; Grünwald et al. 2009; Prospero et al. 2007, 2009). Nursery data included 40 MLGs across 216 samples of NA1 clones. Of these 40, 12 MLGs matched the forest sample and 28 MLGs were unique to the nurseries. The only region that did not contain genotypes that matched those found in nurseries was Pistol River South Fork. When considering those 12 genotypes that were present in both data sets, with the exception of Joe Hall Creek, all genotypes were first isolated from nurseries before discovery in the forest. At the most variable locus, PrMS43, both nursery populations had allele 281 at frequencies of 4.5 and 4.9% for California and Oregon, respectively. This allele was not observed in the forest population. Both populations contained allele 489 at >10% frequency and the California nursery population contained allele 493 at a frequency of 1.4%.

    When nursery genotypes were added to the minimum spanning network, MLGs found at Hunter Creek, previously isolated in the network, connected by only a single MLG from the coast, gained more connections to nursery MLGs. Clustering with Nei’s distance revealed that the nursery isolates from California consistently clustering closest with Hunter Creek isolates in both clone-corrected and uncorrected data sets (Fig. 5). DAPC clustering revealed a decrease in overall assignment-success rate at 78%. The nursery isolates received 74 and 83% assignment success for California and Oregon nurseries, respectively.

    Fig. 5.

    Fig. 5. Unrooted neighbor-joining tree with 10,000 bootstrap replicates of Nei’s genetic distance for Phytophthora ramorum populations defined by region. Tip labels are shaded by region. Branches with bootstrap values greater than 50% are shown. Nursery populations are shown as originating from California (CA) or Oregon (OR).

    Download as PowerPoint

    The assignment successes for the regions before or after inclusion of nursery data changed less than 5% for all regions, except for the coast, which saw a decrease of 17.6% when nursery populations were included. Prediction of sources for nursery genotypes against the forest data revealed that 48% of these nursery isolates were predicted to share membership with the coast at ≥95% probability (Supplementary Figs. S6 and S7). In total, 3.2% of the isolates were predicted to share membership with Hunter Creek at ≥99.9% membership probability. Furthermore, 21.75% of the nursery data could not be assigned to any of the forest populations at >60% probability.

    Because 3.2% of the nursery isolates had a very strong signal for Hunter Creek, we predicted sources for Hunter Creek isolates when considering nursery isolates. This approach determines whether Hunter Creek isolates cluster more readily with nursery or coast populations. Indeed, 92% of the Hunter Creek isolates were predicted to share membership with California nurseries at ≥99% membership probability (Supplementary Fig. S8). No Hunter Creek isolate was predicted to share membership with a forest population at >0.45% membership probability.

    DISCUSSION

    To date, populations monitored from 2001 to 2014 show the presence of only the NA1 clonal lineage observed previously (Prospero et al. 2007, 2009). The fact that individuals belonging to the EU1 and NA2 clonal lineages have not been found in Oregon forests, despite their presence on the west coast from British Columbia to California, is welcome news (Goss et al. 2009, 2011; Grünwald et al. 2012). The lack of EU1 or NA2 isolates provides evidence that monitoring for P. ramorum in nurseries by federal and state agencies is helping avoid emergence of new clones in Oregon’s forests.

    Our analysis provides support for a most parsimonious scenario of two introductions into Curry County from nurseries: one initial introduction into Curry County some time before detection of the first infected tanoak trees in 2001 from California (or, less likely, Oregon) nurseries followed by a second introduction into the Hunter Creek area, again from nurseries. The relative position of the nursery populations in the minimum spanning network and DAPC scatter plot (Figs. 2 and 4) suggest that the introductions from nurseries were rare, though more even sampling and migration models could disprove this hypothesis. Since 2001, the epidemic has spread clonally throughout southwestern Curry County mostly north but also west, toward the coast, and southeast. This clonal spread of the pathogen from the Joe Hall area is supported partially by Mantel tests showing significant levels of isolation by distance in years following 2002 (Table 3) along with significant AMOVA results across regions (Table 4). The populations sampled in 2011 in Hunter Creek (Cape Sebastian) appear to have originated from a new source and cluster into a distinct group based on DAPC (Fig. 4). Based on the minimum spanning network, this population would appear isolated in the epidemic if it were not for MLG 32, which is connected with MLG 33 (found on the coast in 2010) by one mutational step at locus PrMS43. This, in turn, is connected with the other MLGs from Hunter Creek by one mutational step at locus PrMS39. However, when considering clustering via Bruvo’s distance in combination with data from nursery populations, these genotypes from Hunter Creek appear to be more similar to California nursery populations (Fig. 2) than Oregon forest populations. Predictions based on DAPC place samples from Hunter Creek as coming from California nurseries. This, in combination with the observation that purely forest genotypes (i.e., those only found in the forest) are connected to Hunter Creek genotypes through nursery genotypes indicates a possible contribution from nursery populations to the epidemic. This is supported by the observation that all population level clustering, with and without clone correction, places the Hunter Creek isolates adjacent to the nursery isolates from California (Figs. 4 and 5). This appears to be driven by the high frequency of allele 246 at locus PrMS39 which, interestingly, appears to segregate in an east-to-west fashion and is increasing in frequency over time (Supplementary Figs. S9 and S10). This, along with the results from the DAPC clustering and subsequent prediction, provide weak support for a potential third introduction into the coastal region from nurseries some time after the Hunter Creek introduction event.

    An interesting aspect is the observation that there appeared to be more than one cluster of genotypes introduced into the Joe Hall area during the early stages of the epidemic. The two dominant clusters that appeared were the ones that contained MLG 22 and MLG 68. The former has been found in the most recent sampling year, whereas the latter has not been observed since 2005 or beyond the Joe Hall area. This latter group was also the most distantly related group overall, more distant than some nursery genotypes. Although it is clear that the eradication effort has not been entirely successful, there is some evidence that it is having an effect because a major genotype cluster has effectively been eradicated, although disappearance of MLGs could also be explained by being less fit than lineages dominating now.

    The Curry County epidemic is, in many ways, different from the epidemic in California. When introduced into California in the mid 1990s, the causal agent of SOD was unknown, which gave it time to clonally expand and diversify because management strategies in natural forest systems were limited (Rizzo et al. 2002). With the foresight of the epidemic in central California, the Oregon Department of Forestry was able to implement a quarantine effort against the import of hosts as soon as the causal agent was known (A. Kanaskie, personal communication). This quarantine, along with aggressive eradication efforts, has affected the spread of P. ramorum (Mascheretti et al. 2008). Drawing conclusions from previous population studies in California and applying them to the Oregon epidemic should be undertaken with great care given the drastically differing management scenarios (Mascheretti et al. 2008, 2009).

    Our work has some inherent drawbacks. Given the cost of aerial surveys and subsequent ground crew work, and the fact that trees are eradicated once found, populations are not hierarchically sampled across all years. The destructive nature of the management approach means that it was not possible to conduct controlled ecological experiments focusing on effects of climate and rainfall on the spread of disease, as was possible in California trials (Eyre et al. 2013). In addition, most of our work only used five microsatellite loci for genotyping. Ideally, more loci should have been used, as was done in other studies (Croucher et al. 2013). Although only five loci were used, clear patterns of population dynamics in space and time emerged and the MLG accumulation curve supported the fact that loci are informative. Finally, the populations genotyped here are clonal and belong to the NA1 clonal lineage. Thus, much of the analytical power provided by population genetic theory does not apply given that basic assumptions would be violated (Grünwald and Goss 2011). Our work uses appropriate methods to infer patterns that are model free yet informative such as spatial clustering. Thus, we believe that this work provides novel and important insights into the P. ramorum population biology in the Siskiyou forest. Our data indicate that there might have been at least two introductions into Oregon forests from nurseries. The nature of the data does not allow inference of directional migrations given the uneven sampling strategy and moderate number of loci used across all years. We are currently exploring genotyping-by-sequencing (GBS) as a method that could provide further detail on how these populations evolved over space and time (Elshire et al. 2011). GBS can provide richer detail by providing codominant single-nucleotide polymorphism data across several thousand loci sampling the whole genome.

    ACKNOWLEDGMENTS

    This work was supported, in part, by United States Department of Agriculture (USDA) Agricultural Research Service (ARS) grant 5358-22000-039-00D, USDA Animal and Plant Health Inspection Service, the USDA-ARS Floriculture Nursery Initiative, the Oregon Department of Agriculture/Oregon Association of Nurseries (ODA-OAN), and the USDA Forest Service Forest Health Monitoring Program.

    LITERATURE CITED