Skip to main content
Download PDF
  • Research article
  • Open access
  • Published:

Network analysis uncovers putative genes affecting resistance to tick infestation in Braford cattle skin

BMC Genomics volume 20, Article number: 998 (2019) Cite this article

Abstract

Background

Genetic resistance in cattle is considered a suitable way to control tick burden and its consequent losses for livestock production. Exploring tick-resistant (R) and tick-susceptible (S) hosts, we investigated the genetic mechanisms underlying the variation of Braford resistance to tick infestation. Skin biopsies from four-times-artificially infested R (n = 20) and S (n = 19) hosts, obtained before the first and 24 h after the fourth tick infestation were submitted to RNA-Sequencing. Differential gene expression, functional enrichment, and network analysis were performed to identify genetic pathways and transcription factors (TFs) affecting host resistance.

Results

Intergroup comparisons of hosts before (Rpre vs. Spre) and after (Rpost vs. Spost) tick infestation found 51 differentially expressed genes (DEGs), of which almost all presented high variation (TopDEGs), and 38 were redundant genes. Gene expression was consistently different between R and S hosts, suggesting the existence of specific anti-tick mechanisms. In the intragroup comparisons, Rpost vs. Rpre and Spost vs. Spre, we found more than two thousand DEGs in response to tick infestation in both resistance groups. Redundant and non-redundant TopDEGs with potential anti-tick functions suggested a role in the development of different levels of resistance within the same breed. Leukocyte chemotaxis was over-represented in both hosts, whereas skin degradation and remodeling were only found in TopDEGs from R hosts. Also, these genes indicated the participation of cytokines, such as IL6 and IL22, and the activation of Wingless (WNT)-signaling pathway. A central gene of this pathway, WNT7A, was consistently modulated when hosts were compared. Moreover, the findings based on a genome-wide association study (GWAS) corroborate the prediction of the WNT-signaling pathway as a candidate mechanism of resistance. The regulation of immune response was the most relevant pathway predicted for S hosts. Members of Ap1 and NF-kB families were the most relevant TFs predicted for R and S, respectively.

Conclusion

This work provides indications of genetic mechanisms presented by Braford cattle with different levels of resistance in response to tick infestation, contributing to the search of candidate genes for tick resistance in bovine.

Background

Cattle are the preferential hosts of Rhipicephalus microplus, a hard tick that attaches to the host skin and feeds for three weeks. Tick attachment and feeding depend upon numerous saliva components that inhibit host hemostatic responses to the parasite bites [1], a process that is the result of millions of years of evolution [2]. Rhipicephalus microplus is the most important ectoparasite of livestock, especially in tropical and subtropical areas [3], causing severe illness in cattle [4], with annual global costs of around US$ 22–30 billion [5].

Acaricides are currently the most common tick control method. However, significant levels of resistance to the different acaricide classes [6, 7] along with potential contamination of milk, beef, and the environment no longer support their use. Vaccination is an alternative for tick control, and several efforts have been conducted to increase its effectiveness [8,9,10]. Genetic resistance can be a permanent solution to tick control [11]. Bovine resistance to R. microplus infestation is a heritable phenotype, and heritability values around 0.34 were observed across different populations [12, 13]. Tick resistance has been studied in several cattle breeds [14,15,16,17,18], and Bos taurus taurus breeds are more susceptible to tick infestation compared to Bos taurus indicus breeds [11]. Braford, a composite breed of 3/8 Zebu (B. t. indicus) and 5/8 Hereford (B. t. taurus), presents considerable variation in tick resistance.

Previous attempts to understand the genetic mechanisms underlying resistance explored host immune responsiveness. However, differences in experimental design hamper comparison of results [19, 20]. Tick response was primarily compared between zebuine and taurine cattle breeds [14, 21,22,23,24,25,26,27], but differences have been reported between low and high resistance levels within the same breed [28,29,30,31]. Overall, these studies show an essential role of the structural protein-coding genes and cellular immunity through innate and acquired mechanisms, including cytokines, chemokines, T cells, B cells, mast cells, and granulocytes. Although some cellular characterization has been done, bovine gene expression has been mainly assessed by RT-PCR and microarray assays, limiting differential gene expression analyses in terms of the number of genes investigated. Studies from Australia [32] and Brazil [33,34,35,36] have identified Quantitative Trait Loci (QTL) underlying host variation. While the Australian study tested a candidate gene (integrin alpha 11), Brazilian studies were either based on microsatellite technology, resulting in wide confidence intervals for QTL locations, or on TagSNPs. Though these studies brought to light some aspects of genetic influence in tick resistance, they did not investigate, in a more comprehensive manner, which genes and pathways are involved in resistance. Jonsson et al., 2014 [14], reviewing genetic marker research, suggested that cell-mediated immunity, hypersensitivity, local inflammation and structural skin components contribute to host resistance. Porto-Neto et al. [37] validated a positional candidate gene (receptor-interacting serine-threonine kinase 2, RIPK2) for tick burden using a knock-out mice model. Genomic approaches to explore genetic variation affecting tick host resistance have been reviewed [38].

In this report we used high throughput RNA sequencing technology to compare gene expression in tick resistant and susceptible Braford cattle. Functional enrichment and network analyses were employed to uncover genetic mechanisms of host resistance to tick infestation. The mechanisms identified could contribute to the understanding of host immunity against ticks.

Results

Differential gene expression analysis

A comprehensive study of skin transcriptomic profile of genetically divergent hosts regarding anti-tick resistance was performed using RNA-Seq technology. Data from 20 resistant (R) and 19 susceptible (S) animals, previously selected from an original population of 974 animals, were collected prior to first artificial tick infestation (pre) and 24 h after the fourth infestation (post). This strategy was employed to address both innate and acquired immunity before the challenge and those elicited by it. The number of reads after filtering was around 10 M per sample. No differences among groups were observed for the mapping statistics with the reads presenting around 70% concomitant pair alignment rate at the gene level. Read information and mapping statistics are presented in Additional file 1: Table S1. Fig. 1 shows Venn diagrams illustrating the distribution of differentially expressed genes (DEGs) from intergroup: Rpre vs. Spre (1) and Rpost vs. Spost (2) and intragroup: Rpost vs. Rpre (3) and Spost vs. Spre (4) comparisons, including unknown genes (Fig. 1a), only the annotated genes (Fig. 1b), and only the annotated TopDEGs (|log2| FC > 1) (Fig. 1c).

Fig. 1
figure 1

Functional enrichment of differentially expressed genes (DEGs) in Braford skin according to host phenotype and tick infestation. Venn diagrams show the distribution of DEGs from inter- (Rpre vs. Spre (1) and Rpost vs. Spost (2)) and intragroup (Rpost vs. Rpre (3) and Spost vs. Spre (4)) comparisons: (a) all DEGs, (b) only the annotated ones, or (c) only annotated TopDEGs. Functional annotation based on protein classes of DEGs from (d) inter- and (e) intragroup comparisons are shown, represented by a symbol following the Metacore® reference guide (https://portal.genego.com/legends/MetaCoreQuickReferenceGuide.pdf). Actual: number of network objects from the dataset(s) for a given protein class; n: number of network objects in the dataset(s); R: number of network objects of a given protein class in the complete database or background list; N: total number of network objects in the complete database or background list; Expected: mean value for hypergeometric distribution (n*R/N); p-value: probability to have the given value of Actual or higher (or lower for negative z-score); z-score: ((Actual-Expected)/sqrt(variance)); Ratio: connectivity ratio (Actual/Expected); % in Dataset: fraction of network objects with a selected function in the dataset

Full size image

At a false discovery rate (FDR) < 0.05, Rpre vs. Spre (1) showed 56 DEGs, with 35 up- and 21 down-regulated (Table 1 and Fig. 1a), and 43 annotated genes (Table 1 and Fig. 1b). Rpost vs. Spost (2) showed 63 DEGs, with 37 up- and 26 down-regulated (Table 1 and Fig. 1a), and 48 annotated genes (Table 1 and Fig. 1b). Among annotated DEGs, 41 genes are TopDEGs in Rpre vs. Spre and 47 in Rpost vs. Spost (Fig. 1c), with overall fold change ranging from 7.59 to − 9.59 (Table 1). The intergroup comparisons showed 38 redundant TopDEGs (Fig. 1c, intersection between 1 and 2). The correlation (r2) among log2 FC values from these redundant TopDEGs was 0.99.

Table 1 Intergroup comparisons of differentially expressed genes in the skin of tick-resistant and -susceptible Braford cattle
Full size table

Regarding cellular metabolism, many TopDEGs were implicated in cell activation against injury. We observed functions such as cellular signaling, ion-dependent vesicular trafficking and transport, free radical depuration, and detoxification of products of oxidative stress with chaperone for superoxide dismutase (CCS) and glutathione S-transferase mu 1 (GSTM1) genes; cytoskeleton organization with kinesin family member 23 (KIF23); ribosomal processing with U3 small nucleolar ribonucleoprotein (IMP3); cellular growth and differentiation with wingless-type MMTV integration site family member 7A (WNT7A) and neudesin neurotrophic factor (NENF); as well as regulation of gene expression with FBJ murine osteosarcoma viral oncogene homolog (FOS) (Table 1). Genes are listed alphabetically with their Ensembl identification and description.

Concerning response against ticks and immunity, TopDEGs found in both comparisons such as epidermal arachidonate lipoxygenase (ALOX12E), acyl-CoA wax alcohol acyltransferase 1 (AWAT1), serum amyloid A 3 (SAA3) and tachykinin receptor 2 (TACR2), have products acting on inflammation. For instance, AWAT1 showed the greatest modulation among all DEGs, being down-regulated in Rpre compared to Spre and in Rpost compared to Spost (respectively log2 FC − 9.59 and − 9.35), whereas SSA3 (log2 FC 5.05 and 5.55), TACR2 (log2 FC 2.34 and 2.39) and WNT7A (log2 FC 7.10 and 7.59) genes were strongly up-regulated in the same comparisons. Unique TopDEGs as dihydropyrimidinase-like 4 (DPYSL4) and urotensin 2 receptor (UTS2R) were up-regulated only in Rpre with high log2 FC (6.14 and 3.29, respectively), whereas 5,10-methenyltetrahydrofolate synthetase (MTHFS) was down-regulated (log2 FC − 1.17). TopDEGs such as galactose-3-O-sulfotransferase 1 (GAL3ST1), retinol binding protein 1 (RBP1), tetmethylcytosinedioxygenase 1 (TET1), troponin T type 3 (TNNT3), uroplakin 3A (UPK3A) and four incompletely annotated myeloid-associated differentiation marker-like (MYADML) genes were modulated only in Rpost vs. Spost (Table 1).

In the intragroup comparisons, R hosts presented 2523 (FDR < 0.05) DEGs after tick infestation, of these 1807 were up- and 716 were down-regulated (Fig. 1a and Additional file 2: Table S2), with 95.52% (2410) annotated (Fig. 1b and Additional file 2: Table S2). The log2 FC varied from 5.28 to − 3.77, and about 12.5% (316) were TopDEGs with 200 up- and 116 down-regulated genes. Of these, 286 DEGs were annotated (Fig. 1c and Additional file 2: Table S2), with a considerable number of uncharacterized or non-annotated genes remaining among TopDEGs. Susceptible hosts showed 2120 significant (FDR < 0.05) DEGs after tick infestation, with 1442 up- and 678 down-regulated genes (Fig. 1a and Additional file 2: Table S2), from which 2024 (95.47%) were annotated (Fig. 1b and Additional file 2: Table S2). The log2 FC variation was similar to that observed for R, from 5.23 to − 3.65 with 260 TopDEGs (228 annotated), 146 (130) up- and 114 (98) down-regulated genes (Fig. 1c and Additional file 2: Table S2). TopDEGs were distributed over the whole host genome in both phenotypes (except for chromosome 20 in S hosts). We also found unique and redundant DEGs (and TopDEGs) in the intragroup comparisons (Fig. 1). Almost two thousand genes were modulated exclusively in one resistance-group (1177 DEGs, 1121 annotated for R and 775 DEGs, 739 annotated, for S). Nonetheless, more than a thousand DEGs (1337, 1282 annotated) were redundant for both resistance groups after tick infestation with a highly correlated variation on expression (r2 = 0.97). These genes point to mechanisms elicited by tick challenge independent of the host phenotype, such as immune response, coagulation, vascularization, and ion transport.

Enrichment analysis (EA) based on functional ontologies

Functional enrichment based on non-redundant TopDEGs predicted more than 30 gene functions, with several categories activated after infestation (p-values ranging from 5.84E-03 to 1.55E-07 and z-score > 2), and only one inhibited (p-value 3.43E-04 and z-score < − 2) in R hosts (Table 2). Unique TopDEGs were classified according to the functional category and specific function they represent. Most of the activated functions are somehow related to anti-tick response, as cellular activation and migration, inflammation, lipid metabolism, molecule transport and blood vessel formation. Activation and chemotaxis were observed for all immune cell types such as leukocytes, phagocytes, and granulocytes, and mostly represented by cytokines, chemokines, growth factors, and inflammatory genes in R hosts (Table 2).

Table 2 Functional enrichment of non-redundant differentially expressed genes presenting higher variation in the skin of tick-resistant Braford cattle after tick infestation
Full size table

Among S hosts, seven functions were predicted as activated according to the same criteria (Table 3). We also identified anti-tick related functions as cellular movement (except for granulocytes), cell adhesion, leukocyte immune response, and quantity of calcium.

Table 3 Functional enrichment of non-redundant differentially expressed genes presenting higher variation in the skin of tick-susceptible Braford cattle after tick infestation
Full size table

The over-representation of molecular functions evaluated for all annotated DEGs, and the statistical values and symbols which represent the protein classes are presented in Fig. 1d and e. In the intergroup comparisons (Fig. 1d), DEGs from R hosts were enriched (p < 0.05) for two protein classes, enzymes and “other” for proteins that do not belong to any other class listed, both before (Rpre) and after (Rpost) tick infestation. In the intragroup comparisons (Fig. 1e), S hosts showed over-representation of two protein classes after tick infestation (Spost) not represented in R: receptors and transcription factors (TFs).

Next, functional enrichment by ontology was applied to identify over-represented processes and pathways among TopDEGs from intragroup comparisons, looking for potential anti-tick responses differently represented between hosts. Enrichment for TopDEGs included immune and inflammatory processes common to both hosts and others exclusive to R (Additional file 3: Table S3A). Leukocyte chemotaxis, T helper cell 17 (Th17)-signaling and Jak-STAT pathway represent enriched processes for TopDEGs from R and S hosts after tick infestation, whereas proteolysis in connective tissue degradation and extracellular matrix (ECM) remodeling, inflammation by cytokines signaling and platelet-endothelium-leukocyte interactions were enriched only for R.

In the leukocyte chemotaxis process network, C-C motif chemokine ligand 13 (CCL13) gene was up-regulated in both phenotypic classes, although its receptor chemokine C-C motif receptor 3 (CCR3) gene was concomitantly up-regulated only in R hosts. Moreover, only R hosts showed two types of cytoskeletal proteins coding genes, tubulin alpha and actinin alpha, as TopDEGs. The lymphocyte attractant CCL20 gene was strongly down-regulated in R hosts, whereas CCR1 was TopDEG only in S, as well as its ligand, the macrophage inflammatory protein 1 (MIP-1 alpha). Finally, S hosts showed over-expression of the interferon-inducible T-cell alpha chemoattractant (I-TAC) gene after tick infestation.

TopDEGs from both hosts were also committed to the Th17-derived cytokines immune process (Additional file 3: Table S3A) with some differences. R hosts over-expressed matrix metallopeptidase 9 (MMP9), granulocyte-macrophage-colony-stimulating factor (GM-CSF), tissue inhibitor of metalloproteinase 1 (TIMP1), interleukin 6 (IL6) and interleukin 22 (IL22) transcripts that code for proteins involved in response to IL17 as well as in the differentiation of the Th17 cell. Regarding S hosts, we found the over-expression of the cytokine receptor activator of nuclear factor kappa-B ligand (RANKL) transcript and the absence of mRNAs corresponding to critical cytokines coding genes in skin inflammation, such as IL6 and IL22.

The Janus kinase-signal transducers and activators of transcription (Jak-STAT) pathway involves signaling cascades that transmit extracellular signals to target genes in the nucleus. Overall, both R and S hosts showed increased gene expression of external stimuli signaling, mainly represented by cytokines, which leads to the activation of this pathway in response to tick infestation. However, S hosts showed an additional inhibition of ciliary neurotrophic factor (CNTF) and suppressor of cytokine signaling 2 (SOCS2) genes.

Connective tissue degradation and skin ECM remodeling were enriched only for R hosts since genes such as MMP9, TIMP1, kallikrein 1 (KLK1) and KLK6 presented low FC in S. Additionally, target proteins coding genes were down-regulated in R hosts, such as collagen I, collagen III, lumican, and osteonectin. Overall, R hosts showed more TopDEGs (Additional file 3: Table S3A) related to these two processes after tick infestation when compared to S.

Finally, the TopDEGs from R hosts enriched inflammation by IL6-signaling with up-regulation of genes participating in this signaling pathway beyond IL6 itself, such as acute-phase protein coding genes such as SAA3, haptoglobin (HP), and lipopolysaccharide-binding protein (LBP), and the acquired-immunity cytokine IL22 (Additional file 3: Table S3A).

The enrichment based on pathways pointed to ECM remodeling and some cytokines for TopDEGs from both hosts, but indicated the putative participation of T cells only in R (Additional file 3: Table S3B). In fact, through pathway analysis, DEGs from R hosts supported the participation of T CD4+ cell subtypes, including IL22-signaling pathway (Fig. 2). Genes, such as CD4, IL22, and c-Fos, were modulated only in R. Others, as tyrosine kinase 2 (Tyk2), ERK1/2 and SOCS3, were modulated in both hosts, while only S hosts showed modulation of STAT3. All DEGs in this pathway were up-regulated, including the effectors (as c-Fos) and the regulators (as SOCS3).

Fig. 2
figure 2

Interleukin 22 (IL22)-signaling pathways in different cell types. IL22, produced mainly by immune CD4+ T cells, acts upon many different tissue cells, with a key role in the skin, digestive, and respiratory tracts. IL22 acting via its heterodimeric receptor induces activation of many kinases (via MAPs and Jak-STAT pathway) and transcription factors compromised with the expression of pro-inflammatory cytokines. DEGs from inter- (Rpre vs. Spre (1) and Rpost vs. Spost (2)) and intragroup (Rpost vs. Rpre (3) and Spost vs. Spre (4)) comparisons and their respective modulation are shown aside the network objects. Values of expression are represented by colored bars in thermometer-like icons. The red bar represents up-regulation. Network object symbols follow the Metacore® reference guide according to the protein class (https://portal.genego.com/legends/MetaCoreQuickReferenceGuide.pdf), and their names may occasionally differ between tables and the figure due to the MetaBase®

Full size image

Network Analysis. The WNT-signaling pathway, including WNT/β-catenin dependent, was the most relevant network predicted for R hosts (Fig. 3). This pathway is strongly supported by TopDEGs, especially the WNT7A gene (Fig. 3a), and also included DEGs (Additional file 2: Table S2). Among them, the TF-coding gene lymphoid enhancer-binding factor 1 (LEF1) and the anti-apoptotic and inflammatory gene B-cell lymphoma 9 protein (bcl9) were up-regulated in R, whereas plasminogen activator, urokinase (PLAU), PLAU receptor (PLAUR), low-density lipoprotein receptor-related protein 12 (LRP12) and members of the fibroblast growth factor (FGF) family were up-regulated in both hosts. On the other hand, the receptor coding gene LRP2 was down-regulated only in S. There was no evidence of the involvement of WNT-Ca pathway (or non-canonical WNT pathway) in resistance. The complete list of predicted networks for R hosts is shown in Additional file 4: Table S4, and includes skin biology, cell surface receptor signaling, response to an endogenous and external stimulus, inflammatory and defense pathways.

Fig. 3
figure 3

Wingless (WNT)-signaling pathway as a central way to Braford anti-tick resistant (R) response. (a) The canonical WNT-family pathway was predicted for the TopDEGs from R hosts after tick infestation, highlighted in cyan. Red circles represent over-expression: higher the value, stronger the intensity of red. (b) Genes prospected from genome-wide association study also indicate the WNT-signaling pathway as the most relevant in resistance. The red circle represents the prospected gene. Network object symbols follow the Metacore® reference guide according to the protein class (https://portal.genego.com/legends/MetaCoreQuickReferenceGuide.pdf), and their names may occasionally differ between tables and the figure due to the MetaBase®

Full size image

Network analysis was next applied to a subset of genes, based on LD (linkage disequilibrium) information from this population [39], within 200 kb genomic regions from each side of TagSNPs explaining up to 20% of the resistance variation in a GWAS [36], similar to the approach described by Mota et al. (2017) [40]. The canonical WNT network (β-catenin dependent) was predicted as the most relevant pathway for those genes (Fig. 3b and Additional file 5: Table S5), corroborating the involvement of this pathway in resistance and supporting our finding of WNT7A as TopDEG in the intergroup comparisons (Table 1). The complete list of predicted networks for genes prospected from the GWAS is shown in Additional file 5: Table S5.

On the other hand, the most relevant network for S hosts involved MMP1, oncostatin M (OSM), stromelysin-1, CCL2, and IL3 genes (Fig. 4 and Additional file 6: Table S6) and was related mainly to the regulation of immune responses. Other networks also pointed to immune-related processes, but with low connectivity scores.

Fig. 4
figure 4

Regulation of immune response as the most relevant pathway in the infested susceptible (S) skin. Regulation of immune response was predicted for the TopDEGs from S hosts after tick infestation, highlighted in cyan. Red circles represent over-expression: higher the value, stronger the intensity of red. Network object symbols follow the Metacore® reference guide according to the protein class (https://portal.genego.com/legends/MetaCoreQuickReferenceGuide.pdf), and their names may occasionally differ between tables and the figure due to the MetaBase®

Full size image

Finally, the prediction of potential TFs underlying the regulation of TopDEGs from R hosts pointed to Ap1 family members, such as c-Fos, JunB, Fra1, FosB, and Fra2 (Additional file 7: Table S7). Fra1 and Fra2 were also TopDEGs. As the TF lists were manually curated according to the recently published bovine TFs compendium [41], only TFs found in cattle were shown. For S hosts, a nuclear factor kappa B (NF-kB) family member, the NF-kB subunit p65 (RelA), was the most relevant TF predicted. Although some Fos family members also appeared in the list, such as c-Fos and c-Jun, the NF-kB family was represented by three other members. Of which, the most relevant TFs are as follows: NF-kB subunit p50 (NFKB1), cRel and RelB. The TopDEG TFs from S hosts were the nuclear factor, erythroid 2 (NF-E2), another NF-kB family member, CCAAT/enhancer-binding protein (C/EBPE) and paired box 8 (PAX8), which were classified below the 60th position according to their p-value and z-score (Additional file 8: Table S8).

Discussion

Differential gene expression analysis

The mechanisms underlying resistance and susceptibility against ticks of phenotypically divergent Braford hosts were addressed using RNA-Seq in cattle skin. With this approach, the extension of genetic modulation induced by ticks could be better examined in terms of differential gene expression and their collective participation in resistance (and susceptibility) at system biology levels. The gene expression profiles found before challenge (pre) could be partially explained by the immunity previously acquired by hosts, since they were not naive. Additionally, 24 h after the fourth tick infestation (post), we expected to see both innate and acquired immune responses, as already demonstrated for re-infested animals in skin-level studies of gene expression [14, 22, 23, 28].

Through the intergroup comparisons, we investigated the gene expression profiles of R and S hosts, looking for candidate genes that would explain phenotypic differences, according to the infestation challenge: Rpre vs. Spre (1) and Rpost vs. Spost (2). Both differential expression analyses revealed highly redundant DEG lists, suggesting that mechanisms of resistance in Braford cattle were present before challenge. Natural infestations before the first skin sampling, the time of skin sampling after infestation (24 h) or even the absence of the tick attachment site at the biopsy could contribute to these findings. Almost all annotated DEGs were classified as TopDEGs with biological functions correlated to cattle resistance against ticks. Overall, they were involved in hemostatic and immune mechanisms, many of them already associated with tick infestation, such as coagulation, iron metabolism, and inflammation [42,43,44]. It is conceivable that the modulation of these functions could help hosts achieve a better initial cellular response against tick antigens when tick feeding starts.

AWAT1, which had the highest variation in the intergroup comparisons, encodes an enzyme that acts in lipid metabolism and sebum composition in the skin [45]. Lipid metabolism has a role in inflammation control [46] and semiochemistry [47], and probably also in anti-tick resistance [24]. The tachykinin family (represented by TACR2 gene) was widely associated with inflammatory response, playing a role in wound healing and immune cell differentiation [48, 49]. Seric amyloid (SAA) family members were associated with anti-tick resistant response [50], and SAA3 expression was up-regulated in R hosts in addition to the intragroup comparison (Additional file 2: Table S2). Amine oxidase copper containing 3 (SAO) gene product is committed to leukocyte trafficking [51, 52], a function already associated with anti-tick responses [53], as well as in vascularization and tissue organization. Genes, such as bolA family member (BOLA) and tetratricopeptidase repeat domain 36 (TTC36), are related to immune mechanisms, including antigen processing and presentation. Moreover, genes involved in skin organization as collagen type XI alpha 1 (COL11A1), laminin alpha 1 (LAMA1) and agouti signaling protein (ASIP) may also contribute physically to tick bite healing, unsuccessful attachment, and establishment of the feeding site. Finally, the ectonucleotide pyrophosphatase/phosphodiesterase 3 (ENPP3) gene encodes an enzyme committed to the hydrolysis of extracellular nucleotides, also known as CD203c. This enzyme is better described as an activation marker of basophils, a cell type associated with acquired anti-tick resistance [25, 54, 55].

Genes modulated only before infestation (Rpre) are probably involved in primary cellular functions that may prepare R hosts to promptly react to infestation. DPYSL4 gene is described as a target of p53 in the apoptotic response to DNA damage [56]; once up-regulated, it could contribute to the apoptosis in the tick bite site. MTHFS is associated with folate turnover rate and depletion [57], and folate is essential for healthy skin [58]. Finally, over-expression of UTS2R, even in the absence of active infestation, suggests better control of vascular dynamics and osmoregulation in R hosts.

On the other hand, genes modulated only in Rpost may indicate putative mechanisms of protective response induced by the presence of ticks. Although the proteins encoded by genes such as UPK3A, TET1 and RBP1 do not play a direct role in immune cell functions, they indirectly help immune response. The function of UPK3A gene product has been linked to epithelial differentiation [59]. TET1, inhibited in Rpost, encodes a transcription repressor of a subset of genes, including interleukin 1 beta (IL1b), IL8, IL23, intercellular adhesion molecule 1 (ICAM1) and chemokine C-X-C motif ligand 1 (CXCL1) [60]. Most of these genes, modulated by TET1 expression, were up-regulated in Rpost (Additional file 2: Table S2). RBP1, also up-regulated in R hosts, encodes the carrier protein of retinol (vitamin A alcohol) transport, making it accessible in peripheral tissues. Besides the role in epithelial tissues, the hormone-like properties of vitamin A binding to nuclear hormone receptors retinoic acid receptors (RAR) and peroxisome proliferator-activated receptor (PPAR) affect immunity [61, 62]. This gene was up-regulated in Rpost compared to Spost and down-regulated in Spost compared to Spre, suggesting its involvement in tick resistance. Although the MYADML family is still incompletely annotated, some products have been described as modulators of cell spread and migration to the skin [63]. Three members of this family were modulated in Rpost with contrasting modulation (log2 FC -5.88, − 2.4 and 6.33; Table 1), which could indicate antagonistic effects regarding response to tick infestation.

In the intragroup comparisons, Rpost vs. Rpre (3) and Spost vs. Spre (4), we were looking for genes elicited by tick challenge, according to the resistance group studied. They revealed extensive gene lists, with thousands of DEGs broadly distributed over the genome, showing both redundancy and high correlation between hosts DEGs as well as non-redundant DEGs for R or S hosts (Additional file 2: Table S2). Critical functions against infestation, such as immune responses, coagulation, vascularization, and ion transport, were shared between hosts. This indicates that tick infestation caused a comprehensive response in Braford skin, regardless of being an R or an S host. On the other hand, non-redundant DEGs with potential anti-tick functions may represent mechanisms of response that contribute to differentiating hosts regarding resistance. This hypothesis was supported by the activation of critical functions in R hosts with massive participation of all immune cells (Table 2), which had the potential to protect them. Among S, cellular activation was not predicted for granulocytes (Table 3), an important cell type for tick resistance [19, 23, 54]. Overall, the data suggest the relevance of cytokines in R skin reactions. The inhibition of grooming predicted for R hosts is probably due to the time of skin sampling (24 h) when larvae had already been rejected [64]. Although very little is known about this, grooming behavior is considered to have an impact on tick load [15] and probably to ectoparasites other than ticks as well.

More than a hundred immune-related genes encoding cytokines, chemokines, CD markers, acute phase proteins, complement proteins, integrins, and TFs were found (Additional file 2: Table S2). Our results revealed novel DEGs associated with resistance/susceptibility, DEGs that support previous results and predictions, and DEGs that diverge with authors comparing gene expression between taurine and zebuine breeds [14, 21,22,23,24,25], within other cattle breeds [28,29,30,31] or murine models [37, 65]. According to Piper and colleagues [21, 23], CXCL2 and CCL2 were modulated only in S hosts; however, we found them also modulated in R. CCR1 and IL2R were described as modulated in R and CD14 in S [26], in contrast with our results where CCR1 and CD14 were modulated in both hosts, and IL2RG up-regulated only in S. IL8 was up-regulated in both R and S, but was previously described as down-regulated in R [33] and up-regulated in S [23]. IL13RA1 was described as up-regulated in S and down-regulated in R [28], whereas we found IL13RA1 to be up-regulated in S and IL13RA2 to be up-regulated in R. Regarding T cell participation in resistance [14, 22, 23, 25, 27, 30, 31], CD4 was up-regulated only in R hosts, in agreement with previous results. These inconsistencies may be the result of a myriad of factors, including sensibility of different techniques for gene expression analysis, experimental design, epistatic interactions with the genetic background and environmental effects, and indicate the need for extensive research on specific tick-host relationships.

Enrichment analysis (EA) based on functional ontologies

In the intergroup comparisons, similar protein classes were enriched among TopDEGs, although the number of TopDEGs was insufficient for network analysis. In the intragroup comparisons, S hosts showed more enriched protein classes than R. Functional enrichment was applied to TopDEGs (Additional file 3: Table S3) as this analysis is a powerful tool to determine the relevance of a gene in the context of a pathway. Some processes, which can be crucial to managing the interference of tick in the skin immune biology were over-represented in both hosts, including inflammation, chemotaxis, and immune response. Connective tissue degradation, matrix remodeling and immune process networks found only in R hosts likely indicate mechanisms they use to respond to ticks which are not elicited in S.

Even though both hosts had enriched leukocyte chemotaxis, differences between R and S in the modulation of this network may contribute to susceptibility (Additional file 3: Table S3A). For example, CCR3 was a TopDEG up-regulated only in R, although the corresponding ligand CCL13 was up-regulated in both phenotypes, thus suggesting differences in CCL13 effectiveness. CCL13 is involved in chemotaxis of eosinophils, basophils, and Th1 cells, all of which were previously associated with resistance [19, 23, 54]. The modulation of cytoskeletal proteins may also contribute to more efficient chemotaxis. On the other hand, the down-regulation of CCL20 in R hosts could indicate inhibition of the recruitment of lymphocyte subpopulation to the resistant skin [66]. The over-expression of MIP-1 alpha and I-TAC in S hosts are other examples: the former encodes a chemokine involved in inflammation, especially in attracting neutrophils and macrophages, which are known to have their functions impaired by tick saliva [43, 67, 68]; I-TAC, a chemoattractant of activated lymphocytes produced mainly by basal keratinocytes [69], did not find its receptor coding gene C-X-C motif chemokine receptor 3 (CXCR3) modulated, which was expected in activated T cells responding to infestation. Although tick infestation also elicits a cellular response in S hosts, it may attract less efficient cells to their skin.

Additionally, Th17-derived cytokines network was also over-represented among TopDEGs from both hosts (Additional file 3: Table S3A). Th17 cells belong to a distinct subset of CD4+ Th lymphocytes characterized as preferential producers of IL17A, IL17F, IL22, and, in humans, IL26. The receptors for IL17 and IL22 are broadly expressed in various epithelial tissues. The effector cytokines released from Th17 cells affect different cellular populations in inflammatory sites, mediating host defense through the activation of several signaling pathways. The participation of the IL17 family in cattle anti-tick response has already been proposed [25].

Finally, Jak-STAT pathway was also enriched for both hosts (Additional file 3: Table S3A). In this pathway, membrane-associated tyrosine kinases lead to the phosphorylation of signaling proteins and transcription factors such as STATs, regulating the expression of genes involved in the inflammatory response, cell survival, and cell cycle in response to numerous cytokines and growth factors. S hosts showed down-regulation of CNTF whose product, in the nervous system, appears to protect tissue during inflammation, and SOCS2, a Jak-STAT regulator coding gene whose product interacts with major molecules of signaling complexes to block further signal transduction.

Connective tissue degradation and ECM-remodeling networks were over-represented only among TopDEGs from R. These are two particularly important processes which require the over-expression and activation of proteinases able to destroy the connective tissue, allowing non-resident inflammatory cells arrival to the skin [23]. Most proteolytic enzymes are MMPs, and proteins including disintegrin and metallopeptidase domain (ADAMs), KLKs, elastases, trypsin, and others also participate. The MMPs comprise a family of enzymes that can collectively degrade all components of the ECM. They are the main proteases involved in ECM remodeling and are negatively regulated by TIMPs. The ability to modulate the skin immune-biology to impair tick success may represent a great difference between R and S hosts.

Inflammation via IL6-signaling was also over-represented only in TopDEGs from R hosts, in contrast with previous reports in which both R and S hosts show up-regulation of IL6 [25]; it is also reported in mice-infested skin, a model considered susceptible to tick [65]. This cytokine is secreted during inflammation by macrophages, endothelial cells, and fibroblasts. Its primary effect is to induce the secretion of acute-phase proteins, whose role as mediators of inflammatory responses has already been associated with resistance against ticks [24, 50]. IL6 also participates in the differentiation of some cell types, particularly in inducing B lymphocytes to differentiate into plasma cells; it can also lead to anti-apoptotic effects [70].

Given this, the up-regulation of IL6 together with many acute-phase protein coding genes and IL22 may create a favorable scenario for the maintenance of healthy skin barrier, since IL22 plays a critical role in skin immunity, inflammation, and repair [71]. Th22 participation in tick resistance was supported by the pathway analysis using TopDEGs from R hosts. Previous works have demonstrated the massive infiltration of CD4+ cells in the skin of infested animals as well as an increased number of these cells in the bloodstream and lymph nodes [14, 22, 26, 27, 30, 31]. IL22 is a member of the IL10 family of cytokines [72] that acts via a heterodimeric receptor complex (IL22 receptor) consisting of the subunits IL22RA1 and IL10RB (Fig. 2). IL22 binding to IL22 receptor activates the Jak-STAT pathway, inducing the phosphorylation of JAK1, Tyk2, and of the transcription factors STAT1, STAT3, and STAT5. IL22 also activates three major mitogen-activated protein (MAP) kinase pathways, c-Jun NH2-terminal kinase (JNK), p38 MAPK and extracellular signal-regulated kinase (ERK1/2), which are possibly required for activation of c-Jun/c-Fos to elicit the expression of pro-inflammatory cytokine genes [73, 74]. The primary sources of IL22 are immune cells, whereas its receptor is found in non-hematopoietic cells, especially at outer body barriers such as skin and the digestive and respiratory tracts [71, 75, 76]. Many genes in the Th22 pathway were DEGs in R hosts after infestation, including IL22 itself and other genes involved in response to this cytokine in tissue cells, as the tyrosine kinases and c-Fos (Fig. 2). Additional file 2: Table S2 illustrates substantial evidence for the differential activation of IL22-mediated immune response in R hosts, with the up-regulation of IL20, IL22, HP, JUN, FOSL1, IL9R, MMP13, CCR3, colony-stimulating factor 2 (CSF2) and IL6 in R but not in S, suggesting a putative role of this immune mechanism in resistance against ticks in Braford cattle. S hosts also showed some DEGs in this pathway, though neither IL22 was up-regulated, nor IL6, whose product acts in the Th22 cell differentiation [71]. Robbertse and colleagues (2018) [19] suggested the modulation of cytokines, such as IL6, IL13, IL5, chemokines, and of receptors as CCL2 and CCR1, acute phase proteins, cell adhesion molecules, CD4 and CD14 molecules, with the resulting attraction of T and B lymphocytes and granulocytes to the skin, to be associated with resistance. Our results partially agree with the predictions of these authors, since the modulation of some DEGs was redundant between R and S Braford hosts.

Network analysis

This analysis was carried out to identify the central mechanisms and to better explain contrasting host phenotypes based on their differential gene expression. It showed partially overlapping size-limited sub-networks, expanded from the TopDEGs lists, giving preference to objects with more connectivity considering that highly connected genes tend to be involved in similar biological functions.

WNT-signaling, the most relevant pathway predicted for R hosts (Additional file 4: Table S4 and Fig. 3a), is an evolutionarily conserved pathway that regulates the key process of organogenesis during embryonic development, cell fate determination, cell migration and polarity [77]. Recently, this pathway has also been associated with functions in immune response [78, 79] and skin biology as well [80,81,82]. Its effects on inflammation can be both pro- and anti-inflammatory, depending on the stimulus, cell type, activation context, and its crosstalk with other signaling pathways [78, 83]. It becomes inflammatory in infections with pathogenic bacteria and in inflammatory bowel diseases [78, 84]. However, it has also been associated with the induction of tolerogenic dendritic cells (DCs), which participate in immune response regulation mainly through the inhibition of Th1/Th17 cytokines [79]. It also participates in skin wound healing, through its ability to induce self-renewal and proliferation of skin stem cells [80, 82]. Finally, this pathway plays a critical role in many physiological processes of the skin from the earliest stages of development to postnatal control of hair cycling, inter-follicular epidermis maintenance as well as hair and skin pigmentation [81]. Furthermore, the deregulation of this pathway has already been associated with skin diseases [81] and host-tick interface [65]. Hair coat morphology and pigmentation were also associated with tick counts [85, 86]. Additionally, the WNT-signaling pathway was predicted in our analysis of genes mined from TagSNPs associated with resistance to ticks in the same breed [36] (Additional file 5: Table S5 and Fig. 3b). The importance of wound healing and structural proteins for tick resistance has been suggested [21, 37, 65]. Our results suggest the participation of WNT pathway in anti-tick resistance through an inflammatory scenario comprising the regulation of cytokines, receptors, matrix proteinases, and transcription factors and, at the same time, through its role in skin healing and remodeling, which hinders tick success.

On the other hand, the most relevant pathway for TopDEGs from S hosts includes MMP1, OSM, stromelysin-1, CCL2, and IL3 (Additional file 6: Table S6 and Fig. 4), which act on the regulation of immune response. This pathway included molecules related to tissue remodeling, such as OSM, MMPs, and TIMPs, potentially activating Jak-STAT and MAPK signaling cascades [87]. Most of the modulated genes in this pathway were DEGs (not TopDEGs), although Piper and colleagues [23] suggested ticks can more easily modulate the skin of S hosts.

Considering the complexity of host anti-tick resistance, which involves multiple mechanisms of activation and regulation of gene expression to drive the cellular and immune responses, we investigated the role of TFs as key molecules to orchestrate this modulation. Despite the general lack of knowledge concerning the regulation of gene expression in cattle, the TF lists were filtered according to the recently published manually-curated bovine TF compendium [41] to show only bovine TF candidates. For R hosts, Fos family members (c-Fos, c-Jun, Fra1, and FosB) were predicted as the most relevant TFs, with the latest two being TopDEGs (Additional file 7: Table S7), reinforcing the complexity of the regulation of gene expression in R host response against ticks. For S hosts, a nuclear factor kappa B (NF-kB) family member was the main TF based on its g-score (Additional file 8: Table S8), although Fos family members were also identified with lower relevance. Even though we extensively analyzed the data integrating different approaches, it was impossible to cover all aspects of the complexity of the host anti-tick response in this study. Experiments designed to compare different skin biopsies, using co-expression analysis and proteomic approaches could certainly drive the next steps of the investigation of this complex trait.

Conclusions

In the present study, we investigated the response to tick infestation using RNA-Seq data and identified putative genetic mechanisms differentiating R and S Braford hosts. DEGs included genes already associated with resistance using other strategies but were far from being restricted to them. Beyond immune-related genes as cytokines, chemokines, integrins, immunoglobulin superfamily members, acute phase proteins, we found TFs, skin structural, and wound healing related-genes. Collectively viewed, our results provide substantial evidence that R hosts presented an inflammatory response to infestation based on cytokines and the WNT-signaling pathway, potentially impairing tick attachment and feeding success. S hosts seemed unable to assemble an effective anti-tick immunity, even though the infestation also elicited their responses. Our findings shed light on genetic mechanisms and candidate genes underlying tick resistance and susceptibility in Braford cattle.

Methods

Animals and skin sampling

Braford heifers, selected from a 974-cohort group belonging to producers affiliated with the Delta G Connection Breeding Program, were segregated into R (n = 20) and S (n = 19) hosts according to their Estimated Breeding Values (EBVs) and tick counts [88]. Statistic parameter estimates used in animal selection are presented in Additional file 9: Table S9. In addition to having the lower EBVs, the R hosts were always at least one standard deviation (SD) below the cohort average tick count and, when taken the average SD of subsequent counts, they were among the lower 10% of the SD distribution. The same criteria were used to classify as S host, i.e., having the higher EBVs, being at least one SD above the cohort average, and among the 10% lower SD distribution of subsequent counts. These criteria guaranteed R and S hosts had consistently lower and higher tick counts, respectively, within the cohort group. Heifers were transferred from farms to EMBRAPA Pecuária Sul experimental station and maintained for three months in a tick-free pasture receiving an acaricide treatment every two weeks to become non-infested. For the RNA profiling experiment, each animal had four artificial infestations at 14-day intervals, with approximately 20,000 R. microplus larvae (obtained from 1 g of eggs) released on the animal’s dorsal line from the neck to the tail insertion. Skin biopsies were collected from each animal using an 8 mm punch before the first infestation (pre) and 24 h after the fourth infestation (post) without including the tick bite site. The neck area was chosen for sample collection since it was close to the larvae release site and is a frequent place for tick attachment. After the experiment, the sampled heifers were returned to their owner’s farms and reintroduced to their production systems.

Tick counts

Phenotypic data were collected in four successive infestations over five consecutive days (D19 to D23) post-infestation, considering engorged females bigger than 4 mm which were attached on the left side of the host, as described in [88]. The average number of ticks counted in R hosts used in the present study was significantly (P ≤ 0.001) lower (5.48 ticks per body left side) than in S (80.55 per body left side) (Additional file 9: Table S9).

RNA isolation and sequencing

Skin biopsies were mechanically crushed with 7 mm metallic beads in the TissueLyser equipment (Qiagen), and total RNA harvested using the RNeasy Fibrous Tissue kit (Qiagen). Sequencing libraries were prepared individually using TruSeq RNA Sample Preparation kit v2 guide (Illumina, San Diego, CA) with standard protocols, clustered using TruSeq PE Cluster Kit v3-cBot-HS kit (Illumina, San Diego, CA) and sequenced in a paired-end mode (2 × 100 bp) on a HiSeq1500 equipment (Illumina, San Diego, CA). Reads passing the CASAVA software (Illumina, San Diego, CA) quality controls were used for data analysis.

RNA-Seq data processing and analysis

Read quality trimming was performed by using SeqyClean v. 1.3.12 [89] to select reads with average Phred score quality (Q) ≥ 30. Quality control was visualized using FASTQC software v. 0.11.1 (http://www.bioinformatics.babraham.ac.uk/pojects/fastqc). The reads were aligned to the bovine reference genome (UMD 3.1, v 75) using TopHat2 [90] following default parameters. Mapped reads were normalized according to the sequencing depth of their group and counted at the gene level using summarizeOverlaps package (May 2014) from R (v. 3.1.1) through “IntersectionNotEmpty” mode.

Differential gene expression analysis

Differential gene expression analysis was carried out with the edgeR Bioconductor package (v. 3.7.17) [91]. Dispersions were estimated using a generalized linear model (GLM) considering animal, phenotype (R or S) and time (pre or post infestation). Normalized reads were compared for intergroups: Rpre vs. Spre (1) and Rpost vs. Spost (2); and intragroups Rpre vs. Rpost (3) and Spre vs. Spost (4). Differential gene expression was tested using the likelihood ratio test corrected for multiple testing by Benjamini-Hochberg method. Genes with significant p-value after correcting for multiple testing (FDR < 0.05) were classified as DEGs. TopDEGs were established with |log2| FC > 1 and FDR < 0.05. DEGs were annotated based on their Ensembl identifiers using Biomart (v. 2.22.0).

Enrichment and network analysis

Functional enrichment was performed using the Ingenuity Pathway Analysis (IPA - Qiagen) and MetaCoreTM (Thomson-Reuters) softwares. IPA was used to evaluate the activation/inhibition of non-redundant gene functions of R and S hosts (|z-score| > 2), with positive and negative values indicating, respectively, activation and inhibition of the function. Using MetaCoreTM®, the functional enrichment analysis classified genes according to protein categories and identified the over-represented processes (Process Networks) and pathways (Pathway Maps) among TopDEGs. Network analysis was based on proprietary algorithms to build and analyze partially overlapping size-limited sub-networks (100 network objects), expanded from the initial gene list, giving preference to objects with more connectivity. Beyond the p-value and FDR, two other scores were associated. The z-score considers the number of DEGs represented in the network (saturation), and the g-score represents the saturation of the network with expressed genes considering the number of canonical (classical) pathways used to build it. Since the g-score accounts for a more robust parameter to classify networks, it was chosen to represent them. Network analysis was also applied to genes harboring the TagSNPs selected in previous work [16], where 3545 animals from the same population were phenotyped (ticks counted) and genotyped with the Illumina BovineSNP50 BeadChip. Finally, network analysis was applied to predict TFs considering the differential expression of each TF as DEGs in the data sets.

Availability of data and materials

The data supporting the conclusions of this article is available in the European Nucleotide Archive (ENA) repository (EMBL-EBI), under accession number PRJEB33196 [https://www.ebi.ac.uk/ena/browser/home].

Abbreviations

CD:

Cluster of differentiation

DC:

Dendritic cell

DEG:

Differentially expressed gene

EA:

Enrichment analysis

EBV:

Estimated breeding value

ECM:

Extra-cellular matrix

FC:

Fold change

FDR:

False discovery rate

GO:

Gene ontology

GWAS:

Genome-wide association study

ID:

Identification

Jak:

Janus kinase

KLK:

Kallikrein

LD:

Linkage disequilibrium

MMP:

Matrix metallopeptidase

QTL:

Quantitative trait loci

R:

Resistant

RNA-Seq:

RNA- Sequencing

Rpost :

R hosts after infestation

Rpre :

R hosts before infestation

S:

Susceptible

SD:

Standard deviation

SNP:

Single nucleotide polymorphism

Spost :

S hosts after infestation

Spre :

S hosts before infestation

STAT:

Signal Transducers and activators of transcription

TF:

Transcription factor

Th:

T helper cell

TIMP:

Tissue inhibitor of metalloproteinase

TopDEG(s):

DEG(s) with |log2| FC > 1

vs.:

Versus

WNT:

Wingless

References

  1. Francischetti IM, Sa-Nunes A, Mans BJ, Santos IM, Ribeiro JM. The role of saliva in tick feeding. Front Biosci. (Landmark Ed). 2009;14:2051–88.

    Article  CAS  Google Scholar 

  2. Ribeiro JM. Blood-feeding arthropods: live syringes or invertebrate pharmacologists? Infect Agents Dis. 1995;4:143–52.

    CAS  PubMed  Google Scholar 

  3. Jongejan F, Uilenberg G. The global importance of ticks. Parasitology. 2004;129(Suppl):S3–14.

    Article  PubMed  Google Scholar 

  4. de la Fuente J, Estrada-Pena A, Venzal JM, Kocan KM, Sonenshine DE. Overview: Ticks as vectors of pathogens that cause disease in humans and animals. Front Biosci. (Landmark Ed). 2008;13:6938–46.

    Article  Google Scholar 

  5. Lew-Tabor AE, Rodriguez VM. A review of reverse vaccinology approaches for the development of vaccines against ticks and tick borne diseases. Ticks Tick Borne Dis. 2016;7:573–85.

    Article  CAS  PubMed  Google Scholar 

  6. Rodriguez-Vivas RI, Jonsson NN, Bhushan C. Strategies for the control of Rhipicephalus microplus ticks in a world of conventional acaricide and macrocyclic lactone resistance. Parasitol Res. 2018;117:3–29.

    Article  PubMed  Google Scholar 

  7. Abbas RZ, Zaman MA, Colwell DD, Gilleard J, Iqbal Z. Acaricide resistance in cattle ticks and approaches to its management: the state of play. Vet Parasitol. 2014;203:6–20.

    Article  CAS  PubMed  Google Scholar 

  8. Willadsen P, Bird P, Cobon GS, Hungerford J. Commercialisation of a recombinant vaccine against Boophilus microplus. Parasitology. 1995;110(Suppl):S43–50.

    Article  PubMed  Google Scholar 

  9. Willadsen P. Tick control: thoughts on a research agenda. Vet Parasitol. 2006;138:161–8.

    Article  PubMed  Google Scholar 

  10. Parizi LF, Reck JJ, Oldiges DP, Guizzo MG, Seixas A, Logullo C, et al. Multi-antigenic vaccine against the cattle tick Rhipicephalus (Boophilus) microplus: a field evaluation. Vaccine. 2012;30:6912–7.

    Article  CAS  PubMed  Google Scholar 

  11. Frish JE. Towards a permanent solution for controlling cattle ticks. Int J Parasitol. 1999;29:57–71.

    Article  Google Scholar 

  12. Shyma KP, Gupta JP, Singh V. Breeding strategies for tick resistance in tropical cattle: a sustainable approach for tick control. J Parasit Dis. 2015;39:1–6.

    Article  CAS  PubMed  Google Scholar 

  13. Porto Neto LR, Jonsson NN, D’Occhio MJ, Barendse W. Molecular genetic approaches for identifying the basis of variation in resistance to tick infestation in cattle. Vet Parasitol. 2011;180:165–72.

    Article  CAS  PubMed  Google Scholar 

  14. Jonsson NN, Piper EK, Constantinoiu CC. Host resistance in cattle to infestation with the cattle tick Rhipicephalus microplus. Parasite Immunol. 2014;36:551–7.

    Article  CAS  Google Scholar 

  15. Regitano LCA, Prayaga K. Ticks and tick-borne diseases in cattle. In: Bishop SC, Axford RFE, Nicholas FW, Owen JB, editors. Breeding for disease resistance in farm animals. Oxfordshire: CAB International; 2010. p. 295–314.

    Chapter  Google Scholar 

  16. Cardoso FF, Gomes CC, Sollero BP, Oliveira MM, Roso VM, Piccoli ML, et al. Genomic prediction for tick resistance in Braford and Hereford cattle. J Anim Sci. 2015;93:2693–705.

    Article  CAS  PubMed  Google Scholar 

  17. Mapholi NO, Maiwashe A, Matika O, Riggio V, Bishop SC, MacNeil MD, et al. (2016) genome-wide association study of tick resistance in south African Nguni cattle. Ticks Tick Borne Dis. 2016;7:487–97.

    Article  CAS  PubMed  Google Scholar 

  18. Piccoli ML, Brito LF, Braccini J, Cardoso FF, Sargolzaei M, Schenkel FS. Genomic predictions for economically important traits in Brazilian Braford and Hereford beef cattle using true and imputed genotypes. BMC Genet. 2017;18:2–15.

    Article  PubMed  PubMed Central  Google Scholar 

  19. Robbertse L, Richards SA, Maritz-Olivier C. Bovine immune factors underlying tick resistance: integration and future directions. Front Cell Infect Microbiol. 2017;7:522–37.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  20. Tabor AE, Ali A, Rehman G, Rocha Garcia G, Zangirolamo AF, Malardo T, Jonsson NN. Cattle tick Rhipicephalus microplus-host interface: a review of resistant and susceptible host responses. Front Cell Infect Microbiol. 2017;7:506–23.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  21. Piper EK, Jackson LA, Bagnall NH, Kongsuwan KK, Lew AE, Jonsson NN. Gene expression in the skin of Bos taurus and Bos indicus cattle infested with the cattle tick, Rhipicephalus (Boophilus) microplus. Vet Immunol Immunop. 2008;126:110–9.

    Article  CAS  Google Scholar 

  22. Constantinoiu CC, Jackson LA, Jorgensen WK, Lew-Tabor AE, Piper EK, Mayer DG, Venus B, Jonsson NN. Local immune response against larvae of Rhipicephalus (Boophilus) microplus in Bos taurus indicus and Bos taurus taurus cattle. Int J Parasitol. 2010;40:865–75.

    Article  CAS  PubMed  Google Scholar 

  23. Piper EK, Jackson LA, Bielefeldt-Ohmann H, Gondro C, Lew-Tabor AE, Jonsson NN. Tick-susceptible Bos taurus cattle display an increased cellular response at the site of larval Rhipicephalus (Boophilus) microplus attachment, compared with tick-resistant Bos indicus cattle. Int J Parasitol. 2010;40:43–441.

    Article  Google Scholar 

  24. Carvalho WA, Domingues R, de Azevedo-Prata MC, da Silva MV, de Oliveira GC, Guimaraes SE. Microarray analysis of tick-infested skin in resistant and susceptible cattle confirms the role of inflammatory pathways in immune activation and larval rejection. Vet Parasitol. 2014;205:307–17.

    Article  CAS  PubMed  Google Scholar 

  25. Franzin AM, Maruyama SR, Garcia GR, Oliveira RP, Ribeiro JM, Bishop R, et al. Immune and biochemical responses in skin differ between bovine hosts genetically susceptible and resistant to the cattle tick Rhipicephalus microplus. Parasit Vectors. 2017;10:51–74.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  26. Piper EK, Jonsson NN, Gondro C, Lew-Tabor AE, Moolhuijzen P, Vance ME, Jackson LA. Immunological profiles of Bos taurus and Bos indicus cattle infested with the cattle tick, Rhipicephalus (Boophilus) microplus. Clin Vacc Immunol. 2009;16:1074–86.

    Article  CAS  Google Scholar 

  27. Robbertse L, Richards SA, Clift SJ, Barnard AC, Leisewitz A, Crafford JE, Maritz-Olivier C. Comparison of the differential regulation of T and B-lymphocyte subsets in the skin and lymph nodes amongst three cattle breeds as potential mediators of immune-resistance to Rhipicephalus microplus. Ticks Tick Borne Dis. 2018;9:976–87.

    Article  PubMed  Google Scholar 

  28. Wang YH, Reverter A, Kemp D, McWilliam SM, Ingham A, Davis CK, Moore RJ, Lehnert SA. Gene expression profiling of Hereford shorthorn cattle following challenge with Boophilus microplus tick larvae. Aust J Exp Agric. 2007;47:1397–407.

    Article  CAS  Google Scholar 

  29. Kongsuwan K, Josh P, Colgrave ML, Bagnall NH, Gough J, Burns B, Pearson R. Activation of several key components of the epidermal differentiation pathway in cattle following infestation with the cattle tick, Rhipicephalus (Boophilus) microplus. Int J Parasitol. 2010;40:499–507.

    Article  PubMed  Google Scholar 

  30. Piper EK, Jonsson NN, Gondro C, Vance ME, Lew-Tabor A, Jackson LA. Peripheral cellular and humoral responses to infestation with the cattle tick Rhipicephalus microplus in Santa Gertrudis cattle. Parasite Immunol. 2017;39:e12402.

    Article  CAS  Google Scholar 

  31. Constantinoiu CC, Lew-Tabor A, Jackson LA, Jorgensen WK, Piper EK, Mayer DG, Johnson L, Venus B, Jonsson NN. Local immune response to larvae of Rhipicephalus microplus in Santa Gertrudis cattle. Parasite Immunol. 2018;40:e12515.

    Article  CAS  PubMed  Google Scholar 

  32. Porto Neto LR, Bunch RJ, Harrison BE, Prayaga KC, Barendse W. Haplotypes that include the integrin alpha 11 gene are associated with tick burden in cattle. BMC Genet. 2010;11:55–67.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  33. Regitano LCA, Ibelli AM, Gasperin G, Miyata M, Azevedo AL, Coutinho LL, Teodoro RL, Machado MA, Silva MV, Nakata LC, Zaros LG, Sonstegard TS, Silva AM, Alencar MM, Oliveira MC. On the search for markers of tick resistance in bovines. Dev Biol (Basel). 2008;132:225–30.

    CAS  Google Scholar 

  34. Machado MA, Azevedo AL, Teodoro RL, Pires MA, Peixoto MG, de Freitas C, et al. Genome wide scan for quantitative trait loci affecting tick resistance in cattle (Bos taurus x Bos indicus). BMC Genomics. 2010;11:280–90.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  35. Gasparin G, Miyata M, Coutinho LL, Martinez ML, Teodoro RL, Furlong J, et al. Mapping of quantitative trait loci controlling tick [Riphicephalus (Boophilus) microplus] resistance on bovine chromosomes 5, 7 and 14. Anim Genet. 2007;38:453–9.

    Article  CAS  PubMed  Google Scholar 

  36. Sollero BP, Junqueira VS, Gomes CCG, Caetano AR, Cardoso FF. Tag SNP selection for prediction of tick resistance in Brazilian Braford and Hereford cattle breeds using Bayesian methods. Genet Sel Evol. 2017;49:49–63.

    Article  PubMed  PubMed Central  Google Scholar 

  37. Porto-Neto LR, Jonsson NN, Ingham A, Bunch RJ, Harrison BE, Barendse W. Cooperative Research Centre for Beef Genetic Technology. The RIPK2 gene: a positional candidate for tick burden supported by genetic associations in cattle and immunological response of knockout mouse. Immunogenet. 2012;64:379–88.

    Article  CAS  Google Scholar 

  38. Mapholi NO, Marufu MC, Maiwashe A, Banga CB, Muchenje V, MacNeil MD, Chimonyo M, Dzama K. Towards a genomics approach to tick (Acari: Ixodidae) control in cattle: a review. Ticks Tick Borne Dis. 2014;5(5):475–83.

    Article  PubMed  Google Scholar 

  39. Biegelmeyer P, Gulias-Gomes CC, Caetano AR, Steibel JP, Cardoso FF. Linkage disequilibrium, persistence of phase and effective population size estimates in Hereford and Braford cattle. BMC Genet. 2016;17:32–43.

    Article  PubMed  PubMed Central  Google Scholar 

  40. Mota RR, Guimarães SEF, Fortes MRS, Hayes B, Silva FF, Verardo LL, Kelly MJ, de Campos CF, Guimarães JD, Wenceslau RR, Penitente-Filho JM, Garcia JF, Moore S. Genome-wide association study and annotating candidate gene networks affecting age at first calving in Nellore cattle. J Anim Breed Genet. 2017;134:484–92.

    Article  CAS  PubMed  Google Scholar 

  41. de Souza MM, Zerlotini A, Geistlinger L, Tizioto PC, Taylor JF, Rocha MIP, Diniz WJS, Coutinho LL, Regitano LCA. A comprehensive manually-curated compendium of bovine transcription factors. Sci Rep. 2018;8:13747–58.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  42. Chmelar J, Calvo E, Pedra JH, Francischetti IM, Kotsyfakis M. Tick salivary secretion as a source of antihemostatics. J Proteome. 2012;75:3842–54.

    Article  CAS  Google Scholar 

  43. Wikel S. Ticks and tick-borne pathogens at the cutaneous interface: host defenses, tick countermeasures, and a suitable environment for pathogen establishment. Front Microbiol. 2013;4:337–46.

    Article  PubMed  PubMed Central  Google Scholar 

  44. Galay RL, Umemiya-Shirafuji R, Mochizuki M, Fujisaki K, Tanaka T. Iron metabolism in hard ticks (Acari: Ixodidae): the antidote to their toxic diet. Parasitol Int. 2015;64:182–9.

    Article  CAS  PubMed  Google Scholar 

  45. Turkish AR, Henneberry AL, Cromley D, Padamsee M, Oelkers P, Bazzi H, et al. Identification of two novel human acyl-CoA wax alcohol acyltransferases: members of the diacylglycerolacyltransferase 2 (DGAT2) gene superfamily. J Biol Chem. 2005;280:14755–64.

    Article  CAS  PubMed  Google Scholar 

  46. Zhang C, Wang K, Yang L, Liu R, Chu Y, Qin X, et al. Lipid metabolism in inflammation-related diseases. Analyst. 2018;143:4526–36.

    Article  CAS  PubMed  Google Scholar 

  47. Burger BV. Mammalian semiochemicals. Top Cur Chem. 2005;240:231–78.

    Article  CAS  Google Scholar 

  48. Brain SD. Sensory neuropeptides: their role in inflammation and wound healing. Immunopharmacology. 1997;37:133–52.

    Article  CAS  PubMed  Google Scholar 

  49. Kang HS, Trzaska KA, Corcoran K, Chang VT, Rameshwar P. Neurokinin receptors: relevance to the emerging immune system. Arch Immunol Ther Exp. 2004;52:338–47.

    CAS  Google Scholar 

  50. Carvalho WA, Bechara GH, More DD, Ferreira BR, da Silva JS, de Miranda Santos IK. Rhipicephalus (Boophilus) microplus: distinct acute phase proteins vary during infestations according to the genetic composition of the bovine hosts, Bos taurus and Bos indicus. Exp Parasitol. 2008;118:587–91.

    Article  CAS  PubMed  Google Scholar 

  51. Salmi M, Jalkanen S. A 90-kilodalton endothelial cell molecule mediating lymphocyte binding in humans. Science. 1992;257:1407–9.

    Article  CAS  PubMed  Google Scholar 

  52. Arvilommi AM, Salmi M, Kalimo K, Jalkanen S. Lymphocyte binding to vascular endothelium in inflamed skin revisited: a central role for vascular adhesion protein-1 (VAP-1). Eur J Immunol. 1996;26:825–33.

    Article  CAS  PubMed  Google Scholar 

  53. Carvalho WA, Franzin AM, Abatepaulo AR, de Oliveira CJ, More DD, da Silva JS, et al. Modulation of cutaneous inflammation induced by ticks in contrasting phenotypes of infestation in bovines. Vet Parasitol. 2010;167:260–73.

    Article  CAS  PubMed  Google Scholar 

  54. Allen JR. Tick resistance: basophils in skin reactions of resistant Guinea pigs. Int J Parasitol. 1973;3:195–200.

    Article  CAS  PubMed  Google Scholar 

  55. Wada T, Ishiwata K, Koseki H, Ishikura T, Ugajin T, Ohnuma N, Obata K, et al. Selective ablation of basophils in mice reveals their non-redundant role in acquired immunity against ticks. J Clin Invest. 2010;120:2867–75.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  56. Kimura J, Kudoh T, Miki Y, Yoshida K. Identification of dihydropyrimidinase-related protein 4 as a novel target of the p53 tumor suppressor in the apoptotic response to DNA damage. Int J Cancer. 2011;128:1524–31.

    Article  CAS  PubMed  Google Scholar 

  57. Anguera MC, Suh JR, Ghandour H, Nasrallah IM, Selhub J, Stover PJ. Methenyltetrahydrofolatesynthetase regulates folate turnover and accumulation. J Biol Chem. 2003;278:29856–62.

    Article  CAS  PubMed  Google Scholar 

  58. Williams JD, Jacobson EL, Kim H, Kim M, Jacobson MK. Folate in skin cancer prevention. Subcell Biochem. 2012;56:181–97.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  59. Ogawa K, Johansson SL, Cohen SM. Immunohistochemical analysis of uroplakins, urothelial specific proteins, in ovarian Brenner tumors, normal tissues, and benign and neoplastic lesions of the female genital tract. Am J Pathol. 1999;155:1047–50.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  60. Neves-Costa A, Moita LF. TET1 is a negative transcriptional regulator of IL-1beta in the THP-1 cell line. Mol Immunol. 2013;54:264–70.

    Article  CAS  PubMed  Google Scholar 

  61. Mora JR, Iwata M, von Andrian UH. Vitamin effects on the immune system: vitamins a and D take Centre stage. Nat Rev Immunol. 2008;8:685–98.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  62. Pino-Lagos K, Benson MJ, Noelle RJ. Retinoic acid in the immune system. Ann N Y Acad Sci. 2008;1143:170–87.

    Article  CAS  PubMed  Google Scholar 

  63. Aranda JF, Reglero-Real N, Kremer L, Marcos-Ramiro B, Ruiz-Saenz A, Calvo M, et al. MYADM regulates Rac1 targeting to ordered membranes required for cell spreading and migration. Mol Biol Cell. 2011;22:1252–62.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  64. Roberts JA. Resistance of cattle to Boophilus microplus. II. Stages of the life cycle of the parasite against which resistance manifest. J Parasitol. 1968;54:667–73.

    Article  CAS  PubMed  Google Scholar 

  65. Heinze DM, Wikel SK, Thangamani S, Alarcon-Chaidez FJ. Transcriptional profiling of the murine cutaneous response during initial and subsequent infestations with Ixodes scapularis nymphs. Parasit Vectors. 2012;5:26–40.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  66. Kim TG, Byamba D, Wu WH, Lee MG. Statins inhibit chemotactic interaction between CCL20 and CCR6 in vitro: possible relevance to psoriasis treatment. Exp Dermatol. 2011;20:855–7.

    Article  PubMed  CAS  Google Scholar 

  67. Poole NM, Mamidanna G, Smith RA, Coons LB, Cole JA. Prostaglandin E(2) in tick saliva regulates macrophage cell migration and cytokine profile. Parasit Vectors. 2013;6:261–71.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  68. Wang X, Shaw DK, Sakhon OS, Snyder GA, Sundberg EJ, Santambrogio L, et al. The tick protein sialostatin L2 binds to annexin A2 and inhibits NLRC4-mediated inflammasome activation. Infec Immun. 2016;84:1796–805.

    Article  CAS  Google Scholar 

  69. Flier J, Boorsma DM, van Beek PJ, Nieboer C, Stoof TJ, Willemze R, et al. Differential expression of CXCR3 targeting chemokines CXCL10, CXCL9, and CXCL11 in different types of skin inflammation. J Pathol. 2001;194:398–405.

    Article  CAS  PubMed  Google Scholar 

  70. Hunter CA, Jones SA. IL-6 as a keystone cytokine in health and disease. Nat Immunol. 2015;16:448–57.

    Article  CAS  PubMed  Google Scholar 

  71. Sonnenberg GF, Fouser LA, Artis D. Border patrol: regulation of immunity, inflammation and tissue homeostasis at barrier surfaces by IL-22. Nat Immunol. 2011;12:383–90.

    Article  CAS  PubMed  Google Scholar 

  72. Dumoutier L, Louahed J, Renauld JC. Cloning and characterization of IL-10-related T cell-derived inducible factor (IL-TIF), a novel cytokine structurally related to IL-10 and inducible by IL-9. J Immunol. 2000;164:1814–9.

    Article  CAS  PubMed  Google Scholar 

  73. Lejeune D, Dumoutier L, Constantinescu S, Kruijer W, Schuringa JJ, Renauld JC. Interleukin-22 (IL-22) activates the JAK/STAT, ERK, JNK, and p38 MAP kinase pathways in a rat hepatoma cell line. Pathways that are shared with and distinct from IL-10. J Biol Chem. 2002;277:33676–82.

    Article  CAS  PubMed  Google Scholar 

  74. Andoh A, Zhang Z, Inatomi O, Fujino S, Deguchi Y, Araki Y, et al. Interleukin-22, a member of the IL-10 subfamily, induces inflammatory responses in colonic subepithelial myofibroblasts. Gastroenterology. 2005;129:969–84.

    Article  CAS  PubMed  Google Scholar 

  75. Gurney AL. IL-22, a Th1 cytokine that targets the pancreas and select other peripheral tissues. Int Immunopharmacol. 2004;4:669–77.

    Article  CAS  PubMed  Google Scholar 

  76. Wolk K, Sabat R. Interleukin-22: a novel T- and NK-cell derived cytokine that regulates the biology of tissue cells. Cytokine Growth Factor Rev. 2006;17:367–80.

    Article  CAS  PubMed  Google Scholar 

  77. Komiya Y, Habas R. Wnt signal transduction pathways. Organogenesis. 2008;4:68–75.

    Article  PubMed  PubMed Central  Google Scholar 

  78. Silva-Garcia O, Valdez-Alarcon JJ, Baizabal-Aguirre VM. The Wnt/beta-catenin signaling pathway controls the inflammatory response in infections caused by pathogenic bacteria. Mediat Inflamm. 2014;2014:310183–9.

    Article  CAS  Google Scholar 

  79. Swafford D, Manicassamy S. Wnt signaling in dendritic cells: its role in regulation of immunity and tolerance. Discov Med. 2015;19:303–10.

    PubMed  PubMed Central  Google Scholar 

  80. Lim X, Tan SH, Koh WL, Chau RM, Yan KS, Kuo CJ, et al. Interfollicular epidermal stem cells self-renew via autocrine Wnt signaling. Science. 2013;342:1226–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  81. Lim X, Nusse R. Wnt signaling in skin development, homeostasis, and disease. Cold Spring Harb Perspect Biol. 2013;5:1–26.

    Google Scholar 

  82. Houschyar KS, Momeni A, Pyles MN, Maan ZN, Whittam AJ, Siemers F. Wnt signaling induces epithelial differentiation during cutaneous wound healing. Organogenesis. 2015;1:95–104.

    Article  Google Scholar 

  83. Staal FJ, Luis TC, Tiemessen MM. WNT signaling in the immune system: WNT is spreading its wings. Nat Rev Immunol. 2008;8:581–93.

    Article  CAS  PubMed  Google Scholar 

  84. You J, Nguyen AV, Albers CG, Lin F, Holcombe RF. Wnt pathway-related gene expression in inflammatory bowel disease. Dig Dis Sci. 2008;53:1013–9.

    Article  CAS  PubMed  Google Scholar 

  85. Ibelli AM, Ribeiro AR, Giglioti R, Regitano LC, Alencar MM, Chagas AC, et al. Resistance of cattle of various genetic groups to the tick Rhipicephalus microplus and the relationship with coat traits. Vet Parasitol. 2012;186:425–30.

    Article  CAS  PubMed  Google Scholar 

  86. Marufu MC, Qokweni L, Chimonyo M, Dzama K. Relationships between tick counts and coat characteristics in Nguni and Bonsmara cattle reared on semiarid rangelands in South Africa. Ticks Tick Borne Dis. 2011;2:172–7.

    Article  PubMed  Google Scholar 

  87. Li WQ, Dehnade F, Zafarullah M. Oncostatin M-induced matrix metalloproteinase and tissue inhibitor of metalloproteinase-3 genes expression in chondrocytes requires Janus kinase/STAT signaling pathway. J Immunol. 2001;166:3491–8.

    Article  CAS  PubMed  Google Scholar 

  88. Biegelmeyer P, Nizoli LQ, da Silva SS, dos Santos TR, Dionello NJ, Gulias-Gomes CC, et al. Bovine genetic resistance effects on biological traits of Rhipicephalus (Boophilus) microplus. Vet Parasitol. 2015;208:231–7.

    Article  CAS  PubMed  Google Scholar 

  89. Zhbannikov IY, Hunter SS, Settles ML. SEQYCLEAN user manual. New York, NY: Elsevier; 2013.

    Google Scholar 

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

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  91. Robinson MD, McCarthy DJ, Smyth GK. edgeR: a bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010;26:139–40.

    Article  CAS  PubMed  Google Scholar 

Download references

Acknowledgments

The authors acknowledge the Delta G Connection/Gensys breeding program for providing the animals and performance data for this research, Elisângela Guedes for helping in experimental data collection and Rebecca L Coutinho for proofreading the manuscript.

Funding

Research was financially supported by EMBRAPA (Brazilian Agricultural Research Corporation) grants 02.09.07.004 and 02.13.10.002, which provided infrastructure, consumables and personnel for the field experiment and laboratory analysis. Personnel and consumables were also supported by CAPES (Coordination for the Improvement of Higher Level Personnel) grant PNPD 02645/09–2 and FAPESP (São Paulo Research Foundation) grant PD 2012/24418–1. The financial institutions had no participation on experimental design, data collection and analysis, nor contributed in writing the manuscript.

Author information

Authors and Affiliations

  1. EMBRAPA Pecuária Sudeste, São Carlos, São Paulo, Brazil

    Daniela D. Moré, Wilson Malagó-Jr & Luciana C. A. Regitano

  2. EMBRAPA Pecuária Sul, Bagé, Rio Grande do Sul, Brazil

    Fernando F. Cardoso, Claudia C. Gulias-Gomes & Bruna P. Sollero

  3. Federal University of Pelotas, Capão do Leão, Rio Grande do Sul, Brazil

    Fernando F. Cardoso

  4. EMBRAPA Informática Agropecuária, Campinas, São Paulo, Brazil

    Maurício A. Mudadu

  5. EMBRAPA Suínos e Aves, Concórdia, Santa Catarina, Brazil

    Adriana M. G. Ibelli

  6. Luiz de Queiroz College of Agriculture, University of São Paulo, Piracicaba, São Paulo, Brazil

    Luiz L. Coutinho

Authors
  1. Daniela D. Moré
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Fernando F. Cardoso
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Maurício A. Mudadu
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Wilson Malagó-Jr
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Claudia C. Gulias-Gomes
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Bruna P. Sollero
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Adriana M. G. Ibelli
    View author publications

    You can also search for this author in PubMed Google Scholar

  8. Luiz L. Coutinho
    View author publications

    You can also search for this author in PubMed Google Scholar

  9. Luciana C. A. Regitano
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

FFC, CCG-G, AMGI and LCAR conceived and designed the experiments. CCG-G, AMGI and WM-J performed the experiments. DDM, FFC, MAM, and BPS analyzed the data. LLC collaborated with reagents. DDM, FFC, LLC and LCAR wrote the paper. All authors reviewed and approved the final version of the manuscript.

Authors’ information

FFC, LLC, and LCAR are research fellows of the National Council for Scientific and Technological Development (CNPq). DDM was a postdoctoral fellow of FAPESP.

Corresponding author

Correspondence to Luciana C. A. Regitano.

Ethics declarations

Ethics approval and consent to participate

The Committee for Ethics in Animal Experimentation of the Federal University of Pelotas, Rio Grande do Sul, Brazil, approved all experimental procedures involving animals under registration number 6849. The animals used in this study were obtained from Delta G Connection Breeding Program affiliated producers under lease and consent agreements (EMBRAPA SAIC system registration numbers: 21900.11/0003–8 - Carlos Edmundo Eichenberg; 21900.11/0014–5 - Pedro Monteiro Lopes, and 21900.11/0004–6 - Alfeu Medeiros Fleck).

Consent for publication

Not applicable.

Competing interests

The authors declare they have no competing interests.

Additional information

Publisher’s Note

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

Supplementary information

Additional file 1: Table S1.

Sequencing throughput and mapping statistics for each resistance group. The average, minimum and maximum values of numbers (N) of reads passing filtering and aligned pairs, and the percentage (%) of overall read mapping, multiple alignments, and concordant pair alignment rate are shown. Values were compared inter- and intragroups. Different letters (a and b) represent statistical significance. (XLS 28 kb)

Additional file 2: Table S2.

Intragroup comparisons of differentially expressed genes in the skin of tick-resistant and -susceptible Braford cattle. DEGs are listed alphabetically with their Ensembl ID and description. The log2 FC and the FDR for each DEG are shown in the last two columns, respectively, with the results comparing Rpost vs. Rpre represented by white cells and Spost vs. Spre represented by gray cells. (XLS 866 kb)

Additional file 3: Table S3.

Enrichment analysis report of TopDEGs from intragroup comparisons. A) Enrichment by process networks of TopDEGs from resistant (Rpost vs. Rpre) and susceptible (Spost vs. Spre) hosts. Networks were classified according to the p-value, the false discovery rate (FDR), and the number of network objects found in the TopDEGs list (In Data). Total represents the number of total network objects in each network. Unique TopDEGs are shown in bold. B) Enrichment by pathway maps of TopDEGs from resistant (Rpost vs. Rpre) and susceptible (Spost vs. Spre) hosts. Maps were classified following the same scheme described for networks. (XLS 39 kb)

Additional file 4: Table S4.

Network prediction for TopDEGs from resistant response to tick infestation. For each network, the key network objects and the GO process associated with them, followed by the number of total and seed nodes, and pathways are shown. The p-value, z-score, and g-score classified the networks. (XLS 42 kb)

Additional file 5: Table S5

Network prediction for genes mined from the TagSNPs associated with resistance to ticks prospected in a genome-wide association study. For each network, the key network objects and the GO process associated with them, followed by the number of total and seed nodes, and pathways are shown. The p-value, z-score, and g-score classified the networks. (XLS 34 kb)

Additional file 6: Table S6.

Network prediction for TopDEGs from susceptible response to tick infestation. For each network, the key network objects and the GO process associated with them, followed by the number of total and seed nodes, and pathways are shown. The p-value, z-score, and g-score classified the networks. (XLS 40 kb)

Additional file 7: Table S7.

Transcription factors (TFs) prediction for resistant hosts. Network analysis predicted potential TFs underlying the regulation of gene expression, based on the variation of TopDEGs. The output file was manually-curated, and only TFs found in cattle were classified according to their z-score. TopDEGs TFs were also indicated with their respective modulation and statistical significance. Network object name: network object name in MetaBase®; Actual: number of network objects in the activated dataset(s) which interact with the chosen object; n: number of network objects in the activated dataset(s); R: number of network objects in the complete database or background list which interact with the chosen object; N: total number of gene-based objects in the complete database or background list; Expected: mean value for hypergeometric distribution (n*R/N); Ratio: connectivity ratio (Actual/Expected); z-score: ((Actual-Expected)/sqrt(variance)); p-value: probability to have the given value of Actual or higher (or lower for negative z-score); Input IDs: original probe/gene IDs in the activated dataset(s); Signal: variation of gene expression in log2 FC; p-value2: statistical significance in the differential expression analysis. (XLS 52 kb)

Additional file 8: Table S8.

Transcription factors (TFs) prediction for susceptible hosts. Network analysis predicted potential TFs underlying the regulation of gene expression, based on the variation of TopDEGs. The output file was manually-curated, and only the TFs found in cattle were classified according to their z-score. TopDEGs TFs were also indicated with their respective modulation and statistical significance. Network object name: network object name in MetaBase®; Actual: number of network objects in the activated dataset(s) which interact with the chosen object; n: number of network objects in the activated dataset(s); R: number of network objects in the complete database or background list which interact with the chosen object; N: total number of gene-based objects in the complete database or background list; Expected: mean value for hypergeometric distribution (n*R/N); Ratio: connectivity ratio (Actual/Expected); z-score: ((Actual-Expected)/sqrt(variance)); p-value: probability to have the given value of Actual or higher (or lower for negative z-score); Input IDs: original probe/gene IDs in the activated dataset(s); Signal: variation of gene expression in log2 FC; p-value2: statistical significance in the differential expression analysis. (XLS 48 kb)

Additional file 9: Table S9.

Estimates of statistic parameters used to select genetically divergent hosts. Statistic estimates of the mean (mCount), mean standard deviation (mSD), minimum (minCount) and maximum (maxCount) tick counts, and estimated breeding value (EBV) are shown individually for resistant (R) and susceptible (S) hosts. (XLS 32 kb)

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

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Moré, D.D., Cardoso, F.F., Mudadu, M. et al. Network analysis uncovers putative genes affecting resistance to tick infestation in Braford cattle skin. BMC Genomics 20, 998 (2019). https://doi.org/10.1186/s12864-019-6360-3

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12864-019-6360-3

Keywords

BMC Genomics

ISSN: 1471-2164

Contact us