Leaves and seeds of cultivated rice (Oryza sativa) showing symptoms of BSR were collected in 2018 and 2019 in different rice fields from Bagré, Bama, Banfora, Banzon and Sourou ( Figure 1 ; Supplementary Tables S1, S2 ). These sites are on average 365 km apart (minimum distance 59 km and maximum distance 622 km) and represent the main rice-growing areas of Burkina Faso. Diseased leaves and seeds were arranged in Petri dishes containing a double layer of moistened blotting paper (Mathur and Kongsdal, 2003). Incubation was carried out at 25°C for 7 days under alternating 12 h of white light and 12 h of darkness. Using a flame-stretched Pasteur pipette, fungal conidia formed in the lesions were collected and spread on bacto-agar medium. After incubation for 24 h (25°C, 12 h white light and 12 h darkness), a single germinating conidium was transferred with a sterile scalpel onto potato dextrose agar medium (PDA, 39 g/L) and incubated under the same conditions for 7 days.

Genomic DNA was extracted from fresh mycelium with modifications to the cetyl trimethyl ammonium bromide (CTAB) method described by Doyle and Doyle (1987). A 5-day-old fragment of mycelium was transferred onto Corn Meat Agar (CMA, 18 g/L) medium previously covered with a sterile cellophane disc. The plates were incubated at 25°C (12 h photoperiod) for 7 days, the mycelium was scraped off with a scalpel and approximately 30 mg of the mycelium was transferred to a 2 mL tube containing sterile sand and 2 iron beads. These tubes were placed directly into liquid nitrogen and the mycelium was cold ground using a Retsch mill (Retsch MM 301) for 4 min at 15 Hz. After grinding, each tube received 1 mL of CTAB lysis buffer with 1% polyvinyl pyrrolidone (MW. 360 000), 1% Na₂SO₃, 0.2 mg/mL RNase and 1 mg/mL proteinase K. The tubes were vortexed (Vortex-Genie 2) and incubated in the oven at 65°C for 2 h. The tubes were then centrifuged for 3 min at 8,000 rpm and approximately 900 µL of the supernatant was transferred to a new 2 mL tube. One mL of chloroform isoamilalcohol (CIAA, 24/1 v/v) was placed in each tube and mixed by inversion before being centrifuged at 10,000 rpm for 10 min. After centrifugation, 700 µL of the supernatant was transferred to 1.5 mL tubes. Six hundred µL of isopropanol cooled to -20°C and 210 µL of 5M NaCl were added and mixed by inversion. The tubes were kept at -20°C overnight. They were then centrifuged at 4°C for 30 min at 13,000 rpm and the supernatant was removed by inverting the tubes. Five hundred µL of 70° ethanol was added to each tube and centrifuged at 4°C for 5 min at 13,000 rpm. Ethanol was removed by inverting the tubes and draining on absorbent paper. The tubes were dried for 1 h in the fume cupboard. DNA was re-suspended in 100 µL of sterile milliQ water. The quality and quantity of the extracted DNA was assessed by electrophoresis in 1% (w/v) agarose gels and spectrophotometry (Nanodrop™ 2000).

The preparation of GBS libraries was carried out in collaboration with the genotyping platform of the UMR AGAP/CIRAD according to the method described by Elshire et al. (2011). Genomic DNA of the strains were normalized to a concentration of 10 ng/µL and digested with the restriction enzyme ApeKI (5 U) at 75°C for 2 h. The restriction digestion of DNA was followed by ligation with adapters. The adapters included different barcodes for tagging individual samples and common adapters. The ligation was performed using T4 DNA ligase at 22°C for 30 min and the ligase was inactivated by holding at 65°C for 10 min. Two pools were created. To assess the reproducibility of GBS, each pool contained 6 randomly selected and duplicated isolates. Ligation products were purified with the QIAquick purification kit from QIAgen. The library was enriched by PCR to 50 µL of reaction solution including 25 µL of 2X Taq Master mix NEB, 19 µL of water, 2 µL of DNA and 2 μL each of the forward and reverse primers (10 μM). Amplification was done according to the following program: 98 ˚C for 30 s, followed by 18 cycles of 98 ˚C for 10 s, 65 ˚C for 30 s, 72 ˚C for 30 s, ending at 72 ˚C for 5 min. PCR products were purified with QIAgen’s QIAquick Purification Kit before being evaluated with BioAnalyzer 2100 (Agilent Technologies) for fragment size distribution. The GBS libraries were sequenced on an Illumina HiSeq 4000 producing 150 bp paired-end reads.

The bioinformatics tool RattleSNP (https://github.com/sravel/RattleSNP) was used to process the data. The workflow performs demultiplexing, sequence mapping to the reference genome and variant calling. The sequences of B. oryzae isolates were mapped onto the genome of strain ATCC 44560 (Condon et al., 2013), and those of E. rostratum onto the genome of BF9006 sequenced by CIRAD (https://www.ncbi.nlm.nih.gov/biosample/SAMN30522240/). Prior to variant calling, strains with a mapping rate of less than 75% to the reference genome, a library size of less than 100,000 reads and a sequencing depth of less than 3 x were considered low quality data and excluded from the dataset. SNPs from strains meeting the criteria (61 isolates for B. oryzae and 151 isolates of E. rostratum) were filtered for missing rate (–max-missing 0.70), minor allele frequency (–maf 0.05), depth (–minDP 2) and site quality (–minQ 30) using VCFtools version 01.16 (Danecek et al., 2011)

Discriminant Analysis of Principal Components (DAPC) (Jombart et al., 2010) and sNMF (Frichot et al., 2014) analysis were used to investigate the genetic structure of these fungal populations. DAPC was carried out with the adegenet package of the R software. First, the find.clusters function was used to detect the number of clusters in the population (Jombart, 2008; Jombart et al., 2010). The function uses K-means clustering which decomposes the total variance of a variable into two components: inter-group and intra-group. This model maximizes the inter-group variation (Jombart et al., 2010).

sNMF of the LEA R software package was used to estimate the number of discrete genetic clusters (K) and admixture coefficients (Frichot et al., 2014). The optimal number of sNMF clusters was defined as the one corresponding to the lowest cross-entropy value.

To further characterize the genetic structure of the population described by the DAPC and sNMF analyses, a hierarchical analysis of molecular variance (AMOVA) (Excoffier et al., 1992) were performed with the poppr package of R (Kamvar et al., 2014). Pairwise genetic differentiation (F_ST) was calculated using the ‘genet.dist’ function in the hierfstat package of R (Goudet, 2005), following the method described by Weir and Cockerham (1984). RAxML (Randomized Axelerated Maximum Likelihood) version 8.2.12 was used to assess phylogenetic relationships between isolates using the maximum likelihood method (Stamatakis, 2014).

Analysis of genotypic diversity and association indices of fully or partially clonal populations is recommended before and after clone correction (Milgroom, 1996; Grünwald and Hoheisel, 2006; Kamvar et al., 2014). To identify clones, the dissimilarity rate between duplicates of the same isolate was calculated using the bitwise.dis function of the R poppr package (Kamvar et al., 2014). For a given species, the maximum dissimilarity between duplicates was measured and used as a threshold to define clones. Two isolates with a dissimilarity rate lower than the threshold were considered as clones.

To estimate genotypic diversity, the poppr package of R (Kamvar et al., 2014), was used to calculate Shannon-Wiener index (H) (Shannon, 1948), Simpson’s complement index of genotypic diversity index (λ) corrected for sample size (Simpson, 1949; Grünwald et al., 2003) and evenness (E5) (Pielou, 1975; Ludwig and Reynolds, 1988; Grünwald et al., 2003). To adjust for the population size, the raw value of the Simpson’s complement index of genotypic diversity index was multiplied by N/(N-1) with N being the number of isolates in the population (Grünwald et al., 2003).

To assess gene diversity, nucleotide diversity (π) values were calculated with VCFtools across the genome with 10 kb windows (Danecek et al., 2011).

Detection of recombination was carried out by calculating linkage disequilibrium (LD) and clonality was tested with the PHI test from genotyping data. Pairwise LD between loci was calculated based on allelic frequency correlations (r²) over a 10 kb window using the PopLDdecay program (Zhang et al., 2019).

A FASTA file was used to generate a NeighborNet in SPLITSTREE4 and the PHI test was performed to test the null hypothesis of clonality (Huson and Bryant, 2006).

Identification of mating types was done by amplification of MAT1-1 and MAT1-2 genes with primers specific to B. oryzae (Ahmadpour et al., 2018) and E. rostratum (Kusai et al., 2016). Each PCR reaction using the kit Go-Taq (Promega) had a final volume of 25 μL, containing 1 μL of template genomic DNA (20 ng/µL),1 μL forward primer (10 µM), 1 μL reverse primer (10 µM), 1µL of dNTPs (10 µM each), 5 μL of the buffer 5x, 0.25 µL Taq polymerase (5U/µL) and 14.75 μL sterile water. For the amplifications of E. rostratum MAT genes, the amount of water was 13.75 µL and 1 µL MgCl₂ 5X was added. PCR cycling was performed as follows: 94°C for 2 min; followed by 35 cycles of 94°C for 1 min, 52°C for 1 min and 72°C for 1 min; plus a final extension at 72°C for 7 min for B. oryzae. For E. rostratum, initial denaturation at 95°C for 2 min; followed by 35 cycles of denaturation at 95 C for 30 s, annealing at 60°C for 30 s, extension at 72°C for 1 min; and final extension at 72°C for 5 min. The null hypothesis of a 1:1 ratio between mating types was tested using the chi-square goodness-of-fit.

The sexual fertility of B. oryzae and E. rostratum isolates was tested in vitro according to the protocol described by Tsuda and Ueyama (1976). A dry, autoclaved rice leaf of about 4 cm was placed in the center of a petri dish on Sach’s culture medium. Pieces of mycelium from the pre-determined MAT1-1 and MAT1-2 strains were placed on either side of the rice leaf. The dishes were sealed with tape to prevent dehydration of the medium and then placed in an oven at 24°C for 30 days.

The genetic diversity of isolates of B. oryzae (n=61) and E. rostratum (n=151) collected in Burkina Faso, was estimated using genetic markers (SNPs) obtained by genotyping-by-sequencing. After alignment to the B. oryzae ATCC 44560 reference genome composed of 619 scaffolds and filtering, 28,907 SNPs were identified in the B. oryzae samples (61 isolates). Sequences of E. rostratum isolates (n=151) aligned to the BF9006 isolate reference genome (377 scaffolds) yielded 102,464 SNPs.

To study the B. oryzae and E. rostratum populations from Burkina Faso, the isolates were clustered by DAPC according to their multilocus genotypes. Inference of ancestry was also performed with sNMF.

For B. oryzae, the bayesian information criterion (BIC) of DAPC suggested the best clustering for four genetically distinct populations (K = 4; Supplementary Figure S1 , Figure 2A ). Similarly, the sNMF analysis indicated an optimum for four clusters ( Supplementary Figure 2 , Figure 3A ). The phylogenetic tree made with RaxML also highlighted four clades ( Supplementary Figure S3 ). The clustering of isolates according to DAPC, sNMF and RaxML was similar. Isolates from cluster 1 and 2 were found in all five geographical collection areas ( Table 1 , Figure 4 ), while cluster 4 was only found in Banfora. Isolates from cluster 3 were detected in three areas namely Bagré, Bama and Sourou. The sNMF analysis showed a significant admixture between the four clusters, but particularly between cluster 1 and 2 ( Figure 3A ). The F_ST values between the B. oryzae isolates from the different rice-growing areas in Burkina Faso showed little differentiation between sites, ranging from 0.005 to 0.132 ( Table 2 ) and with an average value of 0.08. No correlation was observed between geographic distances of rice growing sites and genetic differentiation index values ( Supplementary Figure S4 ). The AMOVA analysis revealed that the variability was 88.65% within localities and 11.35% among the five rice-growing zones ( Table 3 ).

For E. rostratum, the bayesian information criterion (BIC) of DAPC suggested the best clustering for four genetically distinct populations (K = 4; Supplementary Figure 5 , Figure 2B ). Similarly, the sNMF analysis indicated an optimum of four clusters ( Supplementary Figure S6 ; Figure 3B ). The phylogenetic tree made with RaxML also highlighted four clades ( Supplementary Figure S7 ). The clustering of isolates was identical between methods. Isolates from clusters 1 and 2 were found in all five localities sampled. Isolates from cluster 3 were found in all localities ( Figure 5 ) except Banzon, and isolates cluster 4 was sampled in Bagré, Banzon and Sourou. Very low F_ST values ranging from 0 to 0.129 and with an average value of 0.022 were observed between populations of different rice production areas in Burkina Faso ( Table 4 ). Distances between the different rice production sites had no impact on population structure ( Supplementary Figure S8 ). The AMOVA analysis revealed that the variability was 94.4% within and 5.6% among the five rice-growing areas ( Table 5 ).

Samples were analyzed for possible clones. The genetic distance between duplicates of B. oryzae isolates varied from 0.035 to 0.087. Those of E. rostratum varied between 0.015 and 0.033. Bipolaris oryzae and E. rostratum isolates were considered as clones if their dissimilarity rates were less than or equal to 0.087 and 0.033 respectively. Based on this criterion, 8 and 11 isolates of B. oryzae and E. rostratum respectively were considered as clones of two or more isolates. Clones represented 13% (8/61) and 7% (11/151) of the isolates in each sample of B. oryzae and E. rostratum, respectively ( Table 6 ).

Genotypic and gene diversity indexes are presented in Table 6 . For B. oryzae, the genotypic diversity (H) between isolates from different rice-growing areas ranged from 2.20 to 3.04 with an average value of 2.45. These indices are influenced by population size and do not allow comparison between different localities. The Simpson index corrected for sample size ranged from 0.999 to 1 (mean value 0.999) indicating a high level of genotypic diversity within each population. The regularity index (E.5) was equal to 1 for all populations.

In terms of genetic diversity, nucleotide diversity values ranged from 1.3 x10^-4 to 2.4 x10^-4 with an average of 1.9 x10^-4. The Banzon population also showed the lowest value of nucleotide diversity (1.3 x10^-4). The Sourou population showed the highest diversity (2.4 x10^-4). It was followed by the populations of Bama, Bagré and Banfora.

For E. rostratum, the genotypic diversity (H) between isolates from different rice-growing areas ranged from 1.61 to 4.23 with an average value of 2.98. The Simpson index corrected for sample size ranged from 0.999 to 1 (mean value 0.999) indicating high genotypic diversity within each population. The regularity index (E.5) was equal to 1 for all populations. The values of nucleotide diversity varied between 5.7 x10^-4 and 9.8 x10^-4. The average value was 4.8 x10^-4. The isolates from Bama showed the highest value of nucleotide diversity (9.8 x10^-4), while the lowest value was observed in Bagré (5.7 x10^-4).

In B. oryzae populations from Burkina Faso, Linkage Disequilibrium (LD) pattern showed a rapid decay to less than 1 kb in all populations ( Figure 6 ) suggesting high rate of recombination. The PHI coefficient calculated on the populations with sufficient number of isolates, rejected the hypothesis of clonality (p = 0.0). Finally, a reticulated NeighborNet was observed ( Figure 7A ), also supporting the hypothesis of recombination in populations of B. oryzae from Burkina Faso. Finally, the mating types MAT1-1 and MAT1-2 were observed in balanced proportions ( Table 7 ).

In E. rostratum, LD showed a rapid decay to less than 1 kb in all populations ( Figure 6 ) and PHI coefficient rejected the hypothesis of clonality (p = 0.0), supporting the hypothesis of recombination. Low reticulated NeighborNet was observed ( Figure 7B ). The dominance (79%) of only one mating type (MAT1-2) was observed and the null hypothesis of a 1:1 ratio was rejected by statistical analysis. Ten strains repeatedly amplified markers of the two mating types ( Table 7 ). Potential mixture of strains was discarded by repeating monoconidial isolation and potential contamination was discarded by repeated independent DNA extractions.

To test the sexual fertility of B. oryzae and E. rostratum isolates, crosses between isolates of opposite mating type from the same rice field were performed in vitro. Out of a total of 156 in vitro crosses for B. oryzae, two pairs of isolates (BF9259 and BF9268; BF9294 and BF9296) produced pseudothecia containing ascospores. Concerning E. rostratum, out of 364 crosses made, no pseudothecia were observed.