Land Use Change and Water Quality Use for Irrigation Alters Drylands Soil Fungal Community in the Mezquital Valley, Mexico

Soil fungal communities provide important ecosystem services, however, some soil borne representatives damage agricultural productivity. Composition under land-use change scenarios, especially in drylands, is rarely studied. Here, the soil fungal community composition and diversity of natural shrubland was analyzed and compared with agricultural systems irrigated with different water quality, namely rain, fresh water, dam-stored, and untreated wastewater. Superficial soil samples were collected during the dry and rainy seasons. Amplicon-based sequencing of the ITS2 region was performed on total DNA extractions and used the amplicon sequence variants to predict specific fungal trophic modes with FUNGuild. Additionally, we screened for potential pathogens of crops and humans and assessed potential risks. Fungal diversity and richness were highest in shrubland and least in the wastewater-irrigated soil. Soil moisture together with soil pH and exchangeable sodium were the strongest drivers of the fungal community. The abundance of saprophytic fungi remained constant among the land use systems, while symbiotic and pathogenic fungi of plants and animals had the lowest abundance in soil irrigated with untreated wastewater. We found lineage-specific adaptations to each land use system: fungal families associated to shrubland, rainfed and part of the freshwater were adapted to drought, hence sensitive to exchangeable sodium content and most of them to N and P content. Taxa associated to freshwater, dam wastewater and untreated wastewater irrigated systems show the opposite trend. Additionally, we identified potentially harmful human pathogens that might be a health risk for the population.


INTRODUCTION
Studies of soil fungal communities are increasing due to their ecological importance. Fungi are the major eukaryotic lineage in terms of biomass in soil, surpassing all other soil organisms combined (excluding plant roots) (Orgiazzi et al., 2012). Saprotrophic fungi participate in critical processes, such as the decomposition and mineralization of both recalcitrant and labile compounds of plant and animal origin, such as cellulose, hemicellulose, lignin, and chitin (Treseder and Lennon, 2015). Fungi are also involved in symbiosis with plant roots [arbuscular mycorrhiza (AMF) and ectomycorrhiza (EcM)], which enable plants to survive in water and nutrient limited environments (Bennett et al., 2017). Fungal communities are also known to include soil-born pathogenic representatives (Tedersoo et al., 2014). Despite their importance in maintaining soil quality and ecosystem functioning, the fungal community composition in soils is less studied than the bacterial part (Tedersoo et al., 2014). The structure and diversity of soil fungal communities are driven by soil physicochemical characteristics, such as moisture (controlled by precipitation or irrigation) (Tedersoo et al., 2014;Zheng et al., 2017), organic C availability (Lauber et al., 2008), and nutrient content (He et al., 2016) as well as the plant community (Yang et al., 2017).
Drylands are regions with arid, semi-arid, or dry sub-humid climate (Koohafkan and Stewart, 2008). This type of ecosystem covers 40% of Eartht's surface and it is likely to expand due to climate change (Önder et al., 2009;Maestre et al., 2015). The agricultural productivity in drylands, although limited by soil moisture, supports one third of the global population (United Nations Environmental Management Group [UNEMG], 2011). As freshwater resources are scarce in these regions, it is common to use treated and untreated wastewater for irrigation to enhance the productivity (Siebe, 1998;United Nations Environmental Management Group [UNEMG], 2011;Norton-Brandão et al., 2013;Becerra-Castro et al., 2015).
The land cover change from shrubland to agriculture in drylands often diminishes soil carbon stocks and modifies the quantity and quality of organic matter, especially if it is converted to maize monoculture (Sánchez-González et al., 2017). The use of wastewater for crop irrigation particularly improves the availability of labile organic carbon and nutrients (N and P) in soil (Siebe, 1998), but it also increases soil carbon stocks due to large root residue inputs (Sánchez-González et al., 2017). However, this practice adds soluble salts, heavy metals, pharmaceuticals, and microbial pathogens to the soil, which might affect the microbial soil communities, as well as crop and human health (Siebe, 1995(Siebe, , 1998Hong and Moorman, 2005;Broszat et al., 2014;Dalkmann et al., 2012).
Soil fungal communities affected by land use change and management have been frequently studied in forest and grassland ecosystems, often contrasting between tilled or fertilized systems (Birkhofer et al., 2012;Klabi et al., 2015;He et al., 2016;Wang et al., 2016). Special attention has been given to arbuscular mycorrhizal fungi (AMF), concluding that the transformation to agriculture reduces its richness and diversity (Moora et al., 2014;Xiang et al., 2014;Verbruggen et al., 2015). The minority of soil fungal studies addresses communities in drylands (Vargas-Gastélum et al., 2015), only a few of them have targeted the effect of the water quantity and quality used for crop irrigation on the soil fungal community composition, diversity and function in dry regions (Ortega-Larrocea et al., 2001;Alguacil et al., 2012;Disciglio et al., 2015;García-Orenes et al., 2015).
In this study, and for the first time, the soil fungal community composition is compared between semiarid shrubland and agricultural land irrigated with different water qualities. We determined the soil fungal community composition and diversity of natural dryland vegetation (shrubland), rainfed system, and irrigation agriculture with different water qualities during the dry and rainy seasons using amplicon-based MiSeq sequencing of the internal transcribed spacer 2 (ITS2) region. DNA-based analysis to assess total fungal community was employed and FUNGuild (Nguyen et al., 2016) was used to evaluate specific fungal trophic modes. The study was conducted in the area known as Mezquital Valley, which is located in central Mexico and covers 90,000 ha. This region represents the largest continuous area in the world in where crop irrigation is performed with untreated wastewater (Siebe and Cifuentes, 1995;Siemens et al., 2008). The soil of this area provides two important services: (a) filtering the untreated wastewater recharging the aquifer, which is used by the surrounding population for consumption and recreational activities, and (b) crop production, mainly fodder (lucerne, oat, rye grass, as well as maize) (Lüneberg et al., 2018). In 2014 more than 4 million tons of fodder crops and vegetables were produced in this area (Conagua and Semarnat, 2015).
We stressed the hypotheses (i) that land use change and different irrigation practices in the Mezquital Valley area affect fungal diversity and community composition. It has been shown that fungal community structure was affected by land use change in different ecosystems (Orgiazzi et al., 2012;Lupatini et al., 2013;McGuire et al., 2015;Brinkmann et al., 2019). Irrigation practices and seasonality modify soil moisture, especially in drylands, which is an important driver of fungal activity (Tedersoo et al., 2014;Vargas-Gastélum et al., 2015), and thus, we hypothesize that (ii) seasonality plays a role in structuring fungal communities, and soil moisture is the strongest driver, followed by nutrient availability. Due to previous studies reporting an increased frequency of saprophytic and plant pathogenic fungi in irrigated and fertilized soil (Calderon et al., 2016;Zheng et al., 2017), and a decrease of AMF in wastewater irrigated soil (Ortega-Larrocea et al., 2001;Alguacil et al., 2012), we assume that (iii) saprophytic and pathogenic fungi will increase and AMF will decrease in soils under periodic irrigation with wastewater. This study will help to better understand the impact of environmental manipulation on fungal communities in drylands, especially of those taxa involved in essential functional aspects of these soil ecosystems.

Study Sites and Sample Collection
The sampling area was described by Lüneberg et al. (2018). Briefly, The Mezquital Valley has a semiarid climate with average annual temperatures from 16 to 18 • C and an average annual rainfall from 400 to 600 mm (British Geological Survey [BGS], 1998). The natural vegetation in the area is classified as xerophytic shrubs with mesquite (Prosopis juliflora) as the dominant tree species. The detailed description of the study area and a map showing the sampling locations can be found in the Supplementary Text S1 and Supplementary Figure S1.
We chose sampling plots with different land use systems and with clay to silty clay loam texture, as it is one of the representative soil types of the area (Haplic Phaeozem; IUSS Working Group WRB, 2014). The sampling encompassed four plots × five land use systems × four samplings. The land use systems included natural dryland vegetation (shrubland plots S1-S4), rainfed plots (R1-R4), plots irrigated with freshwater (FW1-FW4), plots irrigated with wastewater coming from the Endhó dam (DWW1-DWW4) and plots irrigated with untreated wastewater (UTWW1-UTWW4) (Supplementary Figure S1). The size of the plots ranged among one to two ha. The sampling was performed twice during the rainy season (June -October) and twice during the dry season (November -May). All the plots used for agriculture were cropped with maize during the rainy season. During the dry season the rainfed plots were abandoned and the irrigated plots were fallow or cropped with lucerne, oat or grass. At each sampling plot, 20 bulk soil cores from the upper 10 cm were sampled in a regular systematic grid. The cores were homogenized and pooled into one composite sample per plot. The sample was frozen in liquid nitrogen and stored at −80 • C until DNA extraction.

Soil Properties
The soil samples were air-dried for 24 h, homogenized and sieved using a metallic mesh (2 mm). The physical and chemical properties were determined using standard procedures (Soon and Hendershot, 1992). To determine soil pH, 10 g of each soil sample were suspended at a soil-to-liquid ratio of 1:2 (soil − 0.01 M CaCl2). Subsequently, the pH in the supernatant was measured with a glass electrode (Conductronic pH 120, Puebla, Mexico). The gravimetric soil water content (%) was calculated from oven-dried subsamples. Soil organic C content and total N was determined using a Perkin Elmer 2400 CHNS/O elemental analyzer (Massachusetts, United States). Available P was determined by spectrophotometry (Genesis 20, Massachusetts, United States) based on the methodology developed by Olsen (van Reeuwijk, 2002). Soil particle size distribution was determined using the Bouyoucos hydrometer method. The electrical conductivity was determined in a 1:2 soil:distilled water suspension with a conductivity meter (Hanna HL4321, Rhode Island, United States). Exchangeable bases (Ca 2+ , Mg 2+ , Na + and K + ) were extracted with 1N NH 4 OAc, and the determination of Ca 2+ and Mg 2+ was done by atomic absorption spectrophotometry (Perkin Elmer 3110, Massachusetts, United States), for Na + and K + flame emission spectrometry (Sherwood Scientific 36, Cambridge, United Kingdom), according to van Reeuwijk (2002). Table 1 shows the main soil properties. Written and graphic description of soil properties can be found in the Supplementary Text S2 and Supplementary Figure S2.

DNA Extraction and Amplification of ITS2
Soil samples were freeze-dried with liquid nitrogen and kept frozen at −80 • C before nucleic acid extraction. Total DNA was extracted by employing the MoBio PowerSoil DNA isolation kit (MoBio Laboratories, Carlsbad, CA, United States) following the manufacturer instructions. DNA concentrations were quantified by using a NanoDrop ND-1000 Spectrophotometer (Thermo Scientific, Schwerte, Germany). Fungal ribosomal internal transcribed spacer (ITS) amplicons were generated using primers ITS3_KYO2 (GATGAAGAACGYAGYRAA) (Toju et al., 2012) and ITS4 (TCCTCCGCTTATTGATATGC) (White et al., 1990) with an Illumina adapter overhang (forward: TCGTCGGCAGCGTCAGATGTGTATAAGAGACAG, reverse: GTCTCGTGGGCTCG GAGATGTGTATAAGAGACAG). For the PCR reaction mixture contained Phusion High Fidelity DNA Polymerase (Thermo Scientific, Schwerte, Germany), 5% DMSO, 0.2 mM of each primer, 1 U DNA Polymerase and 25 ng of isolated DNA as template. Thermal cycling scheme for ITS2 amplicons was as follows: 1 initial min at 98 • C, 25 cycles of 45 s at 98 • C, 45 s at 48 • C, and 30 s at 72 • C, and a final extension at 72 • C for 5 min. The resulting PCR products were checked by agarose gel electrophoresis for appropriate size and purified by magnetic bead clean-up using the MagSi-NGSPrep Plus as recommended by the manufacturer (MagnaMedics Diagnostics B.V., Geleen, Netherlands). Quantification of the PCR products was performed using the Quant-iT dsDNA HS assay kit and a Qubit fluorometer (Invitrogen GmbH, Karlsruhe, Germany) following the manufacturer's instructions. All PCR reactions were performed in triplicate and pooled in equal amounts. Pooled PCR products were used to attach indices and Illumina sequencing adapters using the Nextera XT Index kit (Illumina, San Diego, CA, United States). Index PCR was performed using 5 µl of template PCR product, 2.5 µl of each index primer, 12.5 µl of 2x KAPA HiFi HotStart ReadyMix and 2.5 µl PCR grade water. Thermal cycling scheme was as follows: 95 • C for 3 min, 8 cycles for 30 s at 95 • C, 30 s at 55 • C and 30 s at 72 • C and a final extension at 72 • C for 5 min. Fungal ITS2 region were sequenced using the dual index paired-end 2x 300 bp approach (v3 chemistry) and the Illumina MiSeq platform as recommended by the manufacturer (Illumina).

Sequence Processing and Analyses
Demultiplexing of raw sequences was performed by CASAVA data analysis software (Illumina). Paired-end sequences were merged using PEAR v0.9.10  with default parameters. The resulting sequences that fulfilled at least one of the following criteria were removed with the split_libraries_fastq.py script from QIIME 1.9.1: average quality score lower than 20 and containing unresolved nucleotides . Reverse and forward primer sequences were removed by employing cutadapt v1.12 (Martin, 2011) with default settings. In preparation for the generation of amplicon sequencing variants (ASV), which represents exact sequence variants (Callahan et al., 2017), USEARCH (9.2.64) with the UNOISE2 (Edgar, 2016) algorithm was used to dereplicate (-fastx_uniques), and removed sequences shorter than 140 bp (-sortbylength). After denoising (-unoise2) chimeric sequences were removed using UCHIME2 (Edgar et al., 2011) included in software package USEARCH (9.2.64) in de novo and reference mode (-uchime2_ref) against UNITE database (v7.1) 1 (Koljalg et al., 2013). All quality-filtered sequences were mapped to chimera-free ASVs and an ASV abundance table was created using USEARCH (-usearch_global) with default settings. Taxonomic classification of the picked reference sequences (ASV) was performed with parallel_assign_taxonomy_blast.py against the UNITE database (v7.1) 1 (Koljalg et al., 2013). Extrinsic domain ASVs and unclassified ASV were removed from the data set by employing filter_otu_table.py. Finally, all unidentified fungal ASVs were searched with BLAST (blastn) (Altschul et al., 1990) against the nt database (version from March 2017) to remove non-fungal ASVs and only as fungi classified reads were kept. Sample comparisons were performed at the same surveying effort, utilizing the lowest number of sequences by random selection (total 16,800). Species richness, alpha and beta diversity estimates and rarefaction curves were determined using the QIIME 1.9.1 script alpha_rarefaction.py.

Prediction of Fungal Trophic Modes
Fungal trophic modes were predicted from the resulting ASV table using the software package FUNGuild (Nguyen et al., 2016). Fungal ASVs assigned only to one trophic mode (symbiotrophy, saprotrophy, and pathotrophy) were included in the analyses, as environmental conditions can modify the guild of specific groups (Nguyen et al., 2016). For pathothrops the abundance 1 https://unite.ut.ee/repository.php of plant and animal pathogens was analyzed separately. The abundances of symbiotrophs of arbuscular mycorrhizal (AM) and ectomycorrhizal (EcM) fungi were also analyzed separately.

Identification of Human Pathogens
The rarefied soil and wastewater ASV tables were screened for genera containing potentially human pathogenic members such as Ajellomyces, Alternaria, Apophysomycetes, Aspergillus, Bipolaris, Blastomyces, Candida, Cephaliophora, Cladosporium, Coccidioides, Colletotrichum, Cryptococcus, Cunninghamella, Curvularia, Fonsecaea, Fusarium, Histoplasma, Lasiodiplodia, Lichtheinia, Malassezia, Microsporum, Mucor, Paracoccidioides, Penicillium, Pneumocystis, Rhizomucor, Rhizopus, Rhodotorula, Scedosporium, Schizophyllum, Sporothrix, and Trichophyton based on the information provided by Mycology Online of the University of Adelaide, Australia 2 . This search was made using QIIME 1.9.1 script filter_taxa_from_otu_table.py. The corresponding sequences were double checked against the nucleotide collection of the NCBI using blastn. ASVs that exhibited > 98% identities to known pathogenic members of the genus were further analyzed. The pathogenic nature of the thereby identified taxa was confirmed by literature searches. Taxa for which no pathogenic features were described were discarded.

Statistical Analyses
All statistical analyses were conducted employing R version 3.3.1 (R Core Team, 2016). The results of all statistical tests were regarded significant with p ≤ 0.05. To evaluate the differences among land use systems and between seasons we employed a generalized mix effect model with the land use as a fix effect, followed by Bonferroni's multiple comparison tests (Shannon index, observed species, and fungal trophic modes). Kruskal-Wallis and Dunn's tests were used to determined differences among land use systems for the most abundant phyla and potential human pathogens. To evaluate the correlation of diversity indexes and fungal trophic modes with soil properties Spearman's rank correlation test was employed. To identify the fungal families associated with the different land use systems, an analysis based on the point biserial correlation coefficient was performed using multipatt (package "indicSpecies") (de Cáceres, 2013). For visualization, a network was generated using the land use systems as source nodes, and the fungal families as target nodes. All taxa with significant (p ≤ 0.05) associations were visualized in the networks. Network generation was performed using the edge-weighted spring embedded layout algorithm in Cytoscape (Shannon et al., 2003), with the edge weight corresponding to the association weight of each family with each land use system. To evaluate and graph the most abundant genera we employ ampvis2 (Andersen et al., 2018), unidentified ASVs were removed, color scale of abundance was square rooted to better visualize low abundant taxa. To assess the (dis)similarity of fungal communities between land use systems and season ANOSIM test was performed in QIIME  using Bray-Curtis dissimilarities with the script compare_categories.py. To visualize the multivariate dispersion of community composition a non-metric multidimensional scaling (NMDS) was performed using Bray-Curtis dissimilarities employing the "vegan" package (Oksanen et al., 2015). Standard deviation ellipses by land use system were projected onto the ordination, utilizing the function ordiellipse. The effects of environmental parameters on the fungal community were analyzed using the envfit function and significant correlated parameters were projected into the ordination with arrows. In order to evaluate the multicollinearity among environmental parameters, a Spearman's rank correlation test was performed and the variance inflation factors (VIF) were determined (Supplementary Table S1). As nitrogen and carbon were highly correlated (R 2 > 7 and the VIF value of nitrogen was > 4), nitrogen was removed from the analysis.

Sequence Data Deposition
The ITS2 region sequences were deposited in the National Centre for Biotechnology Information (NCBI) Sequence Read Archive (SRA) under bioproject accession number PRJNA386070.

General Characteristics of the ITS Dataset
To analyze the fungal community structure, DNA was isolated from 80 soil samples derived from all analyzed land use systems.
Despite several attempts with modified conditions, the isolation of high quality DNA from one rainfed sample (R1D1) failed. After removal of low-quality sequences and singletons, amplicon-based analysis of the ITS2 region of the fungal rRNA resulted in 5,645,373 high-quality sequences. The number of sequences per sample ranged from 16,830 to 342,599. After rarefaction analysis with the minimal amount of sequences per sample (16,800), 9,685 ASVs were obtained (773 + 203 per sample) (Supplementary Figure S3). The Good's coverage index of 0.98 (+ 0.004) indicates that the dataset enclose all major fungal groups inhabiting the studied land use systems.

ASVs Richness and Diversity
Fungal diversity (Shannon index H ) and ASVs richness (number of observed ASVs) responded to land use change and water used for irrigation (p ≤ 0.05, Figure 1). The diversity of the fungal community was largest in the shrubland soil (H 6.8), mainly during the dry season (H 7.4) and smallest in the wastewater-irrigated fields (H 6.0). The same was found for the ASVs richness, which were higher in shrubland soil, in average 956 ± 72, compared to 619 ± 18 from the UTWW system. The shrubland ASV richness was especially higher during the dry season 1123 ± 39 compared to 789 ± 112 found during the rainy season (p ≤ 0.05).
The evaluation of soil properties impact on the fungal diversity and ASVs richness ( Table 2) showed a strong negative correlation to soil moisture and Na + content, and a weaker but significant negative correlation to available P, K + , and Mg 2+ content (p ≤ 0.05). , and Untreated wastewater (UTWW) irrigation systems, during the dry and rainy season. Box are extended from the 25th to 75th percentiles, the line in the box is plotted at the median. Whiskers represent the smallest and the largest value. A mix effect model followed by BonferroniŠs multiple comparison tests were used to determine differences among land use systems, and between seasons. Different letters indicate statistical difference among land use systems and * indicate statistical differences between seasons (p ≤ 0.05). Only parameters differing significantly between seasons in each land use system are shown, if they did not differ, samples of both seasons were merged.

Fungal Community Composition
The composition of fungal communities differed among land use systems as indicated by ANOSIM test (p = 0.001). The (dis)similarities seem to be related to whether the land was untilled or used for agriculture. Shrubland samples were separated from the other land use systems, while the communities of rainfed and FW showed certain similarity (Figure 2), and the communities irrigated with wastewater (DWW and UTWW) were highly akin. Based on the ANOSIM test, none of the studied land use systems differed seasonally, however, during the dry season, the agricultural fields were either abandoned, fallow or cropped with oat, grass or lucerne.
To evaluate for multicollinearity of the soil physico-chemical properties a Spearman correlation test was performed and the VIF was determined. Total C and N showed high correlation (R 2 > 7) and N obtained a VIF > 4. Based on these results we removed N from the analysis. After fitting the soil parameters FIGURE 2 | Non-metric multidimensional scaling (NMDS) of the total fungal community composition from land use systems, shrubland, rainfed, freshwater, dam wastewater and untreated wastewater irrigation and dry and rainy season. Soil parameters with strong and significantly correlation (p ≤ 0.05) to the fungal community structure are indicated by arrows.
onto the ordination, we found that exchangeable sodium content exhibited the strongest effect on the fungal community (R 2 = 0.62, p = 0.001), followed by soil pH and moisture (R 2 = 0.55, R 2 = 0.37, respectively, p = 0.001; Figure 2). Interestingly, soil base cations (Mg 2+ , K + , and specially Na + ) played an important role in the composition of the fungal community, being almost as important as soil humidity, and even more important than the nutritional status of soil (P content and C:N ratio).
To identify families significantly associated with one, two, or more land use systems, a correlation-based association analysis was performed. Here, all fungal families were included, however, on average the relative abundance of 52% of fungal families were not significantly different with respect to the land use system. The correlation-based association analysis (Figure 5) was consistent FIGURE 3 | Relative abundances of soil fungal orders. Land use systems: Shrubland (S), Rainfed (R), Freshwater (FW), Dam wastewater (DWW) and Untreated wastewater (UTWW) irrigated, during dry and rainy season. Bacterial orders with average relative abundances >0.5% is visualized; orders contributing ≤0.5% were summarized as rare taxa. with the multivariate analysis (Figure 2). The fungal community from shrubland soil was separated from the agricultural land use systems, representing the system with the highest number of associated taxa.
Most fungal taxa significantly associated to shrubland soil belong to Ascomycota, the families with stronger association were Annulatascaceae, Chaetothyriaceaea, Verrucariaceae, Geastraceae, Massariaceae, Trichomaceae, Didymosphaeriaceae, Lentitheciaceae, Teloschistaceae, and Parmeliaceae (Figure 5). The latter was also associated to rainfed soil, also associated to rainfed soil were Spizellomycetaceae, Pezizaceae, Hyponectriaceae, Filobasidiaceae, Pyronemataceae, Claroideoglomeraceae, and Kickxellaceae. Families associated to the three land use systems free from wastewater (shrubland, rainfed, and FW) were Dothioraceae and Sporomiaceae. The relative abundance of fungal families associated to shrubland, rainfed and partially FW had a negative correlation to soil moisture, Na + content, pH, and available P (Supplementary Table S2). Some of these taxa are Annulatascaceae, Geastraceae, Massariaceae, Trichocomaceae, Didymosphaeriaceae, and Lentitheciaceae. Most of the taxa associated to rainfed soil were additionally negatively correlated to N content, including, Hyponectriaceae, Spizellomycetaceae, Filobasidiaceae, Claroideoglomeraceae, and Kickxellaceae. Two of the families associated to the land use systems under periodic irrigation (FW, DWW, and UTWW) were Cucurbitariaceae and Ceratostomataceae. Among the families that were significantly associated to systems under wastewater irrigation we detected Nectriaceae, Tephromelataceae, Blastocladiaceae, Mortierellaceae, Vibrisseaceae, Pichiaceae, Lasiosphaeriaceae, Heterogastridiaceae, and Sonoraphlyctidaceae (Figure 5). Most families associated to the wastewater irrigated systems were in contrary to shrubland and rainfed systems positively correlated to soil moisture, Na + content and available P (Supplementary Table S2).

Fungal Trophic Modes
Saprotrophy was the predominant trophic mode in all land use systems, showing no significant differences in relative abundances (Figure 6). The most abundant saprotrophs FIGURE 5 | Correlation-based association network at family level of the fungal communities among land use systems: shrubland, rainfed, freshwater, irrigated with dam and untreated wastewater, respectively. Only statistical significant fungal families are visualized (p ≤ 0.05). The size of the nodes is proportional to the taxon relative abundance and the edge width corresponds to the association strength of each taxon with the land use system. belonged to the Nectriaceae, Lasiosphaeriaceae, and Mortierellaceae families. The 20 most abundant saprothrophic genera are depicted in the Supplementary Figure S5. The abundance of saprotrophs was positively correlated to total C and available P in the soil, R 2 = 0.25 and 0.32, respectively (p ≤ 0.05; Table 3). Symbiotrophic fungi showed a decrease in the land use systems under periodic irrigation (FW, DWW, UTWW) in comparison to shrubland and rainfed systems (p ≤ 0.05). AMF had a higher relative abundance in the rainfed soil and lower in the UTWW system. The relative abundance of symbionts and AMF was higher during the dry season in the shrubland soil. The most abundant symbiothrophic genera were Capronia, which was only detected in the Shrubland soil and, Glomus that was detected in the land use systems under periodic irrigation (Supplementary Figure S6). The relative abundance of EcMF did not show differences among land use systems. The most abundant AMF genus was Glomus and the most abundant EcMF genera were Pezisa and Amanita. Symbiotrophs in general were negatively associated to humidity and available P (R 2 = −0.30). EcMF were also negatively associated to base saturation cations, which are increased in wastewater irrigated systems (p ≤ 0.05; Table 3). The abundance of pathotrophs including plant and animal pathogens was not higher in the wastewater-irrigated systems than in the other systems, the same was observed specifically for plant pathogens. The relative abundance of animal pathogens was higher in the FW system and lower in the UTWW one. During the dry season the abundance of plant pathogens was higher in the shrubland soil. Three of the most abundant plant pathogens detected were Bipolaris cynodontis, Lectera colletotrichoides, and Gibellulopsis nigrescens (Supplementary Figure S7). The most abundant animal pathogenic taxa detected were Metarhizium marquandii, M. anisopliae, and Metacordyceps chlamydospori. Pathotrophs in general showed negative correlation to several physico-chemical soil properties enhanced in wastewater irrigation systems, such as humidity, E.C., total N, K + , and Na + content and available P (p ≤ 0.05; Table 3).
To determine whether the wastewater increased the abundance of human pathogens in soil we searched for specific genera known for their pathogenicity in soil and wastewater samples. Pathogenic species of Aspergillus, Candida, Fusarium, Malassezia, and others were found in the wastewater samples (Supplementary Table S3). The shrubland system had a lower abundance of human pathogens, while higher abundances were detected in the wastewater-irrigated systems FIGURE 6 | Relative abundances of fungal trophic modes within each land use system. Shrubland (S), Rainfed (R), Freshwater (FW), Dam wastewater (DWW), and Untreated wastewater (UTWW) during the dry and rainy season. Boxes extend from the 25th to 75th percentile, the line in the box is plotted at the median. Whiskers represent the smallest and the largest value. A mix effect model followed by Bonferroni's multiple comparison tests were used to determine differences among land use systems, and between seasons. Different letters indicate statistical difference among land use systems and * indicate statistical differences between seasons (p ≤ 0.05). Only parameters differing significantly between seasons in each land use system are shown, if they did not differ, samples of both seasons were merged.  than in the rainfed and FW systems. However, this difference was not statistically significant. In the wastewater-irrigated systems the most abundant species were Bipolaris cynodontis, Candida quercitrusa, and Fusarium fujikuroi (p ≤ 0.05) (Supplementary Figure S8).

DISCUSSION
In this study we evaluated the effect of changing soil physical and chemical parameters due to land use change and the quality of water used for irrigation on the diversity and composition as well as the functionality of fungal community from a dryland soil.
The study was conducted in the Mezquital Valley, which provides representative conditions to study soil fungal communities under different land use systems, including preserved natural vegetation (shrubland), rainfed agriculture or irrigated agriculture with three different water qualities: freshwater, untreated wastewater stored in a dam and untreated wastewater.

Irrigation With Wastewater Decreases Fungal Richness and Diversity
Fungal diversity and richness showed a progressive decrease in the land use systems where irrigation is applied and the quality of water is worsened. The decline of fungal richness and diversity in the DWW and UTWW systems could be related specifically to the higher amount of P that the wastewater introduces into the soil (Siebe, 1998). Shrestha et al. (2015) reported a reduction of fungal diversity in fertilized soils in comparison to non-amended soil. The moisture effect in the fungal diversity and richness is also evident as the systems under wastewater irrigation were more humid than shrubland, rainfed and FW. These results are in line with several previous studies reporting a more abundant and diverse fungal community under drought or dry soil conditions, and less abundant and diverse emerging during wetter periods and in wetter soils (Hawkes et al., 2011;Acosta-Martínez et al., 2014;Ma et al., 2015;Zheng et al., 2017). The thriving of fungi in dryer conditions is related to the capability of fungal hyphae to grow in small pores and take up water that may be retained in them, as they do not experience the limitation of movement that is inherent to a unicellular structure (Chowdhury et al., 2011;Mchugh et al., 2014). Shrubland soil was the only system that showed differences between seasons, probably several groups encompassing this community are more sensitive to changes in water content that are not present in the other systems. Interestingly, fungal diversity and richness was also negatively correlated to the content of exchangeable soil cations Na + , K + , and Mg 2+ , which are enhanced in wastewater-irrigated soil. Specifically, the effect of Na + content and soil salinity on fungi is still unclear (Rath and Rousk, 2015). Some studies (Chowdhury et al., 2011;Baumann and Marschner, 2013) have reported a decreased proportion of fungal phospholipid fatty acids (PLFAs) with higher salinity. Sardinha et al. (2003) reported a decrease in the fungal biomarker ergosterol with increasing salinity under acidic soil conditions. Also, salt affected soils have shown a lower fungal to bacteria ratio (Pankhurst et al., 2001). Salinity can also have a negative effect on AM fungi, however, several reports showed improved growth and performance of mycorrhizal plants under salt stress conditions (Porcel et al., 2012). It cannot be excluded that contaminants present in the wastewater not measured in this study have affected the diversity of the fungal community.
In the soil of the Mezquital Valley, bacteria responded differently than fungi. Bacteria had higher diversity in the agricultural system under freshwater irrigation than in the shrubland system, bacterial diversity was correlated positively to pH and Ca 2+ and negatively to P and Na + content (Lüneberg et al., 2018), this suggests that the soil conditions caused by wastewater irrigation affect the diversity of fungi and bacteria, but fungi are affected more severely.

Fungal Community Composition Is Affected by Land Use
The soil fungal community composition of the Mezquital Valley is influenced primarily by whether the land conserved its natural vegetation or if it was use for agriculture, and secondly by the quality of water use for irrigation. The clustering according to land use system shows the importance of interactions between fungal community structure and vegetation diversity, as shrubland fungal community is separated from the monocultural land use systems. Similar trends were recorded during conversion of temperate and tropical forests into plantation (Kasel et al., 2008;Purahong et al., 2014;McGuire et al., 2015) or grassland (Lauber et al., 2008;Lopes Leal et al., 2013), and when natural grassland was transformed to agricultural systems (Bissett et al., 2011;Xiang et al., 2014). Regarding agricultural management, many studies focused on practices such as fertilization, tillage, crop rotation and organic farming (Morris et al., 2013;Moora et al., 2014;Wang et al., 2016), but not on water quality used for irrigation. Here, soil fungal structure shifts radically in soil under wastewater irrigation, similar to the findings of Disciglio et al. (2015) for treated wastewater irrigated soil.
The strong impact of Na + content, as well as K + and Mg 2+ on the composition of the fungal community agrees with reports in agricultural land and Mediterranean and semiarid sites, where base saturation was reported to alter total fungal or AMF community structure (Jansa et al., 2014;Ramalho da Silva et al., 2014;Ciccolini et al., 2015). Additionally, in a dry forest Gonçalves Lisboa et al. (2014) found that fungi to bacteria ratio was related to soil cations. The effect of base saturation can be related to the inhibition of spore germination, hyphal development and the production of arbuscules by salinity, as has been observed for AMF (Miransari, 2010). Compared to the well-established role of pH in the determination of soil bacterial composition, the effect of pH on the soil fungal composition is still under debate. On one hand, studies reported pH as one of the main drivers of total and AM fungal composition (Rousk et al., 2010;Jansa et al., 2014;Moora et al., 2014;Qin et al., 2015), or a pH effect on the abundance of fungi (Tedersoo et al., 2014;Maestre et al., 2015) and the fungi to bacteria ratio (Rousk et al., 2010). On the other hand, studies are available that did not identify soil pH as a driver of the fungal community composition (Lauber et al., 2008;Ramalho da Silva et al., 2014;Ciccolini et al., 2015). Here, pH was one of the strongest drivers of the community given the narrow pH range of the soils (7.1-7.8). In general, fungi show optimal growth under a wide pH range (5-9), although different taxa of soil fungi showed preference for certain pH values (Rousk et al., 2009). However, pH is directly or indirectly related to different soil properties such as nutrient availability, solubility of metals and salinity (Lauber et al., 2009;Rousk et al., 2009).
The composition of bacteria in the Mezquital Valley soils was also determined by Na + content, pH and soil moisture (Lüneberg et al., 2018), which exemplifies the importance of these soil properties when drylands are transformed to agricultural land and are irrigated with untreated wastewater. In the soil of the Mezquital Valley fungi seem to be more sensitive to salinity than bacteria, this is in line with the report of Baumann and Marschner (2013).
The fungal community of this study was dominated by Ascomycota and Basidiomycota as previously reported for most dryland and agricultural soils (Steven et al., 2014;Maestre et al., 2015;Wang et al., 2016). The structural shifts in fungal communities with a decrease of Ascomycota and the increase of Zygomycota and Rozellomycota in the wastewater-irrigated systems indicate lineage-specific adaptations of higher taxonomic groups. The correlation-based analysis of associated taxa showed half of the groups were associated to one or more land use systems. Most fungal taxa significantly associated to shrubland soil belong to the classes Eurotiomycetes and Dothideomycetes of the Ascomycota phylum; this agrees with previous reports that have shown Ascomycota dominating fungal communities of semiarid grasslands and drylands (Porras-Alfaro and Bayman, 2011;Maestre et al., 2015;Vargas-Gastélum et al., 2015). Specifically Eurotiomycetes have been associated to abandoned agricultural soil compared to maize monoculture system (Ciccolini et al., 2015). Among the taxa associated to shrubland, several lichenized families were identified, such as Chaetothyriaceae, Verrucariaceae, Teloschistaceae, and Parmeliaceae, the two latter belong to the class Lecanoromycetes, from where are the most lichen fungi that form symbioses with filamentous Nostoc cyanobacteria (Grube and Wedin, 2016). Cyanobacteria was one of the indicator bacterial taxa of shrubland soil (Lüneberg et al., 2018), as well as the order Rhizobiales which is a lichen-associated group (Grube and Wedin, 2016). Most of the taxa associated to shrubland and rainfed soil were adapted to drought, hence sensitive to Na + content and P available. Fungal families associated to rainfed soil included two lichen families of the order Pezizales. While the taxa associated to wastewater-irrigated soil, such as Nectriaceae, Tephromelataceae, Blastocladiaceae, Mortierellaceae, Vibrisseaceae, Pichiaceae, Heterogastridiaceae, and Trichosporonaceae were adapted to higher soil moisture, P content and concentration of Na + provided by wastewater irrigation (Siebe, 1998). Nectriaceae includes numerous important human and plant pathogens being maize one of its hosts, also several species are used in industrial applications as biodegraders and biocontrol agents (Lombard et al., 2015). The families associated to wastewater-irrigated soil might be also adapted to surfactants and heavy metals present in wastewater (Friedel et al., 2000) not measured here.

Soil Moisture and Available P Define Fungal Trophic Modes
To translate fungal sequence information into ecologically meaningful information we employed FUNguild (Nguyen et al., 2016) and assigned ASV sequences to trophic modes. The lack of differences in the relative abundance of saprothrophs observed in this study can be related to the wider range of heat and drought conditions, saprophytic populations are able to resist (Acosta-Martínez et al., 2014), this trend was previously observed comparing forest and grassland soil (Nguyen et al., 2016). Our results about the saprotrophic fungi abundance differs from the findings of Disciglio et al. (2015), who reported higher colonization of saprophytic species in soil irrigated with treated wastewater compared to groundwater-irrigated soil. The positive correlation of saprotrophic fungi to total C and N supports the importance of saprothrophs as key regulators of nutrient cycling (Crowther et al., 2012). The abundance of symbiothrophs is related to soil moisture, AMF are known to resist a wide range of heat and drought conditions (Acosta-Martínez et al., 2014), and especially AMF seem to be sensitive to water quality. Previous studies showed specifically a reduced diversity of AMF in treated wastewater irrigated soil due to increase in nutrient availability (Alguacil et al., 2012;Bastida et al., 2017). Additionally, other studies have shown that AM fungal richness and diversity are highly sensitive to available P (Bai et al., 2013;Yang et al., 2014;He et al., 2016). The results we obtained agreed with the previous culture-based study of AMF in the soils of the Mezquital Valley by Ortega-Larrocea et al. (2001), who found a reduction of AMF abundance in soil irrigated with wastewater, and Glomus to be the most abundant fungal genus. It is important to notice that the general fungal ITS2 primers used in this study might underrepresent AMF, however, the primers are adequate to evaluate treatment responses (Lekberg et al., 2018). Ectomycorrhizal fungi are primarily associated to trees, its propagation is potentially ceased in agricultural systems, due to irrigation and the use of fertilizers and tillage (Munyanziza et al., 1997;Tedersoo and Smith, 2013;Zheng et al., 2017). The abundance of pathothrophs was not higher in the wastewater irrigated systems as we hypothesized (iii). These results can be possibly related to the higher abundance of bacterial taxa such as Pseudomonas in the soil irrigated with wastewater (Lüneberg et al., 2018). Pseudomonas are known for their antifungal compound production (Ligon et al., 2000). A similar trend was reported for tomato crops irrigated with secondary-treated agro-industrial wastewater, revealing a decrease of the pathogen Fusarium spp. in the soil close to the roots (Disciglio et al., 2015). The most abundant plant pathogens Bipolaris cynodontis, Lectera colletotrichoides and Gibellulopsis nigrescens are known for infections in maize and lucerne, the main crops of the area (Hu et al., 2011;Inderbitzin et al., 2011;Cannon et al., 2012;Manamgoda et al., 2014). The most abundant animal pathogen taxa detected were Metarhizium marquandii, M. anisopliae, and Metacordyceps chlamydospori which are recognized entomopathogens and pathogens of nematodes and snails (Bischoff et al., 2009;Kepler et al., 2012;Stirling, 2014;Garrido-Jurado et al., 2015;Sanjaya et al., 2016). Surprisingly, most of the taxa associated with distinct land use systems were either negatively or positively correlated with available P. Additionally, P content was the only soil property correlated with the abundance of all the different fungal trophic modes. These results demonstrate the importance of soil nutrient status for the fungal community in each land use system (He et al., 2016). The human pathogens found in the wastewater irrigated systems (Bipolaris cynodontis, Candida quercitrusa, and Fusarium fujikuroi) produce opportunistic infections in immunocompromised patients (Guarro, 2013;Manamgoda et al., 2014;Xiao et al., 2014). As we could detect more pathogenic species in the systems where wastewater is used for irrigation of the crops, we assume that most pathogenic fungi introduced by wastewater were unable to tolerate the soil conditions.

CONCLUSION
In line with the general hypothesis (i) dryland transformation to agricultural land under irrigation with different water quality has a significant effect on the soil fungal diversity and composition. Soil under natural semi-arid shrubland has the highest fungal richness and diversity, while periodic wastewater irrigation negatively affects these fungal indices, by the increment of soil moisture, available P and base saturation content. In accordance with our hypothesis, soil moisture in addition with sodium content and soil pH were the strongest drivers of the fungal community structure. The (dis)-similarities in the structure of the fungal community in the soils of the Mezquital Valley are primarily due to whether the land conserved its natural vegetation or if it was used for agriculture, and secondly by the quality of water used for irrigation. In contrast with our hypothesis (ii), the structure of the fungal community did not differ significantly along the seasons in any of the studied land use systems, despite differences in soil moisture and vegetation cover. The observation from the fungal trophic modes oppose our hypothesis on an increase of saprophytic and pathogenic fungi in soils under periodic irrigation with wastewater (iii); interestingly we observed the opposite, saprotrophic and pathotrophic fungi (including plant and animal pathogens) showed no differences among land use systems. However, three of the most abundant plant pathogens identified are known for infections in maize and lucerne, the main crops of the area. Symbiotic fungi were affected by irrigation and AMF specifically by irrigation with untreated wastewater. The results indicate lineage-specific adaptations to soil properties, as wastewater irrigation reduces the abundance of Ascomycota, while the abundance of Zygomycota and Rozellomycota is increased. We found that fungal families associated to shrubland, rainfed and partially FW were adapted to drought, hence sensitive to exchangeable sodium content and most of them to total N and available P. Taxa associated to freshwater, dam wastewater and untreated wastewater irrigated systems showed the opposite trend. Additionally, the identified potentially harmful human pathogens might be a health risk for the population. The present study contributes to better understand the impact of shrubland transformation into agricultural systems irrigated with different water quality on soil fungal communities in managed semiarid ecosystems.

AUTHOR CONTRIBUTIONS
CS, RD, and KL conceived the study. KL carried out the field and laboratory work. DS and KL prepared and analyzed the data. CS, RD, DS, NB, and KL interpreted the results and wrote the manuscript.

FUNDING
The work has been funded by the Mexican Consejo Nacional de Ciencia y Tecnología (CONACyT, Grant No. CB-221789). KL studies in Mexico and her study grant and travel expenses in Germany were supported by CONACyT and the Posgrado en Ciencias Biológicas, UNAM.