Plant roots release considerable amounts of labile exudates and debris into the soil, which results in intense microbial activity in the rhizosphere soil which surrounds roots, and the selection of communities which are structurally and functionally distinct from the bulk soil (Morgan et al., 2005) The rhizosphere microbiome can have major impacts on plant growth and nutrition, which can be both positive and negative, through complex direct and indirect interactions. Although it is widely appreciated that a broad range of microbial groups can inhabit the rhizosphere, most studies have focussed on bacteria and fungi (Mendes et al., 2013; Philippot et al., 2013). For these groups, a range of factors can determine the specific communities which assemble in the rhizosphere, including plant characteristics such as genotype and age, environmental properties including soil type and climate, and in agricultural soils, management interventions such as crop rotation and fertilization (Bennett et al., 2012). Protists are also components of the rhizosphere microbiome (Mendes et al., 2013), and can also have marked impacts on plant growth through direct and indirect pathways (Bonkowski, 2004). However, they are typically not considered in studies of rhizosphere microbiology, largely because culture independent techniques to profile complex protist communities remain limited (Adl et al., 2014), with the result that there is little understanding of the factors which shape protist communities in the rhizosphere.

There is one protist group in particular—plasmodiophorids—that rarely receive attention in rhizosphere ecology studies yet are well known as plant pathogens and virus vectors. Plasmodiophorids (Rhizaria; Endomyxa; class Phytomyxea, Order Plasmodiophorida) are parasites and symbionts of angiosperms and stramenopiles (Bass et al., 2009; Neuhauser et al., 2014). Plasmodiophorid diversity is much greater than the few known parasitic taxa and their broader role(s) in the rhizosphere are of great interest. Plasmodiophorids form obligate associations with their hosts, which are often green plants, but in some instances, they can also infect other parasites including heterotrophic stramenopiles, e.g., Woronina pythii, which infects Pythium spp. (Dylewski and Miller, 1983; Neuhauser et al., 2014). Plasmodiophorids are the causative agents of economically significant diseases of crops including brassicas, potatoes, and grain crops (e.g., maize, rice, wheat, sorghum). Plasmodiophora brassicae is the commercially most important and best studied plasmodiophorid, causing clubroot disease in cruciferous plants such as oilseed rape (OSR). Clubroot has been shown to result in average crop losses of 10-15% on a global scale (Dixon, 2009, 2014; Hwang et al., 2012). Other plasmodiophorids include Spongospora subterranea, which causes powdery scab of potato and can also vector Potato Mop Top Virus (Beuch et al., 2015; Falloon et al., 2015). Polymxa graminis can infect most graminaceous crops where it does not cause disease symptoms itself, but can transmit several viruses, such as soil-borne wheat mosaic virus (SBWMV), which is considered as one of the most important diseases of winter wheat in the Central and Eastern USA (Kanyuka et al., 2003). Polymyxa betae vectors Beet Necrotic Yellow Vein Virus, which causes sugar beet “rhizomania,” resulting in ca. 10% loss or the world production of sugar beet (Lemaire et al., 1988; Desoignies et al., 2014; Hassanzadeh Davarani et al., 2014; Biancardi and Lewellen, 2016).

As obligate biotrophs, plasmodiophorids require specific, living hosts for the completion of their life cycle and to enable them to successfully reproduce. But beyond these “primary” hosts, in which the full life cycle can be completed, some plasmodiophorids are associated with a variety of alternative hosts, in which often only the sporangial part of the life cycle can be completed, which results in the formation of short lived zoospores. For example, Spongospora subterranea, whose primary hosts belong to the Solanaceae, can also cause small sporangial infections in hosts within the plant families Poaceae, Brassicaceae, Leguminosae and Geraniaceae (Qu and Christ, 2006). Polymyxa graminis, which has primary hosts including most Poaceae, can for example also infect Arabidopsis thaliana as alternative host (Smith et al., 2011). Similarly, Polymyxa betae, which was considered to be a specialist pathogen of Chenopodiaceae (Desoignies et al., 2014), has now been found to also infect graminaceaous hosts such as wheat (Smith et al., 2013).

The complexity of the plasmodiophorid life cycle (Bulman and Neuhauser, 2017), coupled with their small size (3–6 μm), makes them difficult to study. Specimen-independent molecular probing and sequencing of microbial diversity in environmental samples, referred to as eDNA (environmental DNA) offers an alternative perspective on elusive and cryptic microbes to more classical organism-centric studies (Bass et al., 2015). A recent eDNA study (Neuhauser et al., 2014) investigated plasmodiophorid biodiversity by analysing root and soil-associated plasmodiophorids in a range of habitats, including a vineyard, flood plain and glacier forefield. Eighty-one new OTUs were discovered from just 6 locations, significantly adding to the 41 known phytomyxid (combined plasmodiophorid and phagomyxid) lineages. This suggests that many lineages remain uncharacterized, the biological function of which is unknown. Given the importance of plasmodiophorids as crop disease agents and viral disease vectors, understanding the recently demonstrated “expanded” diversity and distribution of the group within agricultural systems is important to more fully understand crop health, parasite load, and organismal interactions.

A novel insight into the structure of plasmodiophorid communities and the factors modulating this can be deduced by objectively partitioning the community into core and satellite taxa (van der Gast et al., 2011). By decomposing the relationship between distribution (number of sample communities that taxa occupy) and mean abundance (across those sample communities), core members of the plasmodiophorid communities can be identified as those that are locally abundant and non-randomly distributed. With the rare satellite taxa defined as those that are typically in low abundance and randomly distributed through sample communities (Hanski, 1982; Magurran and Henderson, 2003). This is particularly relevant given that so little is known about the ecological processes acting on plasmodiophorid communities, including immigration and extinction, and competition and niche partitioning, all of which are currently unknown in plasmodiophorids (Ulrich and Zalewski, 2006).

For this study, we designed plasmodophorid-specific PCR primers to generate 18S rDNA amplicons suitable for targeted amplicon high-throughput sequencing (HTS), to determine the local scale drivers of plasmodiophorid distribution and to compare their community assembly in rhizosphere/root and bulk soil samples. The use of HTS approaches for functional ecology studies of micro-eukaryotes lags behind that for diversity and/or phylogenetically oriented diversity. However, the use of HTS for both types of study confers the advantage of being able to detect a phylogenetically defined set of lineages without biases associated with sampling for, and accurately identifying, small and cryptically differentiated microbes. We demonstrate that HTS methods can provide information about host associations of microbial lineages without prior knowledge of the microbes involved or any assumption of specific host-microbe relationships.

From the 160 samples, 196,196 plasmodiophorid sequences remained after quality screening, with an average sequence read length of 404 bp and an average of 1,226 sequences per sample. These resolved into 28 plasmodiophorid OTUs (plus three discarded plant OTUs), whose phylogenetic position in relation to characterized plasmodiophorids with full 18S rDNA sequence overlap with the OTU amplicons is shown in Figure 1. Gray vertical lines to the right of the branch tip labels indicate the five OTUs that were effectively identical (allowing for low levels of PCR/sequencing error) to sequences in public databases. The other 22 OTUs were all novel, although some are likely to correspond with phylotypes shown in Neuhauser et al. (2014), which can't be directly tested as the amplicon regions do not overlap and therefore preclude phylogenetic comparison. Three OTUs (2, 11, 1) were identical to named plasmodiophorids: Polymyxa graminis, P. betae, and Spongospora nasturtii respectively. The OTUs were distributed across the known phylogenetic range of plasmodiophorids. Ten OTUs were closely related to characterized plasmodiophorids (Neuhauser et al., 2014): Sorosphaerula veronicae (associated with the herb Veronica spp.), S. viticola (vitaceous vines), Polymyxa graminis (Poaceae) and P. betae (Chenopodiaceae/Ameranthaceae), Woronina pythii (parasite of the oomycete Pythium spp.), Spongospora subterranea (Solanaceae), and S. nasturtii (watercress). The other OTUs were not clearly related to known plamodiophorids, and in some cases (OTUs 17, 14, 23, 15, 19, 4, and OTUs 13, 22, 9, 6) formed two diverse clades whose only previously known members were from environmental sequencing studies. The first of these clades is particularly interesting as it has a moderately well supported sister relationship with the clubroot pathogen Plasmodiophora brassicae, with OTUs 4 and 19 the most similar to P. brassicae, but OTUs 14, 15, 17, and 23 are more closely related to sequences from the Baltic Sea (FN690466) and a freshwater lake moss pillar in the Antarctic (AB695525), both high latitude habitats. Furthermore, this clade is apparently absent from the environmental survey in Neuhauser et al. (2014); we hereon refer to this clade as PlasX. The second clade includes the freshwater-derived EU910610 and two previously detected environmental sequences from the Volga floodplain (Neuhauser et al., 2014); which we refer to subsequently as the “PlasY” clade.

Since there were only minor effects of sampling time on community composition (Table S1), plasmodiophorid samples were separated into four distinct habitat types from across the 2007 and 2008 seasons: wheat rhizosphere (n = 28), wheat bulk soil (n = 28), OSR rhizosphere (n = 52), and OSR bulk soil (n = 52). Subsequently, the relative abundance and distributions of the plasmodiophorid OTUs were analyzed within a metacommunity framework for each habitat. Distribution abundance relationships (DARs) were plotted to ascertain whether each habitat metacommunity exhibited a significant positive DAR and therefore represented a coherent metacommunity in each instance (Figure 2A). Consistent with this prediction, for each habitat, the abundance of individual OTUs was significantly correlated with the number of sample communities that they occupied. Next, DARs were objectively partitioned into core and satellite OTU groups by decomposing the overall distribution using the ratio of variance to the mean abundance for each OTU. The variance to mean ratio, or index of dispersion, is an index used to model whether taxa follow a Poisson distribution, falling between the 2.5 and 97% confidence limits of the χ² distribution. The indices of dispersion were plotted against sample occupancy for OTUs in each habitat metacommunity (Figure 2B).

Of the 22 OTUs that comprised the wheat rhizosphere metacommunity, eight were non-randomly distributed and classified as core OTU group members. 14 OTUs were randomly distributed across samples, falling below the 2.5% confidence limit line, and were classified as satellite OTUs. The phylogenetic distribution of core and satellite OTUs, and the samples in which they were detected are shown on Figure 1. Of the 24 OTUs in the wheat soil metacommunity, 9 were core and 15 satellite. Within the OSR rhizosphere and OSR soil metacommunities there were 8 and 9 core OTUs and 16 and 19 satellite OTUs, respectively. Further, the core group OTUs accounted for the majority of relative abundance in each habitat: wheat rhizosphere, 97.8%; wheat soil, 96.7%; OSR rhizosphere, 97.8%; and OSR soil, 96.9%.

Plasmodiophorid diversity between habitats was compared using Fisher's alpha index of diversity (Figure 3). For the whole metacommunities, mean sample diversity was not significantly different between the wheat and OSR rhizosphere (P = 0.28) and wheat and OSR bulk soils (P = 0.18), but significantly different at the P < 0.05 level in all other instances (Figure 3A). These patterns of diversity were also reflected between core OTU groups: wheat and OSR rhizosphere, P = 0.76; and wheat and OSR soil, P = 0.14 (Figure 3B). Within the satellite OTU groups, although diversity was more variable, OSR rhizosphere mean diversity was significantly lower when compared to each of the soil satellite groups; P > 0.05 in all instances (Figure 3C).

Analysis of similarities (ANOSIM) tests demonstrated that the whole, core, and satellite groups compared between habitats were highly dissimilar and significantly divergent from each other (P < 0.0001, in all instances; Figure 4), with the exception of wheat and OSR soil OTU groups which were not significantly different (P < 0.05). Similarity of percentages (SIMPER) analysis of the whole metacommunity was used to identify those OTUs that contributed most to the dissimilarity between the four habitat metacommunities. These OTUs are listed in Table 1, with typically, the abundant core OTUs contributing most to the compositional dissimilarity between metacommunities. The impact of habitat type was assessed further at the individual OTU level using volcano plots, by plotting fold-change in relative mean abundance against significance (P) values from two-sample t-tests of differences in relative abundance for each taxon (Figure 5). Minimal significant impact was observed between wheat and OSR soil and wheat and OSR rhizosphere metacommunities, respectively. While more pronounced significant fold-changes in relative abundance were observed between rhizosphere and soil metacommunities, irrespective of whether planted with wheat or OSR. ANOSIM showed that rotation (i.e., crop frequency over the preceding 3 years) had no effect on the whole rhizosphere plasmodiophorid communities in OSR, but did affect the bulk soil community (Table S1A). However, for wheat there was evidence for effects of rotation on composition of both rhizosphere and bulk soil communities. These effects were mediated through effects on the core community, with no rotational effects on community composition for rhizosphere or bulk soil of either crop in the satellite community (Table S1B,C).

Our group-specific PCR approach coupled with HTS showed that plasmodiophorids are common and diverse in both rhizosphere and bulk agricultural soils. Positive relationships between abundance and distribution have been observed at many spatial scales for taxa when classified into different types of ecological organization (for example, guild or community; Guo et al., 2000). Within the current study, we also observed significant positive DARs (Figure 2A). This indicated that plasmodiophorid OTUs, that were widely distributed throughout each habitat metacommunity were more locally abundant than the taxa with a more restricted distribution. Therefore, as has been observed for other ecological communities, the commonness and rarity of plasmodiophorid taxa within the different metacommunities was found to be related to their occupancy in the local communities (van der Gast et al., 2011), and allowed decomposition of DARs to objectively catergorize core and satellite OTUs (Figure 2B).

Plasmodiophorid assemblage differed strongly between the rhizospheres of OSR and wheat and between all rhizosphere and bulk soil samples. The only non-significant community comparison was between the two bulk soil samples. Therefore the plants clearly exerted a selective force modifying the plasmodiophorid community in rhizosphere/root samples. It is known that the exudates of certain roots trigger germination of resting spores of P. brassicae and then can serve as alternative hosts (Friberg et al., 2006; Rashid et al., 2013), and it is likely that root exudates stimulate the germination of other phytomyxid species in a similar way. Plasmodiophorid diversity was lower in the rhizosphere compared to bulk soil, illustrated by significantly lower Fisher's alpha diversity indices, particularly with respect to core OTUs in both wheat and OSR rhizospheres (Figure 3). Further, Plasmodiophorid community composition was significantly different in this respect (Figure 4). In these cases the rhizosphere is selective: recruiting a subset of the pool of diversity from the bulk soil into the rhizosphere microbiota. This trend has been observed in other microbial taxa in the rhizosphere (Morgan et al., 2005). Even if the change of local biodiversity is simply caused by missed detection of some OTUs due to relatively very high representation of others, or if its caused by an imbalance in copy numbers of the genes, it is indicative of a diversity shift between rhizosphere and bulk soil habitats. Significant shifts in the frequency of OTU occurrence in the rhizosphere implies a positive or negative functional relationship between plasmodiophorids and plants, enabling us to infer a functional relationship between an environmental sequence and its hosts potentially offering a tool for the identification of secondary interactions. In some cases these involved characterized plasmodiophorids, which corroborated or extended existing knowledge (Braselton, 2001). However, the majority of the OTUs defined in the study were phylogenetically distant from characterized species.

Bacterial and fungal communities inhabiting the rhizosphere are shaped by a variety of factors, including plant genotype, rotation and plant age (Hilton et al., 2013; Chaparro et al., 2014). In our study, most variation in plasmodiophorid rhizosphere communities was associated with plant species, with rotation having minor effects on the community in wheat but not OSR. In OSR and wheat crops 15 and 11 OTUs respectively were significantly differently distributed between plasmodiophorid metacommunities in rhizosphere compared to soil habitats (Figures 5A,F). However, the community shifts were not the same in each case. In the OSR rhizosphere all but three of the 15 significantly different OTUs decreased in relative mean % abundance compared to OSR soil OTUs (Figure 5F). The three OTUs with increased abundance in the OSR rhizosphere were OTU 1 (the brassicaceae-associated Spongospora nasturtii), OTU 7 (Spongospora-like), and OTU 17, which groups with the environmental PlasX clade. Of the other 12 OTUs which had decreased abundance in rhizosphere relative to soil, 2, 4, 5, 9, 24, and 26 had the greatest abundance shifts and lowest P-values: OTU 2 is the largely cereal-associated Polymyxa graminis which will only form secondary, alternative infections in OSR roots, OTU 4 is a novel member of PlasX, more deeply-branching than OTU 17, OTU 5 is the oomycete parasite Woronina pythii, OTU 9 groups in the PlasY clade, and OTU 24 and 26 are related to the Plantaginaceae-associated Sorosphaerula veronicae (Figure 1). These may be associated with other angiosperm species (e.g., weeds; Veronica spp. are abundant agricultural weeds) present in the bulk soil area but not in the crop plant rhizospheres, or other organisms (e.g. oomycetes) associated with other plant species.

The sequence of OTU 1 was very similar to that of Spongospora nasturtii (AF310901), differing only in one homopolymer region in the amplicon, and therefore possibly representing the same sequence in this region. However, the traditional taxonomic concept of phytomyxids based on the morphology and arrangement of the resting spores is not always well reflected in molecular phylogenies. Different plasmodophorid species can have very similar 18S rDNA sequences (Neuhauser et al., 2014), so this does not prove that OTU 1 is S. nasturtii, but it is certainly very closely related. S. nasturtii causes crooked root disease of water cress (Nasturtium officinale, a brassica) in the UK (Claxton et al., 1998), and also Belgium, France, and the US (CABI, EPPO, 2011). Therefore we now show that this or a closely related lineage is preferentially associated with OSR (another brassica), and also negatively associated with wheat rhizosphere relative to proximal bulk soil. Overall the increase of OTU1 in the OSR rhizosphere is interesting from a biological point of view as it indicates a wider host range within brassicas including OSR of the crook root parasite S. narsturtii or the existence of a closely related lineage capable of infecting OSR.

Equally interesting are OTUs representing currently uncharacterized lineages detected only via our environmental sequencing, which were also responding positively and negatively to the rhizosphere habitat, indicating that these lineages interact with the respective plant or another closely plant-associated organism. It is notable that in addition to OTU 7, which has no close relative in Figure 1 and OTU 17 in the PlasX clade, the related OTUs 14 and 23 were also positively associated with OSR, being detected only with OSR in this study, in both soil and rhizosphere. Whether or not these interactions are those of a fully compatible host-parasite pathosystem or if they are of an “alternative host” type cannot be answered at this stage. But any form of increased interaction of host and parasite will have an ecological role which quickly can translate into productivity changes in the agricultural context.

The wheat rhizosphere showed different plasmodiophorid associations compared to bulk soil from OSR plots. Of the 11 OTUs which were significantly differently distributed between wheat rhizosphere and proximal bulk soil, only two were relatively more abundant in the rhizosphere: OTU 24 (which decreased from OSR soil to rhizosphere), and the closely related OTU 25, both of which are closely related to Sorosphaerula veronicae but almost certainly not the same species. Molecular phylogenies have shown previously that Sorosphaerula, Polymyxa and Ligniera form a well-supported clade with not as well defined borders between the genera (Neuhauser et al., 2014). The known hosts of S. veronicae include different Veronica spp. (Plantaginaceae). It is worth noting that a species called S. radicalis has been described from the root hairs of different grasses from the UK (Ivimey Cook and Schwartz, 1929) as well as Ligniera pilorum which was reported from Poa spp. root hairs (Karling, 1968). However, no DNA sequences of these species are available It is therefore possible that some of the OTUs found here might correspond to already described species without a validated DNA record.

The nine OTUs less frequently detected in wheat rhizosphere than bulk soil (OTUs 1, 2, 4, 5, 6, 7, 8, 9, 11; Figure 5A) include seven (2, 4, 5, 6, 9, 11) that were also more abundant in OSR bulk soil than OSR rhizosphere. However, the other two (OTUs 1 and 7; Spongospora relatives) were more abundant in OSR rhizosphere than soil, converse to the situation in wheat. Directly comparing wheat with OSR rhizosphere samples (Figure 5B), the Sorophaerula relatives OTUs 24, 25, 26, 28, and the PlasX member OTU 4 were significantly less abundant in OSR, whereas OTUs 1, 3, 7, 18 (Spongospora relatives), and 17 (long-branched PlasX) were significantly more abundant in OSR. Therefore the rhizosphere habitat of both crop types positively and negatively selected plasmodiophorid lineages from the surrounding bulk soil. Although these sets of OTUs overlapped, they were not identical, with different recognized genera being significantly shifted in each case. Further, the direction of OTU abundance shifts differed between crop species. It is important to note that different plasmodiophorid taxa are associated with particular plant hosts with which they from fully compatible interactions. But from all the species mentioned above in relation to wheat and OSR rhizospheres it is known that a number of other hosts can be utilized (at least) for a short time. Our results suggest some strong associations between protist and plants that were not previously recognized, and it will be up to future research to indentify the biological basis of these associations.

Three OTUs were found in soil only and not any rhizosphere samples—OTUs 16 and 20 (both branching basally to the Sorosphaerula clade), and OTU 19 (sister to OTU 4, which showed negative abundance shifts in both wheat and OSR rhizospheres relative to soil). These OTUs are not sufficiently phylogenetically close to any characterized plasmodiophorid to infer their ecological roles. It is possible that they are directly excluded by other, positively rhizosphere-associated plasmodiophorids or other organisms. Another possibility is that they may be symbionts of other plant species present at the site in limited abundance (weeds) or that they are present in the form of resting spores, or that these lineages are symbionts of other organisms occurring throughout the soil which are relatively less abundant in the more specialized rhizosphere communities. The only lineages with significantly different distribution between the OSR and wheat bulk soil samples were OTUs 6 and 12, which both group in clades not known to be associated with higher plants: OTU 6 in PlasY (Figure 1), closely related to an environmental sequence from a freshwater aquifer (and therefore possibly a parasite of an aquatic alga or oomycete), and OTU 12 grouping in a clade with Woronina (a parasite of oomycetes).

In the light of the high plasmodiophorid diversity detected in this study the absence of any OTU identical or similar to Plasmodiophora brassicae is notable. The PCR primers used had no mismatches with available P. brassicae sequences. However, there is no history of clubroot disease at the sites studied, so its absence is not unexpected.

We show for the first time that many plasmodiophorids beyond the five well-studied plasmodiophorid agricultural parasites (P. brassicae, Spongospora subterranea, S. nasturtii, Polymyxa graminis, P. betae) are present in significant numbers in agricultural soils and in the rhizosphere, even if their known primary hosts are absent or rare. Generally, it is assumed that the distribution of plasmodiophorids follows that of their hosts, but the fact that phytomyxids can use alternative hosts (Neuhauser et al., 2014) means that predicting their diversity and distribution is non-trivial. Within the wheat rhizosphere samples Polymxya graminis-like OTUs 2 and 24 were predictably dominant, as wheat is a primary host plant of P. graminis (Table 1). On the other hand, the S. nasturtii-like OTU 1 dominated the OSR rhizospheres. The primary host of S. nasturtii is another brassicaceous plant (Nasturtium spp.), therefore the organism represented by OTU1 may interact in a similarly compatible way with a range of brassicas. Spongospora sp. are known vectors for plant viruses (Merz and Falloon, 2009), pointing to an additional interesting aspect of this interaction.

The DNA-based detection used in this study does not itself discriminate between active and dormant forms. However, in this system it is apparently sensitive enough to strongly indicate shifts in interaction dynamics between host and symbionts. We show that functional as well as phylogenetic and distribution information can be inferred from environmental sequencing (eDNA) methods combined with a structured and biologically informed sampling strategy. Our results show that diverse plasmodiophorid lineages were positively associated with rhizosphere/root samples compared to bulk soil, and that the enriched lineages were different in wheat and OSR rhizosphere/roots indicating that selection processes in the rhizosphere/root play a role in the establishment and persistence of plant-associated phytomyxids with the potential to increase or decrease the load of pathogenic species.

DNA sequences: GenBank accessions KX263011- KX263038.

DB, ST, and GB designed the research, ST and SH performed the research, DB, CvdG, and GB performed the analyses, DB, CvdG, SN, and GB wrote the paper.