2.5 | Metagenomic sequencing and analysis
Each representative DNA sample from the sea cucumber intestine for which sufficient volumes were available after 16S typing was used for metagenomic sequencing. Metagenomic DNA paired-end libraries were prepared with an insert size of 350 bp and were quantified using a Qubit Fluorometer and an Agilent 2100 TapeStation system. Sequencing was performed on an Illumina HiSeq PE150 platform.
Raw reads were preprocessed using FasqMcf to exclude adapter sequences and low-quality sequences (Aronesty, 2011), and reads derived from host contamination were filtered using bowtie2.2.4 (Langmead & Salzberg, 2012). Clean reads were assembled and analyzed by SOAPdenovo software (R. Luo et al., 2012). Then, gene prediction was performed on contigs larger than 500 bp by MetaGeneMark software with the default parameter, and gene models with CDS lengths less than 100 bp were filtered out (Nielsen et al., 2014; Qin et al., 2010). A gene catalog was constructed using the gene models predicted from each sample by CD-HIT-EST (version 4.6.6) (Li & Godzik, 2006) with the parameters ‘-c 0.95 -n 10 -G 0 -a S 0.9’, which adopts a greedy incremental clustering algorithm and criteria of identity > 95% and overlap > 90% of the shorter genes.
DIAMOND software (version 0.9.9) (Buchfink, Xie, & Huson, 2015) was used to align unigenes to the sequences of bacteria, fungi, archaea and viruses, which were all extracted from the NCBI nr database (version: 2018-01-02) with the parameter ‘-e 1e-5’. Functional assignments of protein sequences were made on the basis of DIAMOND alignment against the KEGG protein database (version 2018-01-01) (Minoru et al., 2006), eggNOG database (version 4.5) (Sean et al., 2014), and CAZy database (version 20150704) (Cantarel et al., 2009) by using the best hit with an e value < 1e-5. For each sequence’s BLAST result, the best BLAST hit was used for subsequent analysis.
Phylum, class, order, family, genus, species, KEGG orthology (KO), and orthologous group (OG) relative abundances were calculated by summing the abundance of the respective genes belonging to each category per sample based on the taxonomic assignments and KO and OG annotations. The relative gene abundance profile was also summarized into KEGG, eggNOG and CAZy functional profiles for functional analysis.
2. 6 | Quantitative real-time PCR (qPCR) assays
To assess the amounts of intestinal bacterial in different samples, we also quantified the bacteria through qPCR using the primers 341F and 518R, which are universal primers targeting a short fragment (ca. 171 bp) of the bacterial 16S rRNA genes. The qPCR assay was based on the fluorescence intensity of the SYBR Green dye and performed as previously described (F. Sun et al., 2015). Using the vector-targeted primers M13F/M13R, linear fragments were obtained from PCR amplification of circular plasmids (pTZ57R/T vector; Fermentas), which contained inserts of the bacterial 16S rRNA gene fragments. The standard DNA was also quantified using the PicoGreen dsDNA reagent kit (Invitrogen).