Population Genetic Structures of Puccinia striiformis f. sp. tritici in the Gansu-Ningxia Region and Hubei Province, China

Wheat stripe rust, caused by the fungal pathogen Puccinia striiformis f. sp. tritici (Pst), is a destructive wheat disease in China. The Gansu–Ningxia region (GN) is a key area for pathogen over-summering in China, and northwestern Hubei (HB) is an important region for pathogen over-wintering, serving as a source of inoculum in spring epidemic regions. The spatiotemporal population genetic structure of Pst in HB and the pathogen population exchanges between GN and HB are important for estimating the risk of interregional epidemics. Here, 567 isolates from GN and HB were sampled from fall 2016 to spring 2018 and were genotyped using simple sequence repeat markers. The genotypic and genetic diversity of Pst subpopulations in HB varied among seasons and locations. Greater genetic diversification levels were found in the spring compared with fall populations using principal coordinate analysis and Bayesian assignments. In total, there were 17 common genotypes among the 208 determined, as shown by a small overlap of genotypes in the principal coordinate analysis and dissimilar Bayesian assignments in both regions, which revealed the limited genotype exchange between the populations of GN and HB.


Introduction
Wheat stripe rust, caused by the basidiomycete fungus Puccinia striiformis Westend. f. sp. tritici Erikss. (Pst), is a destructive disease worldwide [1][2][3]. China is the largest wheat-producing and wheat-consuming country, and wheat stripe rust annually causes various levels of yield loss [4][5][6]. In 2017, outbreaks occurred in 18 provinces (containing the main wheat-producing areas of China in the northwest, southwest, and east), and the occurrence area reached up to 5.56 million hm 2 [7].
Understanding the relationships among pathogen population structures from different regions is critical to designing a country-wide disease management strategy. Clearly, the movement of pathogens between over-summering and over-wintering regions is a key issue of interregional epidemics. Such movements promote the completion of the pathogen's annual cycle on a national scale, allowing it to maintain itself over generations. The inoculum can significantly accumulate in over-wintering regions and serve as the source initiating interregional epidemics [5].
The wheat stripe rust pathogen in China cannot over-summer in most of the eastern wheat-growing regions [5]. The mountainous area in the northwestern region includes Southern and Eastern Gansu, Eastern Qinghai, Northwestern Sichuan, and Southern Ningxia ( Figure 1) [8]. The pathogen can over-summer on volunteer seedlings and latematuring wheat at high elevations. The over-summered populations can migrate to low elevations to infect fall seedlings in the same regions, allowing the pathogen to accumulate during fall and winter. Such an inoculum can also disperse southward and eastward to infect fall seedlings, allowing it to survive winter under favorable conditions [5]. The main regions in which the pathogen can over-winter include the Sichuan Basin, Yunnan, Guizhou, Northwestern Hubei (HB), and Southern Henan [5,8]. In these regions, a large number of urediniospores are propagated continuously in winter owing to relatively warm temperatures, and this confirms the contribution of the over-wintering regions to countrywide spring epidemics. Thus, the quantity of urediniospores may greatly influence the intensity of spring epidemics in the eastern regions, including Shaanxi, Henan, Shandong, and other provinces ( Figure 1) [8,9]. Gansu and Ningxia Provinces are the most important pathogen over-summering regions [6]. Wheat is planted from late September in the lowlands to late October at various elevations in the mountains. The pathogen is exchanged frequently among subregions in Gansu and Ningxia Provinces [10]. These two regions may be considered a joint area, providing over-summered inocula for other regions through long-distance dispersal, resulting in interregional epidemics [9].
In addition to the Sichuan Basin, the northwestern part of HB has become another important pathogen over-wintering region in China in recent years. HB Province is adjacent to the southeast of Gansu and close to eastern Sichuan Province ( Figure 1). Wheat is planted from late September to mid-October, and stripe rust usually occurs in late November, forming outbreaks in late February in the northwestern area of HB [11]. Northwestern HB, in which Pst can over-winter [9], may serve as an important geographic area linking the main pathogen over-summering regions (northwestern China) to the major spring epidemic regions. Additionally, northwestern HB represents the rough northern borderline from which Pst can over-winter [8] and accumulate as the initial inoculum of epidemics [12]. This feature may significantly affect the country-wide spring wheat stripe rust epidemics [4]. Many questions regarding the role of the over-wintering areas on interregional epidemics still remain. The population genetic structure of Pst in HB is unclear. Previous studies have shown that CYR34 has become the main race in Gansu, Ningxia, Sichuan, Hubei [13], etc., which indicates a certain spread of the pathogen among different regions. It has been proved that the Sichuan Basin could receive the inoculum from Gansu Province in the fall [13]; the genetic relationship between the inocula in the main over-summering regions and HB is still unclear.
Phylogeographical analyses were helpful to study the migration pathway or direction of pathogen dispersal, and molecular markers could greatly improve our ability to infer the dispersal routes of pathogens [3,10]. Regarding molecular markers, although singlenucleotide polymorphism (SNP) has a high level of genetic variation and has been used in studies of many organisms, such as developing SNP markers with the technology of kompetitive allele-specific PCR (KASP), it is still at the initial stage for Pst [14]. Simple sequence repeats (SSRs) have been the most commonly used marker in recent years [3,10,15] because of their character of codominant and high polymorphism.
In this study, we compared the diversity levels and genetic structures of Pst in the Gansu-Ningxia region (GN) in two fall seasons with those in HB in the subsequent spring seasons to infer the possible migration from GN to HB. The objectives of this study were to (i) determine the genetic structure of Pst in different seasons in HB and (ii) identify possible Pst migration events between the GN and HB populations using various phylogenetic analysis methods.

Sampling and Reproduction of the Isolates
Samplings were conducted in fall 2016, spring 2017, fall 2017, and spring 2018. The sampling periods were from September to January and from February to May. In fall 2016, diseased wheat leaves showing Pst sporulation were collected from fields ( Figure 1), covering three districts having six counties in Gansu Province, one district having four counties in Ningxia Province, and three districts having seven counties in HB Province (Table 1). In spring 2017, samples were collected from two districts having three counties in HB. In fall 2017, diseased leaves were sampled from three districts having six counties in Gansu and one district having one county in Ningxia. In spring 2018, samples were collected from three districts having seven counties in HB. In total, 567 isolates were obtained: 249, 96, 113, and 109 collected in fall 2016, spring 2017, fall 2017, and spring 2018, respectively (Table 1).
The pathogens from the samples were propagated on seedlings of the susceptible wheat cultivar 'Mingxian 169'. Approximately eight seedlings were used per 1 L pot. After a sampled leaf was moisturized at 4 • C for 12 h, urediniospores were picked up from a single uredium pustule and rubbed into a drop of water on the leaf surface of the seedlings as the inoculant. The inoculated seedlings were incubated at 8-10 • C with dew for 24 h in the dark to promote infection and then transferred to an incubator at a 17 • C 14 h day:14 • C 10 h night photo cycle. After approximately 13 days, urediniospores were harvested from the lesions. Each isolate was propagated one or two times, as described above, to obtain an adequate amount of urediniospores. The urediniospores were transferred into a 0.5 mL Eppendorf tube, dried in desiccators at 4 • C for 3-4 days, and stored at −20 • C for DNA extraction.

DNA Extraction and Simple Sequence Repeat (SSR) Genotyping
Genomic DNAs of the Pst urediniospores were extracted with the cetyltrimethylammonium bromide (CTAB) method [16]. In total, 567 DNA samples were extracted. UV absorptions at wavelengths of 260 and 280 nm were measured using a Nanodrop 2000 spectrophotometer (Gene Company Limited, Hong Kong, China) to determine the DNA concentration of each sample. The samples were then diluted to 50 ng DNA/µL for amplification. Three isolates representing each subpopulation were selected for the amplification of the internal transcribed spacer (ITS) region of the ribosomal DNA with the primers ITS1 (5 -TCCGTAGGTGAACCTGCGG-3 ) and ITS4 (5 -TCCTCCGCTTATTGATATGC-3 ) and sequenced. All isolates were identified as Puccinia striiformis f. sp. tritici based on the ITS sequences and morphology of the pathogen [9]. In total, 12 pairs of published SSR or microsatellite primers were used to genotype each isolate ( Table 2). All the SSR primers were labeled with FAM, ROX, TAMRA, or HEX fluorescence at their 5 ends ( Table 2). The PCRs for primers CPS08, CPS13, CPS27, CPS34, and WSR44 consisted of 1 µL of 10× reaction buffer, 0.

Data Analysis
GeneMarker 2.2 [21] was used to analyze the allelic configuration of each sample. A genotype accumulation curve was generated to determine the minimum number of loci necessary to discriminate between individuals in the total population using the poppr package in R [22]. The population genotypic diversity was calculated as corrected lambda, which concerned the differences in the sample size among subpopulations using the poppr package in R [22]. Gene diversity was assessed using Nei's unbiased gene diversity [23]. The subsequent calculations were performed on the clone-corrected data (i.e., every different genotype was represented only once in a population) using the GenClone 2.0 program [24], which reduced the effects of over-representation by certain clones owing to the epidemics. Such over-representations may influence genetic analyses. The genetic diversity within a population was assessed using POPGENE version 1.3.1 [25]. A mean Shannon's information index of each locus was calculated using the Nei method [26], and the standard error was calculated by dividing the standard deviation output produced by POPGENE version 1.3.1 by the square root of 12 [24]. Shared genotypes among subpopulations in the same season or between successive seasons were calculated using the GenClone 2.0 program [24].
A nonparametric principal coordinate analysis (PCoA) based on Nei's distances between all the pairs of SSR genotypes was generated using GenALEx 6.5 [27,28]. The distribution of the genotypes among the studied populations was plotted in two dimensions. A model-based Bayesian method implemented in STRUCTURE 2.3.4 [29] was used to identify genetic clusters and to evaluate the extent of admixtures among them. A total of 267 isolates that were clone-corrected from 567 isolates were used in the STRUCTURE software-based analysis. A model allowing admixture and independent allele frequencies among populations was used. All of the Pst isolates were assigned into K clusters ranging from 1 to 10 on the basis of their clone-corrected SSR multi-locus genotypes. For each simulated cluster K, 50 independent runs were repeated, with 100,000 iterations in Monte Carlo Markov chain replications and a burn-in period of 20,000. The best K value was obtained by STRUCTURE HARVESTER (http://taylor0.biology.ucla.edu/structureHarvester) and evaluated as ∆K, as proposed by Evanno et al. [30]. Indi and pop files of 50 runs corresponding to the best K were processed using CLUMPP 1.1.2 [31], and then, a bar plot was generated by Distruct 1.1 [32] using the output of CLUMPP 1.1.2 [31]. The bar plot analysis results were graphed as pie charts for each location on the basis of the membership fraction Q (Q > 0.8) [29]. Thus, individuals with Q > 0.8 in a given cluster were assigned to the cluster, whereas individuals with 0.2 < Q < 0.8 were assigned to the admixture cluster. Genetic distances between isolates was calculated using GenALEx 6.5, and the phylogenetic trees of the 267 isolates were reconstructed by an unweighted pair-group method with the arithmetic means (UPGMA) approach using MEGA 7. Linkage disequilibrium (LD) tests were performed using Multilocus version 3.1 [33]. The significance of the LD was assessed by the standardized index of association (rBarD). The rBarD values ranged from 0 to 1, with 0 indicating a random association of alleles at different loci, namely linkage equilibrium, and 1 indicating a complete association of alleles at different loci (complete LD) [33]. A randomization test, executed 1000 times, was carried out using the hypothesis of rBarD = 0. The statistical significance of the observed rBarD value was assessed by comparison with the 1000 rBarD value generated under the null hypothesis.

Genotypic Diversity
Among the 567 Pst isolates, 208 unique genotypes were identified. During the wheatgrowing year of 2016-2017, the eMLG of 16F-HB-SY, 16F-HB-XY, and 16F-HB-SZ was much lower (eMLG < 5) than that of 16F-GN (eMLG = 9); however, it increased in the spring of 2017 according to the results of the rarefaction curve ( Figure S1, Table 3). In the next wheat-growing year, the level of eMLG among all subpopulations was similar. The same distribution of genotypic diversity was observed in the corrected lambda (Table 3).  N-1). e Multilocus linkage disequilibrium (LD) was assessed using the standardized index of association (rBarD) and estimated for each population, as described previously [33]. The significance of rBarD was tested with 1000 randomizations of the data by comparing the observed value to the expected value under the null hypothesis of rBarD = 0. The null hypothesis of LD was rejected if p < 0.01. ** Significant at p = 0.01.
The level of genotype exchange could be determined by assessing the shared genotypes among subpopulations. In accordance with the relationship between the two successive seasons in HB, only 4 (g10, g17, g19, and g20) of 69 genotypes (Supplementary file 1: Table S1) were shared, and g20 was the dominant genotype ( Figure 2c). Of 208 genotypes, 10 were shared between the GN and HB subpopulations in the 2016-2017 season (Figure 2a

Genetic Diversity
Genetic diversity was estimated using the parameter of Nei's unbiased gene diversity, and the values in subpopulations in 16F-GN and 17F-GN were similar. Nevertheless, the values in 16F-HB were much lower than those in 16F-GN and increased in the spring of 2017. Except for 18S-HB-JM, Nei's unbiased gene diversities of 18S-HB were higher than those of 17F-GN (Table 3).

Population Subdivision on Spatial and Temporal Scales
PCoA showed that the horizontal and vertical coordinates represented 21.72% and 13.69%, respectively, of the total genetic variance in the 2016-2017 season and 19.54% and 16.31%, respectively, in the 2017-2018 season (Figure 3). The 16F-GN sample was located in all four quartiles, whereas the 16F-HB-SY, 16F-HB-XY, and 16F-HB-SZ samples were mostly located in the second and third quartiles (Figure 3a). This result indicated a certain separation of the HB population from the GN subpopulation in 16F. However, the genotypes in 17S-and 18S-HB were mixed with parts of genotypes in 16F-GN and 17F-GN, respectively (Figure 3a,b), suggesting that a small level of population exchange exists between the fall subpopulations in GN and the spring subpopulations in HB. The result of the Bayesian assignment analysis carried out using STRUCTURE revealed three genetic clusters among the 267 isolates ( Figure 4). Individuals in the GN subpopulation were assigned into the G1, G2, G3, and admixture groups, indicating a certain degree of diversification (Supplementary file 2: Table S2, Figure 4a). Three subpopulations in HB in 2016 were mainly assigned into the G3 and admixture groups, revealing a large difference compared to the GN subpopulation in the same season (Figure 4c). Although a number of G1 and G2 isolates appeared in 17S-HB, the composition was more similar to that of 16F-HB than that of 16F-GN. Similarly, in 2017-2018, the proportion of G2 isolates was much lower in 18S-HB than in the previous season in GN. The phylogenetic analysis showed similar results with PCoA and Bayesian assignment analysis. The isolates of 16F-HB were mainly clustered with 17S-HB, and a small number of isolates in 16F-GN had closed genetic relationships with 17S-HB (Figure 5a). Similarly, some isolates of 17F-GN were mixed with 18S-HB (Figure 5b).   (Table 3).

Discussion
The genetic diversity levels and population structures of seasonal subpopulations of Pst in HB were investigated for the first time in this study. Owing to the severe disease epidemics that occurred from 2016 to 2017 in China under the influence of high winter temperatures [7,34], a large number of diseased leaves were collected in fall 2016, which was not common. Previous field investigations indicated that there was no possibility for Pst to over-summer in northwestern HB owing to the high temperatures and the absence of Pst hosts [12,35]. Thus, the primary inoculum in fall seasons most likely migrated from the over-summering regions. Therefore, the occurrence of wheat stripe rust disease in HB relied on external sources of Pst.
There were only 17 of 208 shared genotypes, a limited PCoA overlap, and dissimilar Bayesian assignments, which revealed that the Pst populations of GN and HB underwent a low level of pathogen exchange. Pst migration might occur between GN and HB after over-summering; however, the amount of inoculum in HB received from GN was limited. Therefore, the disease intensity in fall seedlings of GN was not an indicator of a possible risk of wheat stripe rust epidemics in HB in the spring.
Shared genotypes, PCoAs, Bayesian assignments, and phylogenetic analysis help infer the levels of population exchange. Liang et al. [10] found that 40 of 55 Pst genotypes are common between Gansu and Ningxia Provinces, indicating that a frequent pathogen exchange occurs between the two regions. Similarly, Wan et al. [16] found that 7 of 72 Pst genotypes are shared between Gansu and Xinjiang subpopulations, illustrating that Xinjiang and Gansu are not completely isolated. Similarly, in Liang et al. [15], a large number of shared Pst genotypes, large PCoA overlap, and homologous Bayesian assignments confirmed that migration occurs between the fall population in Gansu and the spring population in the Sichuan Basin. In this study, 17 of 208 genotypes were shared between the GN and HB subpopulations; therefore, these two regions are not isolated completely.
The amount of Pst in HB in the early spring could influence the occurrence of wheat stripe rust in further eastern regions because of transfer by western winds [5,36,37]. Therefore, it was crucial to clearly define the origin source of Pst in HB after over-summering. A previous study illustrated that Pst spores can disperse to HB from the northwest (containing GN) through an analysis of the upper air current above 5500 m [38]. However, in this study, apart from the shared genotypes between the subpopulations in HB and GN, there were 44 genotypes in HB in S17 and 30 genotypes in S18 (Supplementary file 1: Table S1) that were not detected in previous subpopulations in GN. That implied that the main source of Pst came from other regions to HB after over-summering, because the pathogen cannot over-summer locally.
Field investigations confirmed that Pst cannot over-summer in northwestern HB because of the high temperature and the absence of wheat in the high-altitude area [31]. Therefore, the inoculum in northwestern HB after summer was most likely from the northwest. In addition to Gansu, Ningxia, Shaanxi, Qinghai, Xinjiang, and Tibet also suffer from severe stripe rust disease ( Figure 1) [16,39]. However, previous studies have illustrated that there is limited pathogen exchange between Xinjiang, Tibet, and the other regions of China [16,39]. Pst found east of Qinghai could over-summer; however, spring epidemics in Qinghai rely on inocula from other regions because Pst cannot survive the winter [40]. Thus, the migration from fall to winter inferred in this study might not have originated in Qinghai. This study illustrated that the frequency of pathogen exchange between GN and HB is low. Thus, Pst populations in the middle area (such as western Shaanxi, Figure 1) could migrate to northwestern HB after the summer. Still, the possible migration from the southeast to HB could not be ruled out. Although Pst could survive year-round in some locations of Yunnan and Guizhou [9], the detailed locations of oversummering and over-wintering areas in Yunnan and Guizhou and the genetic relationships of Pst with subpopulations in other epidemic regions are still not clear. Therefore, future studies should focus on the population genetic structure of Pst among populations in western Shaanxi, northwestern HB, Yunnan, and Guizhou to determine the main Pst inoculum source in northwestern HB after over-summering.
This study revealed the usefulness of phylogenetic approaches in inferring the genetic relationships of wheat stripe rust pathogen among subpopulations. More sampling sites in HB would help reveal distinct variations compared with the GN populations. Future studies of the long-distance dispersal of Pst would benefit from next-generation sequencing (e.g., restriction-site-associated DNA sequencing), which detects more genetic differences than microsatellites. Therefore, this technology would improve assignments that infer the relationships among regions having a high level of gene flow [41,42]. These data will be vital to help determine the geographic range and timing of disease management strategies in different regions.
Interregional disease management strategies rely on fully understanding the interactions among geographic pathogen populations and their possible long-distance migration patterns. Such information significantly enhances the designs of country-wide disease management strategies. For instance, the efficient reduction of the inoculum in source regions would significantly reduce the inoculum in the target regions resulting from longdistance pathogen dispersal. In this study, because HB is considered a region in which Pst accumulates during the winter, significantly reducing the amount of inoculum produced from this region becomes an important strategy to reduce the risk of spring epidemics in the northern regions to which the exogenous inoculum is transferred.

Conclusions
In conclusion, we found various levels of Pst genetic and genotypic diversity in HB in different seasons. Although several common genotypes existed in GN and HB, the small overlap in the PCoA, and the dissimilar Bayesian assignment, revealed that the Pst populations of GN and HB experience a low level of exchange. This study revealed that Pst migration might occur between GN and HB after over-summering; however, the amount of inoculum in HB received from GN was limited. We recommend further investigations on the influence of Pst sources from Shaanxi and the southwestern regions of China (such as Yunnan and Guizhou) on Pst in HB. They may provide additional clues to increase our understanding of the migratory pathways of Pst.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/genes12111712/s1. Figure S1: Rarefaction curve for the number of genotypes of Puccinia striiformis f. sp. tritici that would be expected at the largest shared sample size for different subpopulations. Supplementary file 1: Table S1: The frequencies of all the Puccinia striiformis f. sp. tritici genotypes generated with 12 pairs of SSR primers. Supplementary file 2: Table S2