Skip to main content

Advertisement

Transcriptome analysis unraveled potential mechanisms of resistance to Haemonchus contortus infection in Merino sheep populations bred for parasite resistance

Article metrics

  • 1181 Accesses

  • 2 Citations

Abstract

Haemonchus contortus is one of the most pathogenic gastrointestinal nematodes in small ruminants. To understand molecular mechanisms underlying host resistance to this parasite, we used RNA-sequencing technology to compare the transcriptomic response of the abomasal tissue, the site of the host-parasite interaction, of Merino sheep bred to be either genetically resistant or susceptible to H. contortus infection. Two different selection flocks, the Haemonchus selection flock (HSF) and the Trichostrongylus selection flock (TSF), and each contains a resistant and susceptible line, were studied. The TSF flock was seemingly more responsive to both primary and repeated infections than HSF. A total of 127 and 726 genes displayed a significant difference in abundance between resistant and susceptible animals in response to a primary infection in HSF and TSF, respectively. Among them, 38 genes were significantly affected by infection in both flocks. Gene ontology (GO) enrichment of the differentially expressed genes identified in this study predicted the likely involvement of extracellular exosomes in the immune response to H. contortus infection. While the resistant lines in HSF and TSF relied on different mechanisms for the development of host resistance, adhesion and diapedesis of both agranulocytes and granulocytes, coagulation and complement cascades, and multiple pathways related to tissue repair likely played critical roles in the process. Our results offered a quantitative snapshot of changes in the host transcriptome induced by H. contortus infection and provided novel insights into molecular mechanisms of host resistance.

Introduction

Gastrointestinal nematodes (GIN) represent a major health issue for livestock production systems worldwide [1]. GIN infection causes serious losses to farmers, both in impaired production and in control with anthelmintics. Anthelmintics are rapidly becoming ineffective against economically important GIN, due to the increasing incidence of anthelmintic resistance [2]. Moreover, current reliance on anthelmintics has resulted in public concerns for animal welfare and the contaminating residues in animal products [3]. Together, these factors have driven the development of effective and sustainable control strategies, including the application of novel vaccines that can be used to increase flock resistance and the development of genetically resistant populations [4, 5].

Immunological protection against GIN infection is associated with T-helper (Th) responses as well as tissue repair systems. Compared with susceptible animals, which usually evoke a response that shows characteristics of a Th1 response, resistant animals typically produce a dominant, polarized Th2 immune response [6, 7]. Therefore, a hallmark of immunity to GIN is a strong Th2 immune response [8, 9], characterized by the production of cytokines interleukin 4 (IL4), IL5, and IL13 [10]. In sheep resistant to GIN, the characteristics of a Th2 host response include the production of elevated levels of parasite-specific immunoglobulin A (IgA), IgE and IgG1 antibodies [11], eosinophilia [12], mucosal mastocytosis, and goblet cell hyperplasia [13]. However, the molecular mechanism of host resistance to H. contortus has yet to be fully understood.

In this study, we took advantage of the availability of a valuable sheep model that includes two resource flocks developed by selective breeding for genetic resistance or susceptibility to H. contortus [14] and Trichostrongylus colubriformis [15] infections. We used a high-throughput RNA sequencing (RNA-seq) approach in an attempt to identify features in the host transcriptome using both sets of flocks in response to an experimental H. contortus challenge. As such, the work reported here is considered as a companion to and complementary with previous studies [9, 16, 17] in which the candidate gene approach, microarray technology, and proteomic analysis were carried out using the same model. While these prior studies have led to valuable findings and resulted in some working hypotheses for future research, they are largely based on either a limited number of candidate genes or cross-species comparisons (ovine to bovine for microarrays and ovine to human for proteomics). A holistic analysis of the abomasal transcriptome in these resource populations using high-throughput RNA-seq technology and bioinformatics tools, especially on novel transcriptome features, such as transcript isoforms and alternative splicing, will facilitate our understanding of molecular mechanisms of resistance to H. contortus infection.

Materials and methods

Animals and parasitology

All lambs were obtained from two resource flocks of Merino sheep, the Haemonchus selection flock (HSF) and the Trichostrongylus selection flock (TSF), developed by selection and assortative mating for exploring the genetic basis of host resistance to nematodes [16, 18]. Each selection flock contains a resistant (R) line and a susceptible (S) line. The lambs in this study were reared in an animal housing facility on wooden slats from birth. Prior to the commencement of the challenge experiment by H. contortus, lambs from each line were divided into “innate” (I) and “acquired” (A) experimental infection groups. There were between six and seven lambs in each of these groups (Table 1). Lambs in the innate groups (I) received a primary challenge, whereas those in the acquired groups (A) received tertiary challenges. The acquired groups were challenged twice at 16 and 27 weeks of age, respectively, with 5000 H. contortus L3 larvae each (the Kirby strain). The infections were allowed to last for 6 weeks before being terminated using a single dose of 8.25 mg/kg levamisole hydrochloride treatment. At 37 weeks of age, a final infection (3rd) for the acquired groups and a primary (1st) infection for the innate groups were administrated with 10 000 H. contortus L3 larvae. All animals were euthanized on the third day after this infection. All animal procedures were carried out according to the protocols approved by the New South Wales (NSW) Department of Primary Industries, Director General’s Animal Care and Ethics Committee; and Australian Animal Welfare Standards and Guidelines for sheep were strictly followed (Animal Protocol# 03/75, 04/49 and 05/58).

Table 1 Groups of sheep used during Haemonchus contortus challenge trials

RNA extraction and sequencing using RNA-seq technology

Abomasal tissues were scraped and snap frozen in liquid nitrogen and then stored at −80 °C until use. Total RNA samples from the abomasal scrapings were isolated using Trizol (Invitrogen, Carlsbad, CA, USA) followed by DNase digestion and Qiagen RNeasy column purification (Qiagen, Valencia, CA, USA). RNA concentrations were determined using a NanoDrop ND-1000 spectrophotometer (Thermo Scientific, Wilmington, MA, USA), and the RNA integrity was verified using an Agilent BioAnalyzer 2100 (Agilent, Palo Alto, CA, USA) with a RNA integrity number (RIN) value > 7.5. High-quality RNA was processed using an Illumina TruSeq RNA sample prep kit (Illumina, San Diego, CA, USA), following the manufacturer’s instructions. The concentrations of individual RNA-seq libraries were verified and then pooled according to their respective index barcodes. Pooled RNA-seq libraries were sequenced as paired-end reads at 51 bp/sequence using an Illumina HiSeq 2000 sequencer (Illumina), as previously described [19, 20]. The metadata and raw sequence files were deposited in the National Center of Biotechnology Information (NCBI) Sequence Read Archive (SRA access number PRJNA445172).

Data analysis and bioinformatics

Raw sequence reads (mean ± SD = 65 699 458 ± 18 187 533 per sample, N = 50) were first checked using FastQC (Babraham Institute, Cambridge, UK). Low-quality reads were discarded; and low-quality nucleotides of each raw sequence read were trimmed using Trimmomatic [21] with default parameters. The resultant quality reads (mean reads = 57 921 126 per sample) after cleansing were mapped to the sheep reference genome Oar_v3.1 using TopHat2 (version 2.1.1) [22]. TopHat2 was run with the ‘‘–no coverage search’’ option and transcriptome indices were provided via the ‘‘–transcriptome-index’’ option. All other parameters were default settings. Transcripts were assembled and quantified using Cufflinks (version 2.2.0) [23], which was run with “-u” (multi-read correction) and “-b” (bias-correction) options. Transcript assemblies generated by Cufflinks were merged into a single transcriptome annotation using Cuffmerge (version 2.1.1) [23]. The reference annotation file and genomic sequence were provided to Cuffmerge via the “-g” and “-s” options, respectively. Differential expression analysis was performed using Cuffdiff (version 2.2.1) with the ‘‘-u’’ (multi-read correction) and ‘‘-b’’ (bias-correction) options. Differentially expressed genes (DEG) were extracted with the aid of the cummeRbund (R package), using a false discovery rate (FDR) < 0.05 as a cut-off and including an additional fold change filter (> 1.5-fold). The comparison was made independently between resistant and susceptible lines in each flock in response to different infection protocols. In addition, the raw data were also analyzed independently using the STAR-RSEM-EdgeR pipeline as previously described [24] and the CLC Genomics Workbench (Qiagen). Raw reads were mapped to both the sheep genome assemblies Oar_v3.1 and Oar_v4.0 using the ultrafast STAR aligner [25] and the CLC Genomics Workbench, respectively. DEG were further analyzed with Ingenuity Pathways Analysis (IPA) software (Qiagen), essentially as described previously [26]. Gene ontology (GO) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) assignments were analyzed with Database for Annotation, Visualization and Integrated Discovery (DAVID v6.8) [27]. In addition, gene enrichment analysis was performed using ShinyGO v0.50. Principal component analysis (PCA) was conducted using Primer v7.0 (PRIMER-E, Plymouth, UK).

Validation by real-time reverse transcription polymerase chain reaction (qRT-PCR)

The relative expression of 12 genes (see Additional file 1 for their primer sequences) was determined by qRT-PCR, as previously described [19]. Ovine ribosomal protein L19 gene (RPL19), whose expression remained stable among the experimental samples, was used as an endogenous reference gene for all reactions. cDNA was synthesized from total RNA using an iScript cDNA synthesis kit (BIO-RAD, Hercules, CA, USA), according to the manufacturer’s instructions. All qRT-PCR reactions were carried out in a Mic Real-Time PCR Cycler (Bio Molecular Systems, Australia) and analyzed with micPCR v2.2 software (Bio Molecular Systems, Upper Coomera, QLD, Australia). The reactions were run in triplicate in a total volume of 20 μL containing the following: 2 μL of cDNA (200 ng), 0.5 μL of each primer (forward and reverse, 20 nM each), 10 μL of 2 × SensiFAST No-ROX Mix (Bioline, MA, USA) and 7 μL of nuclease-free water. The amplification reactions were subjected to an initial denaturation at 95 °C for 5 min, followed by 40 cycles of 95 °C for 10 s, 60 °C for 30 s, and 72 °C for 20 s. Melting curves were obtained from 60 to 95 °C. Relative gene expression data was calculated using the 2−ΔΔCT method.

Results

Differences in the abomasal transcriptome of two selection flocks in response to experimental H. contortus infections

In this study, approximately 76.40% of raw reads (± 4.05%, SD) were uniquely mapped to the ovine genome assembly (Oar_v3.1) using TopHat2. PCA analysis shows that the separation between the two flocks, between the two lines (Additional file 2), or among the various experimental groups (Figure 1) was indistinguishable, suggesting that neither long-term selective breeding nor different experimental H. contortus infection protocols had a profound effect on the overall abomasal transcriptome structure and composition. Nevertheless, the extent of the changes in the abomasal transcriptome of different lines in different flocks in response to differential infection protocols differed. The numbers of DEG identified between various groups are shown in Figure 2. We identified 127 and 726 DEGs in the innate infection groups between the R and S lines within the HSF and TSF flocks, respectively, while the number of DEG identified in the acquired groups was 19 and 378 between the R and S lines within the HSF and TSF flocks, respectively. It appeared that greater numbers of DEG were identified in the groups using the innate infection protocols than those using the acquired protocol (Additional files 3, 4, 5, 6), suggesting the abomasal transcriptome of the lambs may be more responsive to a primary infection, in both flocks. Moreover, the abomasal transcriptome of the animals in the TSF flock was seemingly more responsive than those in the HSF flock, regardless of the infection protocols used. The expression of 38 genes was significantly (FDR < 0.05) impacted by H. contortus infection in both flocks when R lines were compared to their respective S lines in the innate infection (Table 2). These DEG likely represented a set of early response genes to H. contortus primary infections. Among them, the expression of 9 genes, such as heat shock protein family A (Hsp70) member 6 (HSPA6), lectin, galactoside-binding, soluble, 15 (LGALS15), and complement factor I (CFI), were significantly enhanced by infection in both R lines when compared to their respective S lines.

Figure 1
figure1

Principal component analysis (PCA) based on normalized raw counts. H: Haemonchus selection flock (HSF); T: Trichostrongylus selection flock (TSF); R: resistant line; S: susceptible line; I: innate infection protocol (primary infection); A: acquired infection protocol (repeated infection).

Figure 2
figure2

Venn diagrams and Volcano plots. A Venn diagram showing the number of genes with significant differences in relative abundance induced by infection point in two sheep flocks when comparing the resistant lines with their respective susceptible counterparts at a false discovery rate (FDR) cutoff < 0.05. B Volcano plots. The red color indicates genes detected as differentially expressed between the resistant groups and susceptible groups at FDR < 0.05. H: HSF flock; T: TSF flock; R: resistant line; S: susceptible line; I: innate infection protocol (primary infection); A: acquired infection protocol (repeated infection).

Table 2 38 genes displayed significant differences (false discovery rate FDR < 0.05) in the relative abundance in response to a primary infection in both sheep Haemonchus selection (HSF) and Trichostrongylus selection (TSF) flocks selected for parasite resistance

Three and 17 genes significantly impacted by infection in response to the innate protocol between the R and S lines are related to extracellular matrix (ECM) within the HSF and TSF flocks, respectively (Additional files 3 and 4). Of note, all of them were significantly upregulated in the R lines when compared to their respective S lines. TFF3, involved in wound healing, mucosal protection, cell proliferation and cell migration, was strongly upregulated (21.54-fold) by infection in the R line of the HSF flock when compared to its S counterparts. Moreover, the mRNA expression of several cytokine receptors and chemokines were significantly impacted by infection (FDR < 0.05). The transcripts of IL15 receptor subunit alpha (IL15RA) and IL22 receptor subunit alpha 1 (IL22RA1) were 2.50- and 3.50-fold higher in R animals than S animals in response to the innate challenge protocol in the TSF flock (Additional file 3). Notably, the abundance of seven chemokines, such as C–C motif chemokine ligand 14 (CCL14), CCL20, CCL25, CCL5, C–X–C motif chemokine ligand 10 (CXCL10), CXCL14, and atypical chemokine receptor 3 (ACKR3), was significantly higher by infection in the R line of the TSF flock (Additional file 3). Similarly, CCL19, CCL22, and CXCL13 were significantly upregulated by a primary infection in the R line of the HSF flock when compared to its S line counterparts (Additional file 4).

Of note, greater than 20% of the genes significantly impacted by infection in both flocks are extracellular exosome-related (Additional file 7). The expression of these extracellular exosome-related genes was predominantly enhanced in the R lines compared to their respective S lines. At least 121 genes were significantly upregulated by infection in the R line of the TSF flock, such as adiponectin, C1Q and collagen domain containing (AdipoQ), annexin A1 (ANXA1) and ANXA 3, various complement related genes (CFI, C1QA, C1QB, C1QC, and C1S), and integrin subunit alpha 1 (ITGA1) and ITGA3. On the other hand, 24 genes were significantly upregulated by infection in the R line of the HSF flock, such as CD19, CD22, CD79B, two complement-related genes CFI and CR2 (complement C3d receptor 2), FGG, and FGA (fibrinogen alpha chain). Nevertheless, the infection was also able to repress extracellular exosome-related genes, such as complement C7 (C7) in the R line of the HSF flock, and fibroblast growth factor 9 (FGF9), cadherin 16 (CDH16), and FGG, in the R line of the TSF flock.

Gene ontology (GO), pathways, and molecular networks related to host resistance

Gene enrichment analysis using ShinyGO identified a total of 14 GO terms shared by the two flocks in response to the innate infection protocols (FDR < 0.0001), including 5 Biological process (BP), 8 Cellular components (CC), and 1 Molecular function terms. For example, GO terms, such as Extracellular space, Fibrinogen complex, Response to external stimulus, and Signaling receptor binding, were significantly enriched in both comparisons, TRI vs TSI and HRI vs HSI. However, only 1 GO term, Extracellular space (GO:0005615), was significantly enriched in both flocks in response to the acquired infection protocol. Moreover, Inflammatory response (GO:0006954), Integrin-mediated signaling pathway (GO:0007229), Cell adhesion (GO:0007155) and cell matrix adhesion (GO:0007160) were seemingly among important biological processes involved in the development of resistance to a primary infection. In addition, the cellular component Extracellular exosome (GO:0070062) may play an important role in the development of host resistance, regardless of the flock origin and infection status (Additional file 7).

A total of 20 KEGG pathways were significantly impacted in the innate immune response to infection between the R and S lines in the TSF flock (FDR < 0.1; see Additional file 8). These pathways, such as leukocyte transendothelial migration, MAPK signaling, TNF signaling, cell adhesion molecules, focal adhesion (Table 3), and ECM-receptor interaction (Table 4), were likely involved in the development of host resistance to a primary infection in TSF. On the other hand, only two pathways, Complement and coagulation cascades, and B cell receptor signaling, were significantly enriched between the R and S lines in the HSF flock. Of note, Complement and coagulation cascades (Figure 3), was also significantly enriched in the TSF flock, suggesting that this pathway likely contributed to the development of host resistance to a primary H. contortus infection. Similarly, IPA gene enrichment analysis identified several canonical pathways that were significantly enriched in both flocks in response to the acquired infection protocol. Among these, coagulation system was significantly enriched in both HSF and TSF flocks in response to the acquired infection protocol, in a good agreement with the DAVID results. On the other hand, Interferon signaling was enriched in the HSF flock while the pathway IL-17A signaling in gastric cells may play an important role in the TSF flock.

Table 3 The genes related to gene ontology (GO) term focal adhesion significantly enriched in the Trichostrongylus selection flock (TSF) in response to a primary Haemonchus infection
Table 4 The genes related to the extracellular matrix (ECM)-receptor interaction pathway significantly enriched in the resistant line of the Trichostrongylus selection flock (TSF) in response to a primary Haemonchus contortus infection
Figure 3
figure3

Biological pathways significantly enriched in response to a primary infection. The Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways Complement and coagulation cascades were significantly enriched at a FDR < 10% in both HSF and TSF flocks in response to a primary Haemonchus contortus infection. *Denotes genes that displayed a significantly different abundance between R and S lines in the TSF flock. ^Denotes genes that displayed a significantly different abundance between R and S lines in the HSF flock.

The resistant lines (R) in the both flocks likely utilized different molecular mechanisms for host resistance to a primary infection, as reflected by the difference in interaction patterns and composition in the inferred molecular networks (Figures 4, 5). In the HSF flock, the resistant line may rely on several B cell related canonical pathways to fight a primary H. contortus infection. At least 5 canonical pathways, such as B Cell development, B Cell receptor signaling, B Cell activating factor signaling, PI3K Signaling in B lymphocytes, and altered T Cell and B Cell signaling in arthritis, were significantly enriched in the R line of the HSF flocks. In addition, both agranulocyte and granulocyte adhesion and diapedesis pathways played important roles in this flock. Moreover, it appeared that Complement system as well as coagulation system also played an important role in resistance development in this flock. In contrast, the pathways such as Histamine (and Dopamine) degradation, fatty acid α-oxidation, and acute phase response signaling were seemingly critical in fighting against a primary infection in the TSF flock.

Figure 4
figure4

Visual representation of molecular interactions in the HSF flock. The molecular interactions for the genes with significantly different abundance between the resistant and susceptible lines in the HSF flock in response to a primary Haemonchus contortus infection for 3 days were displayed. The network was inferred using Ingenuity Pathway Analysis (IPA v01-13). Solid lines imply direct relationships between gene products while dotted lines imply indirect interactions. Red: the genes with a significantly higher abundance in the Resistant line than in the susceptible line. Green: the genes with a significantly lower abundance in the Resistant line than in the susceptible line. The degree of changes in gene abundance is represented by the intensity of the color.

Figure 5
figure5

Visual representation of molecular interactions in the TSF flock. The network was inferred using Ingenuity Pathway Analysis (IPA v01-13). The genes in the network included those with significantly different abundance between the resistant and susceptible lines in the TSF flock in response to a primary Haemonchus contortus infection for 3 days. Solid lines imply direct relationships between gene products while dotted lines imply indirect interactions. Red: the genes with a significantly higher abundance in the Resistant line than in the susceptible line. Green: the genes with a significantly lower abundance in the Resistant line than in the susceptible line. The degree of changes in gene abundance is represented by the intensity of the color.

We also analyzed the dataset using the STAR-RSEM-EdegR pipeline. Under the same stringent threshold cutoff (FDR < 0.05), we identified 277 unique genes that displayed a significantly different abundance in one or more comparisons (Additional file 9). Among them, 185 genes, representing approximately 67% of all DEG identified by the STAR-RSEM-EdgeR pipeline, were also detected using the TopHat2–Cufflink–Cuffdiff pipeline. For example, significant changes in the expression of amphiregulin, C–C motif chemokine 25 (CCL25), galectins (LGALS3BP, LGALS7, and LGALS9), granzyme B, prostaglandin G/H synthase 1 (PTGS1 or COX1) and, prostaglandin-endoperoxide Synthase 2 (PTGS2 or COX2) were identified by both pipelines when comparing the TRI group with the TSI group. Moreover, the TopHat2 pipeline also identified transcript isoforms that displayed a significant difference in relative abundance between R and S lines of both flocks in response to the innate infection protocol (Additional file 10). Amphiregulin is a well-known Th2 cytokine enhancing resistance to gastrointestinal nematodes [28]. On the other hand, both pipelines detected a significant upregulation of cadherin 26 (CDH26) only in the resistant line in the HSF flock in response to the innate protocol.

Real-time RT-PCR confirmation

The RNA-seq results of selected genes were validated in the same sample set by real-time (q) RT-PCR (Figure 6). A correlation coefficient (R) of log2 transformed fold change values between qRT-PCR and RNA-seq platforms was 0.7749 (Figure 7). The expression results of CFI, CXCL13, IFIT1, IL22RA1, LGALS15, LPL, MMP1, and MMP13 obtained using qRT-PCR were in good agreement with the RNA-seq results.

Figure 6
figure6

Real-Time RT-PCR analysis (qPCR) of selected genes. Relative expression levels calculated from standard curves were normalized to the endogenous control RPL19 gene. Numbers represent mean plus standard error. ALDH1L2: aldehyde dehydrogenase 1 family member L2; CCL19: C–C motif chemokine ligand 19; CFI: complement factor I; CXCL13: C–X–C motif chemokine ligand 13; IFIT1: interferon induced protein with tetratricopeptide repeats 1; IL22RA1: interleukin 22 receptor subunit alpha 1; LGALS15: lectin, galactoside-binding, soluble, 15; LPL: lipoprotein lipase; MMP1: matrix metallopeptidase 1; MMP13: matrix metallopeptidase 13; TFF2: trefoil factor 2; TFF3: trefoil factor 3. *P < 0.05; **P < 0.01.

Figure 7
figure7

Linear regression analysis of fold changes calculated from qPCR and RNAseq analysis. Blue dots represent log2 transformed fold change values of a single gene in an infected sample obtained from qPCR (X-axis) and RNAseq analysis (Y-axis). Dashed lines: 99% confidence interval. R: correlation coefficient.

Discussion

It has long been recognized that differences in host resistance and susceptibility to parasitic infection exist in various sheep breeds [7, 29] and that genetics may play an important role in regulating host resistance, which has spurred efforts to control parasitic infection through selective breeding for naturally resistant sheep. Over 40 years, two flocks, HSF and TSF, each containing a resistant and a susceptible line, were established by selective breeding [9, 16, 18]. The two flocks have been shown to be ideal models for elucidating the molecular mechanisms of immune resistance and susceptibility of sheep to gastrointestinal nematodes [18].

In this study, the number of genes differentially expressed between resistant and susceptible lines in the HSF in response to both primary and tertiary Haemonchus infections was significantly fewer than those in the TSF (Figure 2). This may support a notion that long-term selective breeding has an effective impact, at least in part, on the immune response of Merino sheep to H. contortus infection. Moreover, our findings provided evidence that different mediators and pathways underlying resistance and susceptibility existed in the two selection flocks, likely due to the difference in their genetic backgrounds. Furthermore, it appeared that the abomasal transcriptome was more responsive to a primary infection than a repeated infection, regardless of the selection flock, as evidenced by a significantly higher number of DEG identified between the R and S lines in response to the innate infection protocol than the acquired infection in both flocks. Following a primary infection, 127 and 726 genes were identified as differentially expressed in the HSF and TSF flocks, respectively. Of these, 38 were shared in both flocks in response to a primary infection, whereas none were identified in both flocks after the tertiary challenge. A proportion of the shared genes differentially expressed following a primary infection have previously been associated with helminth infection. These genes, such as CFI, IFIT1, IFIT2, and LGALS15, may underlie the basis of resistance. CFI is the major complement inhibitor that degrades the activated C3b and C4b in the presence of specific cofactors. CFI deficiency results in secondary complement C3 deficiency due to uncontrolled spontaneous alternative pathway activation, and CFI-deficient individuals suffer from recurring infections [30, 31]. CFI had significantly more abundant transcripts in the abomasum of resistant sheep, which is also observed in the resistant Canaria Hair breed [20]. LGALS15 has been shown to be absent from the normal gastrointestinal epithelium in sheep but induced and secreted into the mucus by epithelial cells within three days after a primary infection with H. contortus [32]. The mRNA levels of LGALS15 were significantly up-regulated in the abomasum of resistant sheep in both flocks after a primary infection in this study. In addition, following a primary infection, the LGALS3 gene in the HSF flock and the LGALS3BP gene in the TSF flock were significantly induced in the abomasum of resistant sheep. IFIT genes, induced by host innate immune defenses after pathogen infection, have been suggested to function as antiviral effectors and immunomodulators [33]. IFIT proteins inhibits viral infections by suppressing translation initiation, binding uncapped or incompletely capped viral RNA, sequestering viral proteins or RNA in the cytoplasm, and regulating cell-intrinsic and cell-extrinsic immune responses through pathways that remain to be defined and/or corroborated. Interestingly, following a primary infection, IFIT1 and IFIT2 genes were repressed more in the resistant animals than in the susceptible animals in the HSF flock, while all members of IFIT family (IFIT1, IFIT2, IFIT3, and IFIT5) were significantly upregulated in the resistant line compared to the susceptible line in the TSF flock (Table 2). The divergent pattern of IFIT expression may contribute to different mechanisms utilized by different selection flocks in the acquisition of immunity to H. contortus infection. IFIT proteins have not been previously associated with the development of immune responses to parasitic nematode infections. Further study is required to determine the expression profiles in different breeds and uncover their functions in the immune response to nematode infection.

Molecular interactions inferred using IPA provided further evidence that the mechanisms responsible for the development of host resistance in different flocks in response to two different infection protocols differed. The resistant lines in both flocks displayed a stronger and more complex interaction network in response to the innate infection protocol than the acquired protocol. The largest molecular network (Figure 4) in the comparison between HRI and HSI consisted of 20 focus molecules while the largest network in the TRI vs TSI included 29 focus molecules (Figure 5). The major molecular function associated with the former network was humoral immune response and immunological disease while the latter network included organismal injury and abnormalities. Indeed, in the HSF flock in response to the innate infection protocol, J chain and immunoglobulin heavy constant mu (IGHM) seemingly played an important role to the coherence of this network. Mucosal IgA is critical to host immune responses and the development of resistance to Haemonchus and Teladorsagia infections in sheep. It is negatively correlated with adult worm length and fecundity in a resistant sheep breed [34]. J chain regulates the formation of polymeric IgA and IgM and has a high affinity for the polymeric Ig receptor (PIGR), which may be involved in the development of parasite resistance in cattle [35]. In the bovine abomasum, PIGR, Complement C3, and J chain are among the most abundant transcripts [19]. Together, interactions among multiple molecules consisting of this network, such as CD19, CD 22, and CD79A/B, and chemokines including CCL19, CCL22, and CXCL13 may contribute to enhanced immune responses in the resistant line in the HSF flock to the innate infection. In the latter network (Figure 5), protein kinase B (also known as Akt) occupied a focal position in the network. In contrast to the compositional difference in both networks in response the innate infection protocol, in both flocks in response to the acquired infection protocol, cell cycle, cell death and survival, and cell signaling were among the primary functions in the networks identified (data not shown).

There were seemingly strong relationships between immune responses to parasite infections and the regulation of focal adhesion, extracellular matrix remodeling, and immune cell trafficking, which are involved in the reorganization of the epithelial layer and wound repair and play important roles in gut physiology [36]. Indeed, the rate of epithelial turnover has been shown to be associated with accelerated worm expulsion [16, 17, 37].

Intriguingly, approximately 20% of the differentially expressed genes after a primary infection in the resistant sheep in both flocks were extracellular exosome-related (Additional file 7). The majority of these exosome related genes were significantly upregulated in resistant sheep. Host-derived extracellular exosomes have been shown to be used as defense mechanisms and play important roles in antigen presentation during parasitic infection. For example, intestinal epithelial cells increased the release of antimicrobial peptide-containing exosomes in response to Cryptosporidium infection, which is driven by enhanced toll-like receptor 4 signaling following recognition of the protozoan parasite [38]. Mycobacterium tuberculosis is also able to induce exosome release from infected macrophages, which consequently promotes recruitment of lymphocytes through heightened proinflammatory chemokine secretion [39]. Mycobacterium bovis-infected macrophage-derived exosomes can promote dendritic cell activation and generate antibacterial T cells in vivo [40]. Importantly, vaccination of chickens with Eimeria parasite antigen-loaded dendritic cell exosomes has been shown to ameliorate symptoms of avian coccidiosis [41]. However, molecular mechanisms of extracellular exosomes in the development of host resistance to H. contortus infection in sheep, especially their roles in eosinophilia and mucosal mastocytosis, warrants further investigation.

References

  1. 1.

    McLeod RS (1995) Costs of major parasites to the Australian livestock industries. Int J Parasitol 25:1363–1367

  2. 2.

    Sangster NC, Cowling A, Woodgate RG (2018) Ten events that defined anthelmintic resistance research. Trends Parasitol 34:553–563

  3. 3.

    Waller PJ (2006) From discovery to development: current industry perspectives for the development of novel methods of helminth control in livestock. Vet Parasitol 139:1–14

  4. 4.

    Stear MJ, Doligalska M, Donskow-Schmelter K (2007) Alternatives to anthelmintics for the control of nematodes in livestock. Parasitology 134:139–151

  5. 5.

    Ekoja SE, Smith WD (2010) Antibodies from sheep immunized against Haemonchus contortus with H-gal-GP inhibit the haemoglobinase activity of this protease complex. Parasite Immunol 32:731–738

  6. 6.

    Wilkie H, Nicol L, Gossner A, Hopkins J, Mucosal (2016) Expression of T cell gene variants is associated with differential resistance to Teladorsagia circumcincta. PLoS One 11:e0168194

  7. 7.

    Alba-Hurtado F, Munoz-Guzman MA (2013) Immune responses associated with resistance to haemonchosis in sheep. Biomed Res Int 2013:162158

  8. 8.

    Lacroux C, Nguyen TH, Andreoletti O, Prevot F, Grisez C, Bergeaud JP, Gruner L, Brunel JC, Francois D, Dorchies P, Jacquiet P (2006) Haemonchus contortus (Nematoda:Trichostrongylidae) infection in lambs elicits an unequivocal Th2 immune response. Vet Res 37:607–622

  9. 9.

    Andronicos N, Hunt P, Windon R (2010) Expression of genes in gastrointestinal and lymphatic tissues during parasite infection in sheep genetically resistant or susceptible to Trichostrongylus colubriformis and Haemonchus contortus. Int J Parasitol 40:417–429

  10. 10.

    Paul WE, Zhu J (2010) How are T(H)2-type immune responses initiated and amplified? Nat Rev Immunol 10:225–235

  11. 11.

    Pernthaner A, Cole SA, Morrison L, Green R, Shaw RJ, Hein WR (2006) Cytokine and antibody subclass responses in the intestinal lymph of sheep during repeated experimental infections with the nematode parasite Trichostrongylus colubriformis. Vet Immunol Immunopathol 114:135–148

  12. 12.

    Rothwell TL, Windon RG, Horsburgh BA, Anderson BH (1993) Relationship between eosinophilia and responsiveness to infection with Trichostrongylus colubriformis in sheep. Int J Parasitol 23:203–211

  13. 13.

    Douch PG, Green RS, Morris CA, McEewan JC, Windon RG (1996) Phenotypic markers for selection of nematode-resistant sheep. Int J Parasitol 26:899–911

  14. 14.

    Woolaston RR, Barger IA, Piper LR (1990) Response to helminth infection of sheep selected for resistance to Haemonchus contortus. Int J Parasitol 20:1015–1018

  15. 15.

    Dawkins HJ, Windon RG, Outteridge PM, Dineen JK (1988) Cellular and humoral responses of sheep with different levels of resistance to Trichostrongylus colubriformis. Int J Parasitol 18:531–537

  16. 16.

    Ingham A, Reverter A, Windon R, Hunt P, Menzies M (2008) Gastrointestinal nematode challenge induces some conserved gene expression changes in the gut mucosa of genetically resistant sheep. Int J Parasitol 38:431–442

  17. 17.

    Nagaraj SH, Harsha HC, Reverter A, Colgrave ML, Sharma R, Andronicos N, Hunt P, Menzies M, Lees MS, Sekhar NR, Pandey A, Ingham A (2012) Proteomic analysis of the abomasal mucosal response following infection by the nematode, Haemonchus contortus, in genetically resistant and susceptible sheep. J Proteomics 75:2141–2152

  18. 18.

    Windon RG (1996) Genetic control of resistance to helminths in sheep. Vet Immunol Immunopathol 54:245–254

  19. 19.

    Li RW, Rinaldi M, Capuco AV (2011) Characterization of the abomasal transcriptome for mechanisms of resistance to gastrointestinal nematodes in cattle. Vet Res 42:114

  20. 20.

    Guo Z, Gonzalez JF, Hernandez JN, McNeilly TN, Corripio-Miyar Y, Frew D, Morrison T, Yu P, Li RW (2016) Possible mechanisms of host resistance to Haemonchus contortus infection in sheep breeds native to the Canary Islands. Sci Rep 6:26200

  21. 21.

    Bolger AM, Lohse M, Usadel B (2014) Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics 30:2114–2120

  22. 22.

    Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, Salzberg SL (2013) TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol 14:R36

  23. 23.

    Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, Pimentel H, Salzberg SL, Rinn JL, Pachter L (2012) Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat Protoc 7:562–578

  24. 24.

    Li W, Richter RA, Jung Y, Zhu Q, Li RW (2016) Web-based bioinformatics workflows for end-to-end RNA-seq data computation and analysis in agricultural animal species. BMC Genomics 17:761

  25. 25.

    Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, Batut P, Chaisson M, Gingeras TR (2013) STAR: ultrafast universal RNA-seq aligner. Bioinformatics 29:15–21

  26. 26.

    Van Meulder F, Van Coppernolle S, Borloo J, Rinaldi M, Li RW, Chiers K, Van den Broeck W, Vercruysse J, Claerebout E, Geldhof P (2013) Granule exocytosis of granulysin and granzyme B as a potential key mechanism in vaccine-induced immunity in cattle against the nematode Ostertagia ostertagi. Infect Immun 81:1798–1809

  27. 27.

    da Huang W, Sherman BT, Lempicki RA (2009) Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc 4:44–57

  28. 28.

    Zaiss DM, Yang L, Shah PR, Kobie JJ, Urban JF, Mosmann TR (2006) Amphiregulin, a TH2 cytokine enhancing resistance to nematodes. Science 314:1746

  29. 29.

    Stear MJ, Boag B, Cattadori I, Murphy L (2009) Genetic variation in resistance to mixed, predominantly Teladorsagia circumcincta nematode infections of sheep: from heritabilities to gene identification. Parasite Immunol 31:274–282

  30. 30.

    Botto M, Kirschfink M, Macor P, Pickering MC, Wurzner R, Tedesco F (2009) Complement in human diseases: lessons from complement deficiencies. Mol Immunol 46:2774–2783

  31. 31.

    Degn SE, Jensenius JC, Thiel S (2011) Disease-causing mutations in genes of the complement system. Am J Hum Genet 88:689–705

  32. 32.

    Dunphy JL, Balic A, Barcham GJ, Horvath AJ, Nash AD, Meeusen EN (2000) Isolation and characterization of a novel inducible mammalian galectin. J Biol Chem 275:32106–32113

  33. 33.

    Diamond MS, Farzan M (2013) The broad-spectrum antiviral functions of IFIT and IFITM proteins. Nat Rev Immunol 13:46–57

  34. 34.

    Hernandez JN, Hernandez A, Stear MJ, Conde-Felipe M, Rodriguez E, Piedrafita D, Gonzalez JF (2016) Potential role for mucosal IgA in modulating Haemonchus contortus adult worm infection in sheep. Vet Parasitol 223:153–158

  35. 35.

    Li RW, Sonstegard TS, Van Tassell CP, Gasbarre LC (2007) Local inflammation as a possible mechanism of resistance to gastrointestinal nematodes in Angus heifers. Vet Parasitol 145:100–107

  36. 36.

    Anderson JM, Van Itallie CM (2009) Physiology and function of the tight junction. Cold Spring Harb Perspect Biol 1:a002584

  37. 37.

    Cliffe LJ, Humphreys NE, Lane TE, Potten CS, Booth C, Grencis RK (2005) Accelerated intestinal epithelial cell turnover: a new mechanism of parasite expulsion. Science 308:1463–1465

  38. 38.

    Hu G, Gong AY, Roth AL, Huang BQ, Ward HD, Zhu G, Larusso NF, Hanson ND, Chen XM (2013) Release of luminal exosomes contributes to TLR4-mediated epithelial antimicrobial defense. PLoS Pathog 9:e1003261

  39. 39.

    Singh PP, Smith VL, Karakousis PC, Schorey JS (2012) Exosomes isolated from mycobacteria-infected mice or cultured macrophages can recruit and activate immune cells in vitro and in vivo. J Immunol 189:777–785

  40. 40.

    Giri PK, Schorey JS (2008) Exosomes derived from M. Bovis BCG infected macrophages activate antigen-specific CD4+ and CD8+ T cells in vitro and in vivo. PLoS One 3:e2461

  41. 41.

    del Cacho E, Gallego M, Lee SH, Lillehoj HS, Quilez J, Lillehoj EP, Sanchez-Acedo C (2012) Induction of protective immunity against Eimeria tenella, Eimeria maxima, and Eimeria acervulina infections using dendritic cell-derived exosomes. Infect Immun 80:1909–1916

Download references

Competing interests

The authors declare that they have no competing interests.

Authors’ contributions

RWL, PH, and AI conceived and designed the experiments. RZ, FL, CL, PH, and RWL performed experiments. RZ, FL, LZ, and RWL analyzed the data. RZ, AI, and RWL drafted the manuscript. All authors read and approved the final manuscript.

Acknowledgements

Mention of trade names or commercial products in this publication is solely for the purpose of providing specific information and does not imply recommendation or endorsement by the U.S. Department of Agriculture (USDA). The USDA is an equal opportunity provider and employer. RZ is supported by China Scholarship Council (No. [2014]5050).

Availability of data and materials

The raw RNA-seq data are available in NCBI Sequence Read Archive (SRA) (BioProject: PRJNA445172).

Ethics approval and consent to participate

The animal experiment was approved by and performed under the guidelines of the F. D. McMaster Animal Ethics Committee, CSIRO Livestock Industries.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Author information

Correspondence to Aaron Ingham or Robert W. Li.

Additional files

13567_2019_622_MOESM1_ESM.xlsx

Additional file 1. Primer sequences used for qPCR validation.

13567_2019_622_MOESM2_ESM.docx

Additional file 2. Principal component analysis (PCA) by flocks (A) and by lines (B).

13567_2019_622_MOESM3_ESM.xlsx

Additional file 3. Differentially expressed genes (DEGs) between the Resistant and Susceptible lines in the TSF flock. Genes with significantly different abundance between the Resistant and Susceptible lines in response to a primary Haemonchus contortus infection in the Trichostrongylus Selection Flock (TSF) at a False Discovery Rate (FDR) < 0.05%.

13567_2019_622_MOESM4_ESM.xlsx

Additional file 4. DEGs between the Resistant and Susceptible lines in the HSF flock. Genes with significantly different abundance between the Resistant and Susceptible lines in response to a primary Haemonchus contortus infection in the Haemonchus Selection Flock (HSF) at a False Discovery Rate (FDR) < 0.05%.

13567_2019_622_MOESM5_ESM.xlsx

Additional file 5. DEGs between the two lines in the HSF flock in response to a tertiary infection. Genes with significantly different abundance between the Resistant and Susceptible lines in response to a Haemonchus contortus tertiary infection in the Haemonchus Selection Flock (HSF) at a False Discovery Rate (FDR) < 0.05%.

13567_2019_622_MOESM6_ESM.xlsx

Additional file 6. DEGs between the two lines in the TSF flock in response to a tertiary infection. Genes with significantly different abundance between the Resistant and Susceptible lines in response to a Haemonchus contortus tertiary infection in the Trichostrongylus Selection Flock (TSF) at a False Discovery Rate (FDR) < 0.05%.

13567_2019_622_MOESM7_ESM.xlsx

Additional file 7. Extracellular Exosome related genes significantly enriched in both flocks in response to Haemonchus contortus infections.

13567_2019_622_MOESM8_ESM.xlsx

Additional file 8. Significantly enriched Database for Annotation, Visualization, and Integrated Discovery (DAVID) terms. DAVID terms at a False Discovery Rate (FDR) < 10% in the Trichostrongylus Selection Flock (TSF) in response to a primary Haemonchus contortus infection were shown.

13567_2019_622_MOESM9_ESM.xlsx

Additional file 9. 277 unique differentially expressed genes detected using the STAR-EdgeR pipeline. 185 of the 277 genes, approximately, 67% of the all DEGs identified by the STAR pipeline are also detected using the same stringency cutoff by the Tophat2-Cufflink-Cuffdiff pipeline.

13567_2019_622_MOESM10_ESM.xlsx

Additional file 10. Significant transcript isoforms between resistant and susceptible lines in response to Haemonchus contortus infection. These isoforms were detected using the TopHat2-Cuffdiff pipeline.

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark