Experiment design
Details on the experimental design (Table 1) and clinical findings and L. intracellularis infection characteristics including pathology and shedding results unrelated to this investigation have been published [24]. In brief, 3-week-old cross-bred pigs (the paternal line was a Hampshire decent i.e. PIC®327, while the maternal line was based on Landrace/Large White descents i.e. FAST 276) were randomly assigned to one of six treatment groups with 10-20 pigs each. The pigs were housed in four separate BSL-2 rooms with 1–3 pens (approximately 10 m2 each) of 10 pigs each equipped with a nipple waterer and a self-feeder. Specifically, the T01-LAW, T02-LAW and T03-LAW pigs were housed in a single room but in three separate pens approximately 2 m apart from each other. Similarly, and the POS-CONTROL group was also housed in a single room in two separate pens. The NEG-CONTROL and VAC-LAW pigs were housed in separate single rooms with a single pen each [24]. At 3 weeks of age, the pigs were either vaccinated with a commercial oral vaccine against L. intracellularis (VAC-LAW, n=10), were supplied feed supplemented with one of the three Bacillus probiotics (T01-LAW, B. amyloliquefaciens, n=10; T02-LAW, B. licheniformis, n=10; T03-LAW, B. pumilus, n=10), or remained non-treated (NEG-CONTROL, n=10; POS-CONTROL, n=20). At 7 weeks of age all groups, excluding the NEG-CONTROL pigs, were challenged with gut homogenate containing L. intracellularis by gastric gavage. The gut homogenate, obtained from a commercial vendor, was produced by collecting intestines from a field pig suffering from ileitis. The affected mucosa was collected and processed (i.e. grinded and diluted), and the resulting homogenate was passaged in experimental specific pathogen free pigs to obtain mucosal scrapings with a high L. intracellularis load free of other common pig pathogens. A homogenized stock was used for inoculation of all pigs to guarantee that each pigs was equally exposed to microbial loads that potentially could have been present in the inoculum. All pigs were euthanized 16 days post L. intracellularis challenge, at approximately 9 weeks old [24].
Table 1. Experimental design
|
|
|
Key events
|
Treatment group
|
Number of pigs
|
Feed information (administered for the entire study duration)
|
Vaccination
3 weeks of age
-28 dpcb
|
Challenge
7 weeks of age
0 dpc
|
Necropsy
9 weeks of age
16 dpc
|
T01-LAW
|
10
|
Base diet + Bacillus amyloliquefaciensa
|
None
|
Lawsonia intracellularis challenge using a gut homogenate
|
Ileum content collection for this study
|
T02-LAW
|
10
|
Base diet + Bacillus licheniformisa
|
None
|
T03-LAW
|
10
|
Base diet + Bacillus pumilusaa
|
None
|
VAC-LAW
|
10
|
Base diet only
|
Commercial vaccine administrationc
|
POS-CONTROL
|
20
|
Base diet only
|
None
|
NEG-CONTROL
|
10
|
Base diet only
|
None
|
None
|
a Supplemented in the feed mill with 1 × 1012 colony forming units (CFU) of the respective Bacillus strain.
b Day post Lawsonia intracellularis challenge.
c Enterisol® Ileitis, Boehringer Ingelheim, serial number 3040187B, via the oral route by drenching 2 ml of the vaccine, reconstituted as per manufacturer’s instructions, into the mouth of each pig.
Measuring L. intracellularis infection kinetics
L. intracellularis shedding was evaluated by measuring the bacterial load in rectal swabs. Rectal swabs were collected at dpc 0, 2, 4, 6, 8, 10, 12 and 15 using polyester swabs and stored in 5 mL plastic tubes containing 1 mL of sterile saline solution at −80°C until testing. Testing was done using a quantitative real-time PCR specific for L. intracellularis as described [24]. Prior to analysis, genomic copy numbers were log10 transformed to normalize them. These values were used to create the total area under the curve (AUC) for each group, the outcome utilized in this analysis. L. intracellularis associated histopathology was described using a combined ileum lesion score based on microscopic lesions and amount of L. intracellularis antigen within lesions. In brief, individual scores of crypt enterocyte hyperplasia (range from 0=normal to 3=marked with or without crypt herniation into the submucosa), inflammation (range from 0=normal to 3=marked cellular infiltrate with or without submucosal infiltrate) and amount of L. intracellularis antigen as determined by immunohistochemistry (range 0=no signal to 3=most crypts in most sessions with marked apical enterocyte labelling) were combined for a maximal score of 9 [24].
Sample collection and sequencing
At necropsy, luminal ileum content was collected from each pig approximately 10 cm anterior to the ileocecal junction and stored at -80°C until use. In addition, the inoculum stock, used to infect the pigs, was also processed and sequenced. Approximately 1 gram (wet weight) of the luminal ileum content from each pig or 1 g of the challenge material was used to extract genomic DNA using PowerMag DNA Isolation Kit (MO BIO Laboratories, Inc. Carlsbad, CA, USA) following the manufacturer’s instructions. The quality of the DNA was assessed by Qubit 4 Fluorometer (ThermoFisher, Waltham, MA, USA). The 16S rRNA amplicon library was prepared as described [25] with minor modifications. Briefly, the V4 region of the 16S rRNA gene was sequenced using the Illumina MiSeq sequencing platform with v2 MiSeq cartridges to produce 2 × 250 bp paired end reads. DNA was amplified by using the 515f/806r primer set: Forward V4 (GTGCCAGCMGCCGCGGTAA), Reverse V4 (GGACTACHVGGGTWTCTAAT). The primers for library preparation were designed based on the V4 primers, multiplex indices and Illumina adapter sequences. Extracted genomic DNA was PCR amplified with Schloss lab indices and AccuPrime™ Pfx SuperMix (ThermoFisher, Waltham, MA, USA) using the following cycling conditions: initial denaturation at 95°C for 2 min followed by 30 cycles of denaturation at 95°C for 20 sec, annealing at 55°C for 15 sec, extension at 72°C for 5 min, and a final extension at 72°C for 10 min. The prepared library was purified using Agencourt AMPure XP beads, (Beckman Coulter®, Brea, CA, USA). Libraries were quantified using Qubit dsDNA HS assay kit (ThermoFisher, Waltham, MA, USA). The pooled library was fed into the sequencing workflow with 15% PhiX control library at a final concentration of 4pM, using MiSeq Reagent Nano kit v2 (Illumina, San Diego, CA, USA) for a 500 cycle (2×250 bp) run. Custom sequencing primers designed on the V4 region were utilized during the sequencing procedure, together with MiSeq sequencing primers (Illumina, San Diego, CA, USA) as follows: read 1 primer (TATGGTAATTGTGTGCCAGCMGCCGCGGTAA), read 2 primer (AGTCAGTCAGCCGGACTACHVGGGTWTCTAAT) and index primer (ATTAGAWACCCBDGTAGTCCGGCTGACTGACT). Raw data were de-multiplexed based on dual indices to generate two fastq files, R1 and R2, for each sample.
Sequence processing using Qiime-2
The paired end reads from the sequencing step above were processed using the Quantitative Insights Into Microbial Ecology-2™ (QIIME-2) software [26]. Combining the raw reads and their corresponding experimental metadata using a manifest generated a Qiime-2 analysis file. The reads were then de-replicated and chimeric sequences removed before the denoising was done in DADA2 [26]. The resultant feature tables were used to identify the appropriate sequence depth at which the subsequent analyses would be performed. In addition, output data from a rarefaction analysis step was used to identify and remove spurious sequences at the appropriate sequence depth. Alpha and beta diversity indices as well as the operational taxonomic unit (OTU) phylogenetic were obtained by following the QIIME-2 workflow [26]. OTU taxonomic classification was done using a naïve Bayes classifier trained on the SILVA database (vMarch 2020) at 97% similarity [27]. First, a training dataset using the sequencing primers was extracted. The classifier was then used after the training to assign taxonomy to our OTU dataset.
Ileal microbiota structure analysis
The analysis of the microbiota community structure was based on characteristics of the alpha and beta diversities. In brief, data from QIIME-2 was exported into the R statistical package (v3.5.1). Four files, including the taxonomic map, the OTU biome table, the experimental metadata and the rooted phylogenetic tree, were combined to generate a phyloseq object [28] which was used for downstream analysis. To analyze the alpha diversity, we first normalized the taxa abundance using the “relative” option in the normalize_data function from the “microbiome” R package [29]. Then the plot_anova-diversity function from the same package was used with the adjusted level of significance set at ≤ 0.05 for species richness (observed OTU count), the Simpson index and the Shannon index [30]. Pairwise comparison between treatment and each alpha index was done using the non-parametric Wilcoxon test with a significant p-value set at < 0.05 [31]. We also explored four beta diversity indices, i.e. Bray-Curtis dissimilarity, Jensen-Shannon divergence (JSD), double principal coordinate analysis (DPCoA), and weighted UniFrac distances, to characterize the impact of treatment on microbiota structure [32-34]. Each of these beta diversity indices exploits slightly different aspects of the microbiota structure and a comparative examination of these tools allowed us to determine the extent to which community structural variations are explained by treatment or sex.
Treatment or sex associated microbial structural changes
To assess multivariate associations of experimental factors and microbial structural changes, a permutation multivariate analysis of variance (PerMANOVA) model with the Adonis function (9,999 permutations) was used in the “microbiome” package in R. Here a model was run for each of the beta diversity indices with the aim to assess the individual and cumulative variation explained by the factors (sex, treatment) used in this study.
Taxonomic profiles of ileal microbiota
The “metacoder” package (v0.3.2) in R [35] was used to profile the ileal microbial taxonomy, and taxonomic assignments (pan-phylogenetic) and their abundances (heat tree) were plotted based on treatment. This allowed for visualization of changes in abundance of each phylogenetic tree branches in each treatment group. In addition, mean abundances per family were computed and visualized using the “ComplexHeatmap” package (v2.5.1) as described [36]. The assumption for this analysis was if the hierarchical clustering of taxa recapitulates the treatment groups, then there should be a strong relationship between the taxonomic profiles and treatment groups [37].
Comparison between L. intracellularis clinical characteristics and ileal microbiota
To identify any relationship between the microbiota characteristics and the clinical signs of infection, i.e. histopathology and fecal shedding of L. intracellularis [24], we used a correlation analysis and the multinomial logistic regression to fit a model using the “ggbupr” package (v0.4.0) in R. First, we assessed the correlation between shedding and histopathology, then examined how each is associated with primary microbiota parameters such as abundance, diversity and richness as well as secondary parameters such as clustering. Parameters and characteristics were defined as directly or inversely correlated and a statistical significance of p<0.05 determined the significance of the association. This allowed us to identify treatment-mediated parameters associated with reduced shedding and histopathology