Illumina MiSeq sequencing

Purified amplicons were pooled in equimolar concentrations and paired-end sequenced (2 × 300) using an Illumina MiSeq platform (Illumina, San Diego, CA, USA) according to standard protocols.

Raw fastq files were demultiplexed, quality-filtered using Trimmomatic and merged using FLASH with the following criteria. (i) Reads were truncated at any site with an average quality score <20 over a 50-base pair (bp) sliding window. (ii) Primers were exact matches, allowing 2 nucleotide mismatches, and reads containing ambiguous bases were removed. (iii) Sequences whose overlap was longer than 10 bp were merged according to their overlap.

Operational taxonomic units (OTUs) were clustered with a 97% similarity cutoff using UPARSE (http://drive5.com/uparse/) and chimeric sequences were identified and removed using UCHIME. The taxonomic assignments for each 16S rRNA gene sequence were obtained by the RDP Classifier algorithm (version 2.2: http://sourceforge.net/projects/rdp-classifier/) against the SILVA (Release119: http://www.arb-silva.de) 16S rRNA gene database using a confidence threshold of 70% (Li et al., 2016). A total of 1,159,668 valid sequences and 511,762,953 bp were obtained by sequence filtration and double-stitched splicing, and the average sequence length was 441 bp. After subsampling each sample to an equal sequencing depth (23,863 reads per sample) and clustering, 3,132 OTUs at 97% identity were obtained.