Exploring the Gastrointestinal “Nemabiome”: Deep Amplicon Sequencing to Quantify the Species Composition of Parasitic Nematode Communities

<div><p>Parasitic helminth infections have a considerable impact on global human health as well as animal welfare and production. Although co-infection with multiple parasite species within a host is common, there is a dearth of tools with which to study the composition of these complex parasite communities. Helminth species vary in their pathogenicity, epidemiology and drug sensitivity and the interactions that occur between co-infecting species and their hosts are poorly understood. We describe the first application of deep amplicon sequencing to study parasitic nematode communities as well as introduce the concept of the gastro-intestinal “nemabiome”. The approach is analogous to 16S rDNA deep sequencing used to explore microbial communities, but utilizes the nematode ITS-2 rDNA locus instead. Gastro-intestinal parasites of cattle were used to develop the concept, as this host has many well-defined gastro-intestinal nematode species that commonly occur as complex co-infections. Further, the availability of pure mono-parasite populations from experimentally infected cattle allowed us to prepare mock parasite communities to determine, and correct for, species representation biases in the sequence data. We demonstrate that, once these biases have been corrected, accurate relative quantitation of gastro-intestinal parasitic nematode communities in cattle fecal samples can be achieved. We have validated the accuracy of the method applied to field-samples by comparing the results of detailed morphological examination of L3 larvae populations with those of the sequencing assay. The results illustrate the insights that can be gained into the species composition of parasite communities, using grazing cattle in the mid-west USA as an example. However, both the technical approach and the concept of the ‘nemabiome’ have a wide range of potential applications in human and veterinary medicine. These include investigations of host-parasite and parasite-parasite interactions during co-infection, parasite epidemiology, parasite ecology and the response of parasite populations to both drug treatments and control programs.</p></div

The determination of species proportions of a series of pairwise parasite species combinations using the amplicon sequencing assay.

<p>DNA lysates were prepared from mock pools of 2000 L3 larvae comprising varying proportions of two different nematode species and the amplicon sequencing assay applied. Panel 3A, <i>O</i>. <i>ostertagi</i> and <i>H</i>. <i>placei</i>; Panel 3B, <i>C</i>. <i>oncophora</i> and <i>O</i>. <i>ostertagi</i>; Panel 3C, <i>O</i>. <i>ostertagi</i> and <i>C</i>. <i>punctata</i>; Panel 3D, <i>C</i>. <i>oncophora</i> and <i>C</i>. <i>punctata</i>. The left hand chart of each pair shows the expected results based on the known numbers of larvae added to each pool and the right hand chart of each pair shows the results of the amplicon sequencing assay (after the application of the appropriate correction factors). The ratios by which the two species of each pair vary from left to right on each chart are as follows: Column 1, 99:1; Column 2 90:10; Column 3; 70:30, Column 4, 50:50; Column 5; 30:70, Column 6; 10:90 and Column 7; 1:99.</p

Neighbor-joining tree of ITS-2 rDNA sequences obtained from 8 different parasitic nematode species from cattle.

<p>ITS-2 rDNA amplicons were directly sequenced by Sanger sequencing from 15 individual L3 larvae obtained from experimentally passaged strains for each of eight cattle parasite species. The neighbor-joining tree (Jukes-Cantor model) was computed with 2000 bootstrap replicates and rooted on <i>C</i>. <i>elegans</i> ITS-2 (Accession: JN636100) using Geneious version 7.1.5 created by Biomatters. Available from <a href="http://geneious.com/" target="_blank">http://geneious.com/</a>. Bootstrap values >70 are shown, values >80 are considered to be significant.</p

Comparison of the amplicon sequencing assay with visual morphological identification of larval populations from field samples.

<p>Fecal samples were collected from 39 individual calves entering feedlots from pasture at different locations in Oklahoma, Arkansas and Nebraska. 100 L3 larvae harvested from each coproculture were fixed and identified to the species level on the basis of morphological features on microscopy. The amplicon sequencing assay was applied to several hundred to several thousand L3 larvae from each sample. <b>Fig 6A</b>. Two bars are shown for each sample on the chart. The left bar of each sample shows the species proportions as determined by larval morphology. The right bar of each sample shows the proportions of each species as estimated by the amplicon sequencing assay (values shown are after application of the appropriate correction factor). The Y-axis shows the percentage proportions of each species as determined by morphology or the sequencing assay respectively. <b>Fig 6B</b>. Linear regression analysis was performed to assess the correlation between the percentage composition of each parasite species in each of 39 field samples as determined by the amplicon sequencing assay and morphological identification. The linear regression plots were produced using SPSS Statistics (IBM Corp. Released 2012. IBM SPSS Statistics for Macintosh, Version 21.0. Armonk, NY: IBM Corp), by plotting the percentage representation of each species from the morphological data (Y-axis) against the percentage representation of each species from the amplicon sequencing data (X-axis) for each sample. The values for R² (coefficient of determination), b (y-intercept) and m (slope) are shown adjacent to each plot. If values are the same the y-intercept (b) will be zero and the slope (m) will be 1.</p

Detection threshold of the amplicon sequencing assay.

<p>Separate mock pools of ~2000 larvae were created consisting of equal proportions of 6 parasite species. The seventh species was then added in increasing numbers of 0, 1, 5 and 20 larvae respectively. The charts show the species proportions estimated from the results of the amplicon sequencing assay (following application of the appropriate correction factor). In each case, the Y-axis shows the percentage proportions of each species and the X-axis indicates the number of <i>H</i>. <i>placei</i> or <i>O</i>. <i>ostertagi</i> larvae added as well as the numbers of reads assigned to <i>H</i>. <i>placei</i> and <i>O</i>. <i>ostertagi</i> for <b>Fig 4A</b> and <b>4B</b> respectively. <b>Fig 4A</b>. Results when <i>H</i>. <i>placei</i> was added in increasing numbers to a pool comprising an equal number of <i>N</i>. <i>helvetianus</i>, <i>C</i>. <i>oncophora</i>, <i>T</i>. <i>axei</i>, <i>T</i>. <i>colubriformis</i>, <i>C</i>. <i>punctata</i> and <i>O</i>. <i>ostertagi</i>. <b>Fig 4B</b>. Results when <i>O</i>. <i>ostertagi</i> was added in increasing numbers to a pool comprising an equal number of shows <i>N</i>. <i>helvetianus</i>, <i>C</i>. <i>oncophora</i>, <i>T</i>. <i>axei</i>, <i>C</i>. <i>punctata</i>, <i>T</i>. <i>colubriformis</i> and <i>H</i>. <i>placei</i>.</p

Repeatability of the amplicon sequencing assay when applied to field samples.

<p>The amplicon sequencing assay was applied to DNA lysates made from populations of L3 larvae collected from field samples from six Canadian cattle farms (Field samples 1–6). <b>Fig 5A</b>: The results of three technical replicates of the amplicon sequencing assay applied to the same single DNA lysate made from larval populations of each farm. Each DNA lysate was made from 300–2000 larvae (Field sample 1 = 1000 L3; Field sample 2 = 400 L3; Field sample 3 = 1000 L3; Field sample 4 = 2000 L3; Field sample 5 = 300 L3 and Field sample 6 = 2000 L3. <b>Fig 5B</b>: The results of three technical replicates of the amplicon sequencing assay applied to three independent DNA lysates made from separate batches of 1,000 and 2,000 larvae from field samples 1 and 4 respectively. The Y-axis of both charts shows the percentage species proportions as determined by the amplicon sequencing assay (25 cycles of amplification and following application of the appropriate correction factor). The field sample number is indicated on the X-axis.</p

Differential expression of 18 fasting response genes*, in the 1 µg/ml (1.1 µM) IVM exposure microarray experiment.

<p>In general genes that were shown to be induced by fasting were also induced following exposure to IVM and fasting repressed genes were also repressed by IVM exposure. *Taken from van Gilst <i>et al.</i> 2005 <a href="http://www.plosone.org/article/info:doi/10.1371/journal.pone.0031367#pone.0031367-VanGilst1" target="_blank">[24]</a>.</p

GFP expression patterns of IVM responsive genes.

<p><i>mtl-1</i>, <i>scl-2</i>, <i>cyp-37B1</i> and <i>ilys-3</i> are expressed in the intestine and C23G10.11 is expressed in the hypodermal cells. These are both proposed sites of fat storage.</p

Top 10 up-regulated probesets based on fold change following 4 hr exposure of DA1316 to 100 ng/ml (0.1 µM) IVM.

<p>Top 10 up-regulated probesets based on fold change following 4 hr exposure of DA1316 to 100 ng/ml (0.1 µM) IVM.</p

Pharyngeal pumping rate following 4 hr exposure of DA1316 and N2 to 0–1000 ng/ml (1.14 µM) IVM.

<p>Whilst strain DA1316 is more resistant to IVM induced pharyngeal paralysis, at concentrations greater than 100 ng/ml (0.1 µM) the pharyngeal pumping rate is significantly reduced. The pharyngeal pumping rates of five worms were counted at each concentration of IVM and the error bars represent the standard deviation.</p