Skip to main content
Download PDF

Extensive differential DNA methylation between tuberculosis skin test positive and skin test negative cattle

BMC Genomics volume 25, Article number: 762 (2024) Cite this article

Abstract

Bovine tuberculosis (bTB), caused by Mycobacterium bovis (M. bovis), represents a significant problem for the agriculture industry as well as posing a risk for human health. Current diagnostic tests for bTB target the cell-mediated immune (CMI) response to infection with M. bovis, primarily through screening of animals with the tuberculin skin test. Epigenetic modifications have been shown to alter the course of the immune response and differentially methylated regions (DMRs) might also influence the outcome of the skin test in cattle. Whole Genome Bisulphite Sequencing (WGBS) was used to profile DNA methylation levels from peripheral blood of a group of cattle identified as test positive for M. bovis (positive for the single intradermal comparative tuberculin test (SICTT) and/or the interferon-γ release assay compared to a test negative control group [n = 8/group, total of 16 WGBS libraries]. Although global methylation profiles were similar for both groups across the genome, 223 DMRs and 159 Differentially Promoter Methylated Genes (DPMGs) were identified between groups with an excess of hypermethylated sites in SICTT positive cattle (threshold > 15% differential methylation). Genes located within these DMRs included the Interleukin 1 receptor (IL1R1) and MHC related genes (BOLA and BOLA-DQB). KEGG pathway analysis identified enrichment of genes involved in Calcium and MAPK signalling, as well as metabolism pathways. Analysis of DMRs in a subset of SICTT negative cattle that were IFN-γ positive showed differential methylation of genes including Interleukin 10 Receptor, alpha (IL10RA), Interleukin 17 F (IL17F) and host defence peptides (DEFB and BDEF109). This study has identified a number of immune gene loci at which differential methylation is associated with SICTT test results and the degree of methylation could influence effective host immune responses.

Peer Review reports

Background

Infection by Mycobacterium bovis (M. bovis) continues to cause significant animal health issues globally with reservoirs of infection in both livestock and wildlife populations [1]. In addition to the significant costs associated with control and eradication programmes, M. bovis is a zoonotic pathogen and therefore has consequential trade and human health impacts [2]. Although the overall proportion of cattle herds testing positive for bTB across the EU remains low, certain countries are experiencing increasing incidences in recent years. As a result, intensive eradication efforts are in place in countries across the EU, including Ireland [3].

The current eradication programme for bTB in Ireland is underpinned by annual screening of cattle over six weeks of age with the Single Intradermal Comparative Tuberculin Test (SICTT), an in vivo skin test that measures the activation of a Delayed Type Hypersensitivity (DTH) reaction within 72 h after the administration of both bovine and avian mycobacterial tuberculin antigens (M. bovis/M. avium purified protein derivative – PPDb/PPDa) [4]. PPDa is used to control for the background sensitisation of animals with cross-reactive environmental mycobacteria and is applied to increase test specificity. A positive standard ‘reactor’ animal is disclosed when the DTH reaction of the response to PPDb is at least 4 mm (B-A > 4 mm) greater than the response to PPDa. More severe interpretations (e.g., B-A > 2 mm) can help to increase the sensitivity of the SICTT. In herds where bTB reactor animals are identified, the interferon- gamma assay (IFN-γ) is often used as an ancillary test, interpreted in parallel with the SICTT to improve the sensitivity of diagnostic testing. The principle of the IFN-γ assay is to detect and quantify release of the IFN-γ cytokine when heparinised whole blood is incubated with PPDb and PPDa [5]. All reactor animals are euthanised soon after disclosure and subjected to post-mortem examination at the abattoir for the presence of lesions consistent with M. bovis infection. However, lesion detection in reactor animals is relatively low: several studies have found that 50–80% of reactor animals had no visible lesions [6] and in a study carried out in Northern Ireland, only 43% of reactors had visible lesions detected at slaughter [7]. Nevertheless, the failure to detect visible lesions at post-mortem does not necessarily indicate absence of infection. The SICTT has a median specificity approaching 100% based on studies of bovine tuberculosis (bTB)-free populations from several countries [5, 6, 8, 9]. Therefore, when multiple SICTT reactors are identified in an exposed herd by the SICTT and / or the IFN-γ test, the probability that they are false-positive reactors is very low. They are deemed to be infected with M. bovis, even in the absence of confirmation by pathological examination. Likewise, in low-risk of infection cohorts, two negative tests (SICTT and the IFN-γ test) in any animal indicate a low probability of bTB infection. Animals that are positive in the IFN-γ test are also considered at higher risk of infection than IFN-γ test negatives [10,11,12].

In an effort to improve on the performance characteristics of current available diagnostic tests, the advent of new technologies, including transcriptomics, have also served to increase our understanding of the immune response in bTB infected cattle [13, 14]. These studies have revealed a complex array of gene expression patterns that govern the specific responses to infection with M. bovis [15, 16]. An emerging area of research has also shown that these responses are underpinned by epigenetic processes that modify chromatin structure to facilitate the coordinated transcription of genes involved in the specific immune responses to infection [17, 18].

Epigenetics refers to chemical modifications to DNA and proteins which control access by the transcriptional machinery to the underlying genetic sequence and is relevant to all phenotypic traits in verterbrates [19, 20]. As each step in the immune response to bTB relies on gene activation prior to the translation of effector proteins, epigenetic modifications are involved in regulation of the IFN-γ and ultimately DTH responses [21]. It has also been shown that infectious organisms, including mycobacteria have evolved means to manipulate the epigenome of their host to facilitate their survival [22, 23]. One major epigenetic mechanism is methylation of DNA through the addition of a methyl group to the cytosine residue of DNA, often in the context of a CpG island [24]. However, other forms of asymmetrical methylation have also been identified in mammals, where methylation occurs to the cytosine in an alternative context - namely CHG or CHH (where H correspond to A, T or C) [25, 26]. The effects of methylation on the immune response will depend on its location: in the promoter region of some genes, for example, methylation can repress the binding of transcription factors potentially downregulating gene expression. Studies of cases with human tuberculosis (caused by M. tuberculosis) have identified specific loci at which hypermethylation results in a dampening of the immune response, including at the IFNGR1 and IL12BR2 genes [27]. Both IFN-γ and IL-12 cytokines are associated with the development of protective immunity to mycobacteria [28] and therefore hypermethylation of these and/or other genes could impact on the outcome of infection as well as the performance of diagnostic tests based on host responses.

In this study, we set out to determine if cattle stratified by routine field interpretation of their SICTT and IFN-γ test responses exhibit differentially methylated immune genes and pathways which underpin their divergent immune responses profiles.

Results & discussion

Summary of whole genome bisulphite sequencing (WGBS) data

Raw WGBS data from 16 whole blood samples representing bTB SICTT positive (n = 8) and negative (n = 8) cattle, each comprising of eight Holstein cows were assessed (Fig. 1 and Supplementary Table S1). After quality control, a total of 5.01 × 109 clean reads were obtained with an average of 3.13 × 108 (average 85.48 Gb per sample) for each sample, indicating a sequencing depth of 30X. GC content was 22.68% on average and the bisulphite conversion rate was 99.72% on average (Supplementary Table S2). On average, 78.88% of non-duplicated clean reads uniquely aligned to the bovine reference genome (ARS-UCD1.2) with an average duplication rate of 24%.

Fig. 1
figure 1

Experimental Design. Four cohorts of samples were identified based on diagnoses using the SICTT test results and IFN-γ assay results from naturally bTB + cattle or test negative controls. Eight SICTT negative and eight SICTT positive cattle were selected and further classified as either SICTT negative and IFN-γ negative (Group 1, n = 4); SICTT negative and IFN-γ positive (Group 2, n = 4); SICTT positive and IFN-γ positive (Group 3, n = 4); SICTT positive and IFN-γ negative (Group 4, n = 4)

Full size image

DNA methylation levels and patterns across the bovine genome

On average, 3.38% of all cytosines in the genome were methylated and 64% of methylated cytosines occurred in a CG context (Supplementary Table S2). As expected, and previously documented in the literature, CG methylation representing the majority of methylation changes present, at 98% of methylated cytosines on average. Asymmetric methylation (CHG and CHH - where H is A, C, or T) accounted for a fraction of overall methylation changes, representing an average of 0.51% and 1.56% of methylated cytosines, respectively. No statistical differences in the proportion of CG, CHG or CHH methylated cytosines was detected between groups (P > 0.05). Across genomic features, clear, reproducible differences in methylation profiles was apparent across all samples with lowest CG methylation in the 5’UTR region of genes (20% on average) rising to maximal levels in exons (71%), introns (81%) and 3’ UTR regions (81%) (Supplementary Fig. S1A and S1B and Supplementary Table S2). The average levels of methylation across genomic features is remarkably similar to what has been reported in cattle recently, with 30% in the promoter and 5’ untranslated region (UTR), 68% in exons, 72% in introns, and 73% in the 3’ UTR [29]. Asymmetrical methylation patterns mirrored the CG profile in general but specific peaks of occurrences were apparent. At the gene level, the gene body accounted for the majority of CG methylation with peaks for CHG and CHH methylation immediately downstream of genes (Supplementary Fig. S2). Low levels of asymmetric methylation levels are also in agreement with what has previously been found in blood DNA samples from pigs [30]. The relevance of methylation in intronic and repeat regions of DNA have less clear relevance to disease processes, whereas methylation has a known role in promoter regions through the regulation of transcription factor binding and ultimately on gene expression [31]. Furthermore, methylation within gene exons is emerging as an important mechanism regulating alternative splicing [32] showing established functional relevance of methylation in these particular genomic regions.

Comparison 1 - Whole genome methylation profile in SICTT positive and SICTT negative cattle reveals differential methylated regions and elevated asymmetrical methylation in SICTT positive cattle

The methylation patterns of SICTT positive and negative cattle followed a similar pattern although differences were apparent and located widely across the genome (Fig. 2A). While overall levels of CG methylation were not statistically different between groups (Fig. 2B), the asymmetrical methylation patterns were quite distinct. Despite accounting for a minor fraction of the overall genomic methylation, both CHG and CHH methylation levels were elevated in the SICTT positive cattle DNA samples relative to the SICTT negative controls (Fig. 2C and D, respectively). Most obvious divergence in both CHG and CHH methylation levels was apparent in the CGI shores, gene promoter, exons and intronic regions. The role of asymmetric methylation is not fully understood [33], particularly in non-model organisms and while the levels reported in the bovine genome are low, the differential pattern in diseased cattle might be an important avenue for further investigation.

Fig. 2
figure 2

Differential methylation across the genome of bTB + cattle samples relative to test negative control samples. (A) Circos density plot of 5-methylcytosine density and difference across the genome. Chromosome numbers and scales are indicated on the periphery, a dark to light colour indicates a low to high level of methylation. Red indicates the greatest degree of methylation difference; (B) Differential CG methylation profile for bTB + cattle (shown in blue) and test negative controls (shown in red) across various genomic features; (C) Differential CHG methylation profile for bTB + cattle (shown in blue) and test negative controls (shown in red) across various genomic features; (D) Differential CHH methylation profile for bTB + cattle (shown in blue) and test negative controls (shown in red) across various genomic features. Both methylation types C and D are referred to as asymmetrical methylation

Full size image

Considering the direction of methylation differences across genomic features, hypermethylated and hypomethylated sites are spread across the genome (Fig. 3A). However, there is an overabundance of hypomethylated CG sites particularly in exons and introns (Fig. 3B). Although CHG and CHH methylation is a lot less abundant (Fig. 3C and E), an excess of hypomethylated sites is apparent in genomic regions including the promoter (Fig. 3D) and is most evident for CHH methylation patterns (Fig. 3F). From all identified DMRs, 8,410 mapped uniquely and were different between SICTT positive and negative controls across the genome (full list provided in Supplementary Table S3). Many of these DMRs represent extensive sequential runs of methylated cytosines spanning multiple genomic features (promoter, exon and introns) of the same gene, indicating extensive regulation of these genes by methylation (Supplementary Table S3). Various levels of stringency can be applied to increase the detection of DMRs with a potential functional impact including the areaStat and differential methylation values. The areaStat is the sum of the test statistics of all CpG sites within a DMR [34] and values in this study varied from − 1779 to 10,426 with a higher number associated with longer DMRs and are thought to indicate a more reliably detected DMR. However, various combinations of these parameters have been used in the literature to increase the stringency of DMR identification and we have used a similar approach to that previously adopted [35, 36] and focused on DMRs with a differential methylation level > 15% between groups for discussion, a level above which a biological relevance in a disease context is more likely.

Fig. 3
figure 3

Ratio of hypermethylation to hypomethylation across the genome of bTB + cattle samples relative to test negative control samples. A circos density plot illustrates the degree of hypermethylation (shown in red) and hypomethylation (shown in blue) for (A) 5-methylcytosine (CG) methylation; (C) CHG methylation and (E) CHH methylation. Chromosome numbers and scales are indicated on the periphery. The ratio of hypermethylation to hypomethylation in each of the CG, CHG and CHH sites is shown for each genomic feature in (B), (D) and (F) respectively

Full size image

A total of 1507 DMRs were identified in exons, which is reduced to 223, when the minimum 15% differential methylation cut off is applied (Table 1A). The range in differential methylation varied between − 0.49 and 0.46 and consecutive CG sites varied from a low of 4 to a maximum of 420 (median of 18 CpGs). In addition, a total of 833 DPMGs were identified, again reduced to 159 with > 15% differential methylation. Differential methylation levels in gene promoters were higher than in exons, varying from − 0.54 to 0.5 and the CpG sites were shorter than those found in exons, with a maximum of 119 CpGs in any DMR (median 12 sites). In both exons and gene promoters, the numbers of hypermethylated sites exceeded that of hypomethylated sites (124 and 90, compared to 99 and 69 respectively, threshold > 15%) in bTB + cattle compared to test negative controls (Table 1A).

The number of genes identified within all DMRs is shown in Fig. 4A and B. The majority of these genes (3268 DMR genes and 749 DPMGs) are CG methylated and the top ranked differentially methylated DMR genes and DPMGs in bTB + cattle relative to test negative controls are shown in Table 1B and Table 1C. Amongst these hypermethylated genes is the Leukocyte Associated Immunoglobulin Like Receptor (LAIR1) gene, which is a regulator of hypersensitivity reactions [37], and shows 26% higher methylation in SICTT positive cattle samples. Immunoglobulin superfamily 6 (IGSF6), which is part of panel of genes proposed to differentiate between human TB patients at different stages of infection [38] shows 21% higher methylation. Immune responses to mycobacterial infection requires the presentation of antigen by innate immune cells via the major histocompatibility complex (MHC) [39] and 39% hypermethylation of the BoLa class II histocompatibility antigen (BOLA-DQB) gene may lead to sub-optimal CD4+ T cell activation (Table 1B). Mycobacterial peptides are also presented via MHC class I, and the 23% lower methylation of the BOLA class I histocompatibility antigen, alpha chain (Table 1C) in the SICTT positive cattle may preferentially lead to a CD8+ T cell response [40]. The gene showing the lowest levels of methylation is the Interleukin 1 receptor type 1 (IL1R1) with 46% lower methylation and of note, Interleukin-1 receptor-like 2 (IL1RL2) is also hypomethylated (19%) showing important regulation of the IL-1 signalling pathway by methylation in SICTT positive cattle. A smaller number of genes within the DMRs exhibit either exclusively CHH and CHG methylation alone or in combination with CG methylation. Three genes, namely Capping actin protein (CAPZB), the Retinoic acid receptor beta (RARB) and Kinesin family member 1 A (KIF1A) show all three types of methylation simultaneously. Interestingly, evidence has just emerged that the retinoic acid receptor pathway is exploited by mycobacteria for their survival [41] but the implications of multiple forms of cytosine methylation for gene expression remains unknown. A full list of all genes is given in Supplementary Table S4.

Fig. 4
figure 4

Numbers of genes identified as differentially methylated (CG, CHG or CHH) in (A) genes or (B) in the promoter region of genes of bTB + cattle samples relative to test negative control samples. The degree of overlap shows the numbers of genes with multiple forms of methylation. The identity of these genes is given in Supplementary Table S4

Full size image

In order to identify pathways enriched by differentially methylated genes, GO and KEGG pathway analysis was conducted on genes identified in DMR and on DPMGs. No GO terms survived the P value adjustment (Supplementary Table S5), but 123 and 22 KEGG pathways were significantly enriched by genes in DMRs and DPMGs, respectively (FDR-P < 0.05). The top enriched pathways are shown in Table 2. The MAPK signalling pathway (bta04010) was the most significantly enriched by DMR genes (FDR-P = 4.29 × 107), with Calcium signalling pathway (bta04020) significantly enriched by both DMR genes and DPMGs. Inhibition of calcium signalling has been identified as an immune evasive mechanism of Mycobacterium tuberculosis [42] and MAP kinases are also critical mediators of Interferon-γ production [43]. Metabolic pathways (bta01100) were also enriched by the greatest number of DMR genes (FDR-P = 2.4 × 105) (Fig. 5). The enrichment of genes within the metabolic pathways including Inositol phosphate (bta00562), Sphingolipid (00600) and the mTOR signalling pathway (bta04150) also points towards important changes in the methylation of genes involved in nutrient partitioning which will ultimately impact on the efficacy of anti-mycobacterial immune responses [44]. A full list of the enriched KEGG pathways is shown in Supplementary Table S6.

Fig. 5
figure 5

Scatterplot showing significantly enriched pathways from (A) DMRs and (B) DPMGs identified using KEGG. FDR corrected P values are indicated by the colour and the numbers of genes represented in the pathway is indicated by the size of the circle. The Rich factor is the ratio of differentially expressed gene numbers annotated in this pathway term relative to all gene numbers annotated in this pathway term

Full size image
Table 1 Differentially methylated genes in bTB + cows relative to test negative controls. (A) summary information on DMRs across exons and gene promoters (all data after line 2 uses a stringent Diff Meth minimum of 15%); (B) top annotated genes (ranked by Diff Meth) hypermethylated in bTB + cattle for DMR and DPMG datasets; (C) top annotated genes (ranked by Diff Meth) hypomethylated in bTB + cattle for DMR and DPMG datasets. Number of CpG sites (nCG); % Differential methylation (% Diff Meth) and areaStat (the sum of the test statistic of all CG sites within the DMR)
Full size table
Table 2 The top-ranked biological pathways enriched for genes within differentially methylated regions (DMRs) and differentially promoter methylated genes (DPMGs) associated with M. bovis infection status assessed using KEGG
Full size table

The number of differentially methylated genes identified in SICTT positive cattle in this study is similar to that recently reported for cattle with Johne’s disease [29]. Johne’s disease is caused by a related Mycobacterium sp. (Myobacterium avium subsp. paratuberculosis) and a total of 3,911, 4,336, and 4,094 DMGs were detected in clinical vs. subclinical, clinical vs. healthy, and subclinical vs. healthy groups, respectively. The same study reported some similarities in terms of differentially methylated genes to this current work (IGF2 and IGF1R), as well as enriched Calcium signalling pathways which suggests that the two related diseases may share similar patho-epigenetic mechanisms [29].

Comparison 2 - Differential methylation in SICTT negative cattle displaying antigen specific IFN-γ responses (IFN-γ positive)

We then considered the methylation differences between Group 2 (SICTT-/IFN-γ+) samples relative to Group 3 (SICTT+/IFN-γ+) samples (n = 4/group). As IFN-γ positive cattle are considered at higher risk of bTB infection [12], we hypothesised that differential methylation at specific immune gene loci might influence the immune processes in these cattle. Considering the direction of methylation differences across genomic features, Group 2 samples showed an abundance of CG hypermethylation at most genomic features (Fig. 6A), a profile that was reduced or absent when examining CHG or CHH methylation patterns (Fig. 6B and Figure C, respectively). In addition, this profile contrasts with that shown when all SICTT positive and test negative controls are compared as shown in Fig. 3.

Fig. 6
figure 6

Numbers of genes identified as differentially methylated (CG, CHG or CHH) in (A) genes or (B) in the promoter region of genes of Group 2 (SICTT-/IFN-γ+) samples relative to Group 3 (SICTT+/IFG-γ+) samples. The degree of overlap shows the numbers of genes with multiple forms of methylation. The identity of these genes is given in Supplementary Table S8. The ratio of hypermethylation to hypomethylation in each of the CG, CHG and CHH sites is shown for each genomic feature in (B), (D) and (E) respectively

Full size image

A total of 4,600 unique DMRs were identified between Group 2 (SICTT-/IFN-γ+) samples relative to Group 3 (SICTT+/IFN-γ+) samples across the genome (full list given in Supplementary Table S7). Within exons, 528 of which had a differential methylation value > 15%. These DMRs had a minimum length of 4 CG nucleotides and a maximum of 947 (median 29). The range in differential methylation across these sites was (-0.58 to 0.54) and the areaStat values varied from − 1384 to 9721. In gene promoters, differential methylation levels are higher (range − 0.61 to 0.65) with 823 DMRs in total, 358 with differential methylation value > 15%. Summary details for both exonic and promoter DMRs is given in Table 3A. In both exons and gene promoters, the numbers of hypermethylated sites exceeded that of hypomethylated sites (332 and 203, compared to 196 and 155 respectively, threshold > 15%) in Group 2 (SICTT-/IFN-γ+) samples relative to Group 3 (SICTT+/IFN-γ+) cattle (Table 3A). This hypermethylation is particularly evident in promoters and exonic genomic features in a CG context (Fig. 6A-C).

Table 3 Differentially methylated regions (DMRs) and differentially promoter methylated genes (DPMGs) in Group 2 (SICTT-/IFN-γ+) cows relative to Group 3 (SICTT+/IFN-γ+). (A) summary information on differentially methylated regions across exons and gene promoters (all data after line 2 uses a stringent Diff Meth minimum of 15%); (B) top 10 annotated genes (ranked by Diff Meth) located within DMRs that are differentially methylated in Group 2 cattle relative to Group 3; (C) top 10 annotated genes (ranked by Diff Meth) located within DPMGs that are differentially methylated in Group 2 cattle relative to Group 3. Number of CpG sites (nCG); % Differential methylation (% Diff Meth) and areaStat (the sum of the test statistic of all CG sites within the DMR)
Full size table

The number of genes identified within the DMRs and DPMGs is shown in Fig. 6D and E. As previously, the majority of these genes (2722 exonic genes and 727 promoter genes) are CG methylated the top ranked differentially methylated genes in DMRs and DPMGs shown in Table 3B and Table 3C. Amongst these is the Interleukin 1 receptor (IL1R1) with 43% higher methylation in Group 2 cattle which may affect inflammatory signalling in these cattle. In addition, BOLA class I histocompatibility antigen, alpha chain (BOLA) is also hypermethylated (25%) in Group 2 samples, which may alter CD8+ T cell responses. Interestingly, the gene encoding the DEP domain containing MTOR interacting protein (DEPTOR) exhibits 32% higher methylation in Group 2 cattle, which may negatively regulate kinase activity [45]. The Interleukin-10 receptor subunit alpha gene (IL10RA) is also hypermethylated by 22% indicating a potential reduction in the regulatory effects of this cytokine. In contrast, the gene encoding Interleukin 17F (IL17F) is hypomethylated in Group 2 samples by 30% (Table 3C). IL-17 is a crucial cytokine with potent inflammatory and anti-mycobacterial functions [46] and hypomethylation suggests increased expression in Group 2 cattle. Alongside IL-17A, IL-17F has been shown to be potently induced in response to mycobacterial antigens in M. bovis infected cattle [47] and is associated with a protective response, including post-vaccination. Host defence genes with well characterised antimicrobial and immunomodulatory roles during mycobacterial infection [48] including defensin genes (BDEF109) and (DEFB), which were hypomethylated by 31% and 18% respectively.

Short, stable, microRNAs (miRNA) have also been investigated as potential biomarkers for multiple mycobacterial infections in cattle and they are proposed to have stage specificity for disease diagnosis [49]. DMRs identified in this study encompass multiple regions regulating miRNA expression which may offer insights into potential diagnostics to identify similar Group 2 cattle. The following miRNA are hypermethylated in Group 2 cattle (bta-mir-126 (17%), bta-mir-2447 (22%), bta-mir-2285ab (34%), bta-mir-11,993 (15%), bta-mir-2385 (16%), bta-mir-2309 (20%), bta-mir-7861 (24%). In contrast, others are hypomethylated (bta-mir-671 (16%), bta-mir-2887-2 (39%), and bta-mir-331 is hypomethylated by 43% (Table 3C and Supplementary Table S7). A smaller number of genes within the DMRs exhibit either exclusively CHH and CHG methylation alone or in combination with CG methylation (Supplementary Table S8).

To investigate the potential associations of the DMGs and DPMGs with IFN-γ status, functional enrichment analysis was performed. While only a single GO term survived the P value adjustment for multiple testing, identifying enrichment of the molecular function ‘Binding’ (GO:0005488, FDR-P = 9.1 × 103) (Supplementary Table S9). In contrast, 132 and 2 KEGG pathways were enriched by identified DMR genes and DPMGs including the Calcium signalling pathway (FDR-P = 2.86 × 10− 9) and the MAPK signalling pathway (FDR-P = 2.11 × 10− 7) (Table 4) which are of major relevance to mycobacterial infection. The Metabolic signalling pathway contained an input of 205 genes and was significantly enriched (FDR-P = 5.23 × 10− 8). Other significantly enriched pathways include the NFkB signalling pathway (bta04064, FDR-P = 3.7 × 10− 3) and Th17 cell differentiation pathway (bta04659, FDR-P = 2.9 × 10− 3). A full list of the enriched KEGG pathways is shown in Supplementary Table S10. For reference, data for all pairwise comparisons of the experimental groups shown in Fig. 1 is supplied in Supplementary Figures S2–S6.

Table 4 The top-ranked biological pathways enriched for genes within differentially methylated regions (DMRs) and differentially promoter methylated genes (DPMGs) associated with Group 2 (SICTT-/IFN-γ+) cows relative to Group 3 (SICTT+/IFN-γ+) cows assessed using KEGG.
Full size table

Conclusions

Deciphering methylation patterns across the entire genome is complex, and consequently the understanding of the epigenome in livestock species is in its early stages. Current diagnostics for bTB rely largely on tuberculin skin test responses as a measure of the probability of infection and it is plausible that atypical immune responses in some cattle may contribute to the imperfect sensitivity and specificity of these tests. Under natural infection conditions where physiological demands on cattle are acute and concurrent infections often occur, aberrant methylation of DNA might affect immune responses, diagnostics and disease outcomes. In addition, mycobacteria share a long evolutionary history with their respective hosts, so it is not unexpected that differential methylation may be manipulated by pathogens to enhance their survival [22, 23]. In this study we have identified extensive differential methylation of host DNA, specifically of genes regulating MAPK, Calcium and Metabolism signalling pathways in cattle with different SICTT test results and identified specific immune gene loci that could potentially affect immune responses and diagnostic test performance. This study provides a new layer of understanding to the complexities of host-pathogen interaction during infection with Mycobacterium bovis [50]. However, consensus on the optimal thresholds for filtering whole genome methylation data has yet to be agreed [51] and as the sample size in the current study is limited, further detailed validation of disease-associated methylation changes in target genes in larger populations will be required to identify specific loci associated with infection which may affect diagnostic test performance.

Materials and methods

Cattle tuberculin skin test, IFN-γ assay and DNA extraction

The SICTT was carried out by intradermal injection of cattle with 0.1 mL PPD-bovine and PPD-avian (Thermo-Fisher Scientific) at sites 12 cm apart in the mid-neck region using a McLintock tuberculin syringe. Skin thicknesses were measured in mm at both sites before the intradermal injection and after 72 h in accordance with Council Directive 64/432/EEC (2015) and OIE (2009). Based on the results of the SICTT, the animal was defined as a standard reactor if the bovine reaction was both positive and exceeded the avian reaction by > 4 mm; as a standard inconclusive reactor if the bovine reaction was either positive or inconclusive, > 1–4 mm above the avian reaction. Peripheral blood was collected via the coccygeal vein in heparin-treated vacutainers from known bTB infected herds maintained by the Irish Department of Agriculture, Food and the Marine (DAFM) as part of the national bTB control programme and were not euthanised as part of this study. Farmer consent was therefore not required. Blood samples were stimulated with PPDb and PPDa tuberculin (Thermo-Fisher Scientific) within 8 h of blood collection. Gamma-Interferon (IFN-γ) production was measured in duplicate samples by ELISA using a semi-quantitative commercial diagnostic kit according to kit instructions (Bovigam, Thermo-Fisher Scientific). Absorbance values at OD450 nm were converted to OD units using the formula, OD450 × 1000. A sample was positive when the OD450 of the PPDb stimulated sample exceeded 100 OD units (PPDb > 100), was greater than the nil un-stimulated sample by 50 OD units (PPDb - Nil sample > 50) and was greater than the PPDa stimulated sample (PPDb - PPDa > 80). These interpretation criteria are used in the Irish national bTB eradication programme [10].

Sixteen animals were chosen as either positive (n = 8) or negative (n = 8) on the single intradermal comparative tuberculin test (SICTT). Among these, four animals in each group were also positive on the IFN-γ ELISA (Supplementary Table S1). In follow-up of animals after the study was completed, all but two of the animals had been euthanised. Post-mortem examination at the abattoir revealed that three animals (2 x SICTT positive / IFN-γ positive, 1 x SICTT positive / IFN-γ negative) disclosed caseous lesions in the retropharyngeal lymph nodes, consistent with M. bovis infection. Ethical approval for their use in this study was provided by the UCD ethics committee under licence number AREC-E-22-33-Meade. Genomic DNA was subsequently extracted from frozen whole blood representing 16 Friesian cows as shown in Fig. 1. A full list of test results for each animal is shown in Supplementary Table S1. DNA was extracted using the QIAamp Blood Mini Kit (Qiagen) and the purity of the extracted DNA was assessed using the 260/280 ratio.

Bisulfite treatment and WGBS library preparation

The EZ DNA Methylation Gold Kit (Cat No. D5005) Zymo Research was used to perform the bisulfite treatment of the unmethylated cytosines present in extracted DNA samples. Cluster generation and high-throughput sequencing of post-bisulfite adaptor tagging (PBAT) libraries were performed using an Illumina® NovaSeq™ 6000 sequencing system on an S4 flow cell with paired-end 150 bp reads using the Illumina S4 reagent kit v1.5, by Novogene (UK) Company Ltd., Cambridge. Uniquely indexed samples were pooled for partial lane sequencing with demultiplexing and initially run to generate one Gigabase (Gb) of sequencing data (approximately 3.3 M PE150 reads) to assess the cytosine conversion ratio. Subsequently, samples with unique indices were pooled for further partial lane sequencing with demultiplexing until approximately 90 Gb (300 M PE150 reads) of data per sample was generated.

Data quality control, estimation of methylation levels and identification of differential methylated regions (DMRs)

FastQC (fastqc_v0.11.5) was used to perform basic statistics assessing the quality of the raw reads. Read sequences produced by the Illumina pipeline in FASTQ format were subsequently pre-processed through Trimmomatic (Trimmomatic-0.36) software using the parameter (SLIDINGWINDOW: 4:15; LEADING:3, TRAILING:3; ILLUMINACLIP: adapter.fa: 2:30:10; MINLEN:36). Reads that passed all the filtering steps (clean reads) were used for all subsequent analyses. FastQC was used to perform basic quality statistics on the clean data reads. Prior to the analysis, reference data for the bovine, were assembled which included the reference sequence fasta file, the annotation file in gtf format, the GO annotation file, the description file and the gene region file in bed format. For the bed files, repeat regions were predicted using RepeatMasker, followed by using a CGI track from a genome using cpgIslandExt. Bismark software (version 0.16.3; [52]) was used to perform alignments of bisulfite-treated reads to the bovine reference genome (-X 700 --dovetail). The reference genome was firstly transformed into a bisulfite-converted version (C-to-T and G-to-A converted) and then indexed using bowtie2 software [53]. Sequence reads were also transformed into bisulfite-converted versions (C-to-T and G-to-A converted) before they were aligned to similarly converted versions of the bovine genome in a directional manner. Sequence reads that produced a unique best alignment from the two alignment processes (top and bottom strand) were then compared to the normal genomic sequence and the methylation state of all cytosine positions in the read was inferred. The same reads that aligned to the same regions of genome were regarded as duplicated. The sequencing depth and coverage were summarized using de-duplicated reads. Results from methylation extractor (bismark_methylation_extractor, -- no_overlap) were transformed into bigWig format for visualization using the IGV browser. The sodium bisulfite non-conversion rate was calculated as the percentage of cytosine sequenced at cytosine reference positions in the lambda genome. To identify the methylation site, the sum Mc of methylated counts was modelled as a binomial (Bin) random variable with methylation rate. In order to calculate the methylation level, the sequence was divided into multiple bins (bin size of 10 kb) and the sum of methylated and unmethylated read counts in each window calculated. Methylation level (ML) for each window or C site represents the fraction of methylated Cs and was defined as: ML(C) = reads (mC)/reads (mC) + reads (C). Differentially methylated regions (DMRs) were identified using DSS software [54,55,56] and the core of DSS is a new dispersion shrinkage method for estimating the dispersion parameter from Gamma-Poisson or Beta-Binomial distributions. DSS possesses three characteristics to detect DMRs: firstly, spatial correlation and proper utilization of the information from neighbouring cytosine sites improves the estimation of ML at each C site, and hence improve DMR detection; secondly, the read depth of the C sites provides information on precision used to improve statistical tests for DMR detection; and thirdly, the variance among biological replicates provides information necessary for a valid statistical test to detect DMRs. According to the distribution of DMRs through the genome, genes related to DMRs were defined as genes whose gene body region (from TSS to TES) or promoter region (upstream 2 kb from the TSS) have an overlap with the DMRs.

Gene ontology and pathway analysis

Gene Ontology (GO) enrichment analysis of genes related to DMRs was implemented by the GOseq R package [55], in which gene length bias was corrected and GO terms with corrected FDR P-value < 0.05 were considered significantly enriched by DMR-related genes. KOBAS software [57] was used to test the statistical enrichment of DMR-related genes using KEGG [58]. Bioinformatic analysis was performed by the Novogene bioinformatics team and the described methods herein represent a common data analysis pipeline described similarly in other studies and which were provided by Novogene (www.novogene.com).

Data availability

The data sets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found at: www.ebi.ac.uk/biostudies/arrayexpress/studies/E-MTAB-13487.

References

  1. Fitzgerald SD, Kaneene JB. Wildlife reservoirs of bovine tuberculosis worldwide: hosts, pathology, surveillance, and control. Vet Pathol. 2013;50(3):488–99.

    Article  CAS  PubMed  Google Scholar 

  2. Allen AR, Skuce RA, Byrne AW. Bovine tuberculosis in Britain and Ireland - A Perfect Storm? The confluence of potential ecological and epidemiological impediments to Controlling a chronic infectious disease. Front Vet Sci. 2018;5:109.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  3. Reviriego Gordejo FJ, Vermeersch JP. Towards eradication of bovine tuberculosis in the European Union. Vet Microbiol. 2006;112(2–4):101–9.

    Article  CAS  PubMed  Google Scholar 

  4. Pollock JM, McNair J, Bassett H, Cassidy JP, Costello E, Aggerbeck H, Rosenkrands I, Andersen P. Specific delayed-type hypersensitivity responses to ESAT-6 identify tuberculosis-infected cattle. J Clin Microbiol. 2003;41(5):1856–60.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Gormley E, Doyle MB, Fitzsimons T, McGill K, Collins JD. Diagnosis of Mycobacterium bovis infection in cattle by use of the gamma-interferon (Bovigam) assay. Vet Microbiol. 2006;112(2–4):171–9.

    Article  CAS  PubMed  Google Scholar 

  6. de la Rua-Domenech R, Goodchild AT, Vordermeier HM, Hewinson RG, Christiansen KH, Clifton-Hadley RS. Ante mortem diagnosis of tuberculosis in cattle: a review of the tuberculin tests, gamma-interferon assay and other ancillary diagnostic techniques. Res Vet Sci. 2006;81(2):190–210.

    Article  PubMed  Google Scholar 

  7. O’Hagan MJ, Courcier EA, Drewe JA, Gordon AW, McNair J, Abernethy DA. Risk factors for visible lesions or positive laboratory tests in bovine tuberculosis reactor cattle in Northern Ireland. Prev Vet Med. 2015;120(3–4):283–90.

    Article  PubMed  Google Scholar 

  8. Clegg TA, Good M, Hayes M, Duignan A, McGrath G, More SJ. Trends and predictors of large tuberculosis episodes in cattle herds in Ireland. Front Vet Sci. 2018;5:86.

    Article  PubMed  PubMed Central  Google Scholar 

  9. Clegg TA, Good M, More SJ. Future risk of bovine tuberculosis recurrence among higher risk herds in Ireland. Prev Vet Med. 2015;118(1):71–9.

    Article  CAS  PubMed  Google Scholar 

  10. Madden JM, O’Donovan J, Casey-Bryars M, Sweeney J, Messam LL, McAloon CG, More SJ, Kenny K, Ryan E, Gormley E. The impact of changing the cut-off threshold of the interferon-gamma (IFN-gamma) assay for diagnosing bovine tuberculosis in Ireland. Prev Vet Med. 2024;224:106129.

    Article  PubMed  Google Scholar 

  11. Clegg TA, Doyle M, Ryan E, More SJ, Gormley E. Characteristics of Mycobacterium bovis infected herds tested with the interferon-gamma assay. Prev Vet Med. 2019;168:52–9.

    Article  CAS  PubMed  Google Scholar 

  12. Lahuerta-Marin A, Gallagher M, McBride S, Skuce R, Menzies F, McNair J, McDowell SW, Byrne AW. Should they stay, or should they go? Relative future risk of bovine tuberculosis for interferon-gamma test-positive cattle left on farms. Vet Res. 2015;46(1):90.

    Article  PubMed  PubMed Central  Google Scholar 

  13. McLoughlin KE, Nalpas NC, Rue-Albrecht K, Browne JA, Magee DA, Killick KE, Park SD, Hokamp K, Meade KG, O’Farrelly C, et al. RNA-seq transcriptional profiling of peripheral blood leukocytes from cattle infected with Mycobacterium bovis. Front Immunol. 2014;5:396.

    Article  PubMed  PubMed Central  Google Scholar 

  14. McLoughlin KE, Correia CN, Browne JA, Magee DA, Nalpas NC, Rue-Albrecht K, Whelan AO, Villarreal-Ramos B, Vordermeier HM, Gormley E, et al. RNA-Seq transcriptome analysis of peripheral blood from cattle infected with Mycobacterium bovis across an experimental time course. Front Vet Sci. 2021;8:662002.

    Article  PubMed  PubMed Central  Google Scholar 

  15. Bhat SA, Elnaggar M, Hall TJ, McHugo GP, Reid C, MacHugh DE, Meade KG. Preferential differential gene expression within the WC1.1(+) gammadelta T cell compartment in cattle naturally infected with Mycobacterium bovis. Front Immunol. 2023;14:1265038.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  16. Hall TJ, McHugo GP, Mullen MP, Ward JA, Killick KE, Browne JA, Gordon SV, MacHugh DE. Integrative and comparative genomic analyses of mammalian macrophage responses to intracellular mycobacterial pathogens. Tuberculosis (Edinb) 2023:102453.

  17. Hall TJ, Vernimmen D, Browne JA, Mullen MP, Gordon SV, MacHugh DE. O’Doherty AM: alveolar macrophage chromatin is modified to Orchestrate Host Response to Mycobacterium bovis infection. Front Genet. 2019;10:1386.

    Article  CAS  PubMed  Google Scholar 

  18. O’Doherty AM, Rue-Albrecht KC, Magee DA, Ahting S, Irwin RE, Hall TJ, Browne JA, Nalpas NC, Walsh CP, Gordon SV, et al. The bovine alveolar macrophage DNA methylome is resilient to infection with Mycobacterium bovis. Sci Rep. 2019;9(1):1510.

    Article  PubMed  PubMed Central  Google Scholar 

  19. Doherty R, C OF, Meade KG. Comparative epigenetics: relevance to the regulation of production and health traits in cattle. Anim Genet. 2014;45(Suppl 1):3–14.

    Article  PubMed  Google Scholar 

  20. Zhang J, Sheng H, Hu C, Li F, Cai B, Ma Y, Wang Y, Ma Y. Effects of DNA methylation on Gene expression and phenotypic traits in cattle: a review. Int J Mol Sci. 2023;24(15):11882.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. de Araujo-Souza PS, Hanschke SC, Viola JP. Epigenetic control of interferon-gamma expression in CD8 T cells. J Immunol Res 2015, 2015:849573.

  22. Gauba K, Gupta S, Shekhawat J, Sharma P, Yadav D, Banerjee M. Immunomodulation by epigenome alterations in Mycobacterium tuberculosis infection. Tuberculosis (Edinb). 2021;128:102077.

    Article  CAS  PubMed  Google Scholar 

  23. Shell SS, Prestwich EG, Baek SH, Shah RR, Sassetti CM, Dedon PC, Fortune SM. DNA methylation impacts gene expression and ensures hypoxic survival of Mycobacterium tuberculosis. PLoS Pathog. 2013;9(7):e1003419.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Villicana S, Bell JT. Genetic impacts on DNA methylation: research findings and future perspectives. Genome Biol. 2021;22(1):127.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Yan J, Zierath JR, Barres R. Evidence for non-CpG methylation in mammals. Exp Cell Res. 2011;317(18):2555–61.

    Article  CAS  PubMed  Google Scholar 

  26. He XJ, Chen T, Zhu JK. Regulation and function of DNA methylation in plants and animals. Cell Res. 2011;21(3):442–65.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. DiNardo AR, Rajapakshe K, Nishiguchi T, Grimm SL, Mtetwa G, Dlamini Q, Kahari J, Mahapatra S, Kay A, Maphalala G, et al. DNA hypermethylation during tuberculosis dampens host immune responsiveness. J Clin Invest. 2020;130(6):3113–23.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. Cooper AM, Magram J, Ferrante J, Orme IM. Interleukin 12 (IL-12) is crucial to the development of protective immunity in mice intravenously infected with mycobacterium tuberculosis. J Exp Med. 1997;186(1):39–45.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Ibeagha-Awemu EM, Bissonnette N, Bhattarai S, Wang M, Dudemaine PL, McKay S, Zhao X. Whole genome methylation analysis reveals role of DNA methylation in cow’s Ileal and Ileal Lymph node responses to Mycobacterium avium subsp. paratuberculosis infection. Front Genet. 2021;12:797490.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  30. Wang M, Feng S, Ma G, Miao Y, Zuo B, Ruan J, Zhao S, Wang H, Du X, Liu X. Whole-genome methylation analysis reveals epigenetic variation in Cloned and Donor pigs. Front Genet. 2020;11:23.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Heberle E, Bardet AF. Sensitivity of transcription factors to DNA methylation. Essays Biochem. 2019;63(6):727–41.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  32. Li S, Zhang J, Huang S, He X. Genome-wide analysis reveals that exon methylation facilitates its selective usage in the human transcriptome. Brief Bioinform. 2018;19(5):754–64.

    Article  CAS  PubMed  Google Scholar 

  33. He Y, Ecker JR. Non-CG methylation in the Human Genome. Annu Rev Genomics Hum Genet. 2015;16:55–77.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Feng H, Wu H. Differential methylation analysis for bisulfite sequencing using DSS. Quant Biol. 2019;7(4):327–34.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Corso-Diaz X, Gentry J, Rebernick R, Jaeger C, Brooks MJ, van Asten F, Kooragayala K, Gieser L, Nellissery J, Covian R, et al. Genome-wide profiling identifies DNA methylation signatures of aging in Rod Photoreceptors Associated with alterations in Energy Metabolism. Cell Rep. 2020;31(3):107525.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Zhang J, Han B, Zheng W, Lin S, Li H, Gao Y, Sun D. Genome-wide DNA methylation Profile in Jejunum reveals the potential genes Associated with paratuberculosis in dairy cattle. Front Genet. 2021;12:735147.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  37. Omiya R, Tsushima F, Narazaki H, Sakoda Y, Kuramasu A, Kim Y, Xu H, Tamura H, Zhu G, Chen L, et al. Leucocyte-associated immunoglobulin-like receptor-1 is an inhibitory regulator of contact hypersensitivity. Immunology. 2009;128(4):543–55.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Mistry R, Cliff JM, Clayton CL, Beyers N, Mohamed YS, Wilson PA, Dockrell HM, Wallace DM, van Helden PD, Duncan K, et al. Gene-expression patterns in whole blood identify subjects at risk for recurrent tuberculosis. J Infect Dis. 2007;195(3):357–65.

    Article  CAS  PubMed  Google Scholar 

  39. Blischak JD, Tailleux L, Mitrano A, Barreiro LB, Gilad Y. Mycobacterial infection induces a specific human innate immune response. Sci Rep. 2015;5:16882.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  40. Harding CV, Boom WH. Regulation of antigen presentation by Mycobacterium tuberculosis: a role for toll-like receptors. Nat Rev Microbiol. 2010;8(4):296–307.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Tenorio de Menezes YK, Eto C, de Oliveira J, Larson EC, Mendes D, Dias GBM, Delgobo M, Gubernat AK, Gleim JL, Munari EL, et al. The endogenous retinoic acid receptor pathway is exploited by Mycobacterium tuberculosis during infection, both in Vitro and in vivo. J Immunol. 2023;211(4):601–11.

    Article  CAS  PubMed  Google Scholar 

  42. Malik ZA, Denning GM, Kusner DJ. Inhibition of ca(2+) signaling by Mycobacterium tuberculosis is associated with reduced phagosome-lysosome fusion and increased survival within human macrophages. J Exp Med. 2000;191(2):287–302.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Pasquinelli V, Rovetta AI, Alvarez IB, Jurado JO, Musella RM, Palmero DJ, Malbran A, Samten B, Barnes PF, Garcia VE. Phosphorylation of mitogen-activated protein kinases contributes to interferon gamma production in response to Mycobacterium tuberculosis. J Infect Dis. 2013;207(2):340–50.

    Article  CAS  PubMed  Google Scholar 

  44. Howard NC, Khader SA. Immunometabolism during Mycobacterium tuberculosis infection. Trends Microbiol. 2020;28(10):832–50.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  45. Catena V, Fanciulli M. Deptor: not only a mTOR inhibitor. J Exp Clin Cancer Res. 2017;36(1):12.

    Article  PubMed  PubMed Central  Google Scholar 

  46. Shen H, Chen ZW. The crucial roles of Th17-related cytokines/signal pathways in M. Tuberculosis infection. Cell Mol Immunol. 2018;15(3):216–25.

    Article  CAS  PubMed  Google Scholar 

  47. Waters WR, Maggioli MF, Palmer MV, Thacker TC, McGill JL, Vordermeier HM, Berney-Meyer L, Jacobs WR Jr., Larsen MH. Interleukin-17A as a biomarker for bovine tuberculosis. Clin Vaccine Immunol. 2016;23(2):168–80.

    Article  PubMed  PubMed Central  Google Scholar 

  48. Dong H, Lv Y, Zhao D, Barrow P, Zhou X. Defensins: The Case for Their Use against Mycobacterial Infections. J Immunol Res 2016, 2016:7515687.

  49. Choi SW, Kim S, Park HT, Park HE, Choi JS, Yoo HS. MicroRNA profiling in bovine serum according to the stage of Mycobacterium avium subsp. paratuberculosis infection. PLoS ONE. 2021;16(11):e0259539.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Blanco FC, Gravisaco MJ, Bigi MM, Garcia EA, Marquez C, McNeil M, Jackson M, Bigi F. Identifying bacterial and Host Factors Involved in the Interaction of Mycobacterium bovis with the bovine Innate Immune cells. Front Immunol. 2021;12:674643.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  51. Seiler Vellame D, Castanho I, Dahir A, Mill J, Hannon E. Characterizing the properties of bisulfite sequencing data: maximizing power and sensitivity to identify between-group differences in DNA methylation. BMC Genomics. 2021;22(1):446.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Krueger F, Andrews SR. Bismark: a flexible aligner and methylation caller for Bisulfite-Seq applications. Bioinformatics. 2011;27(11):1571–2.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012;9(4):357–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Feng H, Conneely KN, Wu H. A bayesian hierarchical model to detect differentially methylated loci from single nucleotide resolution sequencing data. Nucleic Acids Res. 2014;42(8):e69.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  55. Young MD, Wakefield MJ, Smyth GK, Oshlack A. Gene ontology analysis for RNA-seq: accounting for selection bias. Genome Biol. 2010;11(2):R14.

    Article  PubMed  PubMed Central  Google Scholar 

  56. Park Y, Wu H. Differential methylation analysis for BS-seq data under general experimental design. Bioinformatics. 2016;32(10):1446–53.

    Article  CAS  PubMed  Google Scholar 

  57. Mao X, Cai T, Olyarchuk JG, Wei L. Automated genome annotation and pathway identification using the KEGG Orthology (KO) as a controlled vocabulary. Bioinformatics. 2005;21(19):3787–93.

    Article  CAS  PubMed  Google Scholar 

  58. Kanehisa M, Araki M, Goto S, Hattori M, Hirakawa M, Itoh M, Katayama T, Kawashima S, Okuda S, Tokimatsu T, et al. KEGG for linking genomes to life and the environment. Nucleic Acids Res. 2008;36(Database issue):D480–484.

    CAS  PubMed  Google Scholar 

Download references

Acknowledgements

The authors would like to express their gratitude to Dr John Browne, UCD for assistance with DNA extraction.

Funding

This project was funded by Science Foundation Ireland (SFI) Grant 17/CDA/4717 and AP was funded via the Higher Education Authority on behalf of the Department of Further and Higher Education, Research Innovation and Science, and the Shared Island Unit at the Department of the Taoiseach as part of the North-South Research Programme on the Vit-TB project.

Author information

Authors and Affiliations

  1. UCD School of Agriculture and Food Science, University College Dublin, Belfield, Dublin, D04 V1W8, Ireland

    Sajad A. Bhat, Alia Parveen & Kieran G. Meade

  2. UCD School of Veterinary Medicine, University College Dublin, Belfield, Dublin, D04 V1W8, Ireland

    Eamonn Gormley

  3. UCD Conway Institute of Biomolecular and Biomedical Research, University College Dublin, Belfield, Dublin, D04 V1W8, Ireland

    Kieran G. Meade

  4. UCD Institute of Food and Health, University College Dublin, Belfield, Dublin, C15 PW93, Ireland

    Kieran G. Meade

Authors
  1. Sajad A. Bhat
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Alia Parveen
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Eamonn Gormley
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Kieran G. Meade
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

Conceived the study and provided samples: KGM and EG. Performed experiments and interpreted data: KGM, SB and AP. Data analysis: AP and KGM. Wrote the manuscript: KGM, SB, EG and AP. All authors reviewed and approved the final manuscript.

Corresponding author

Correspondence to Kieran G. Meade.

Ethics declarations

Ethical approval

Ethical approval was provided by the UCD ethics committee under licence number AREC-E-22-33-Meade and reporting herein is presented in line with ARRIVE 2.0 guidelines.

Consent for publication

Not applicable.

Conflict of interest

All authors declare that they have no competing interests, and that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Additional information

Publisher’s Note

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

Electronic supplementary material

Below is the link to the electronic supplementary material.

Supplementary Material 1

Supplementary Material 2 - Supplementary Table S1:

Single intradermal comparative tuberculin test (SICTT), IFN-γ release assay results, bTB classification and sampling information on all cattle sampled

Supplementary Material 3 - Supplementary Table S2:

Bisulphite conversion data, genome mapping statistics, % methylation by sample and by genomic feature for 16 WGBS libraries from M. bovis-infected and non-infected control cattle

Supplementary Material 4 (A and B)- Supplementary Figure S1 (A and B):

Methylation profiles for all 16 WGBS samples across genomic features (A) and genes (B). For group information on sample ID, refer to Table S1

Supplementary Material 5 - Supplementary Table S3:

Differentially expressed regions (DMRs) across all genomic features from M. bovis-infected and non-infected control animals is shown in Tab 1. Subsequent tabs are subsets of this data including DMRs in exons and promoters with >15% differential methylation

Supplementary Material 6 - Supplementary Table S4:

Differentially expressed promoter region genes (DMPGs) in various genomic features from M. bovis-infected and non-infected control animals. Genes identified with CH, CHH or CHG methylation types are listed with Ensembl ID and Gene Symbols. The overlap tab lists genes with multiple methylation types

Supplementary Material 7 - Supplementary Table S5:

GO enrichment analysis showing enriched GO terms (biological processes, cellular component and molecular function) for both DMR and DPMG genes from M. bovis-infected and non-infected control animals

Supplementary Material 8 - Supplementary Table S6:

KEGG analysis showing enriched canonical pathways for both DMR and DPMG genes from M. bovis-infected and non-infected control animals

Supplementary Material 9 - Supplementary Table S7:

Differentially expressed regions (DMRs) across all genomic features from Group 2 (SICTT-/IFNG+) samples relative to Group 3 (SICTT+/IFNG+) samples is shown in Tab 1. Subsequent tabs are subsets of this data including DMRs in exons and promoters with >15% differential methylation

Supplementary Material 10 - Supplementary Table S8:

Differentially expressed promoter region genes (DMPGs) in various genomic features from Group 2 (SICTT-/IFNG+) samples relative to Group 3 (SICTT+/IFNG+) samples. Genes identified with CH, CHH or CHG methylation types are listed with Ensembl ID and Gene Symbols. The overlap tab lists genes with multiple methylation types

Supplementary Material 11 - Supplementary Table S9:

GO enrichment analysis showing enriched GO terms (biological processes, cellular component and molecular function) for DMRs and DPMGs for Group 2 (SICTT-/IFNG+) samples relative to Group 3 (SICTT+/IFNG+) samples (as shown in Figure 1)

Supplementary Material 12 - Supplementary Table S10:

KEGG analysis showing enriched canonical pathways for DMRs and DPMGs for Group 2 (SICTT-/IFNG+) samples relative to Group 3 (SICTT+/IFNG+) samples (as shown in Figure 1)

Supplementary Material 13 - Figure S2:

Differential methylation profiles for all 16 WGBS samples, divided according to experimental groups (as shown in Figure 1). Methylation plots are shown for combined levels of CG, CHG and CHH methylation types as well as each type individually for pairwise each group comparison

Supplementary Material 14 - Figure S3:

Ratio of hypermethylation to hypomethylation profiles for all 16 WGBS samples, divided according to experimental groups (as shown in Figure 1) across each genomic feature. Histograms show ratio for CG, CHG and CHH methylation types for each pairwise group comparison

Supplementary Material 15 - Figure S4:

Numbers of genes identified as differentially methylated (CG, CHG or CHH) in DMR regions and in the promoter region of genes for all 16 WGBS samples, divided according to experimental groups (as shown in Figure 1). The degree of overlap shows the numbers of genes with multiple forms of methylation

Supplementary Material 16 - Figure S5:

GO plot enrichment of biological processes, cellular component and molecular function for all 16 WGBS samples, divided according to experimental groups (as shown in Figure 1). Significantly enriched GO categories are identified with an asterix

Supplementary Material 17 - Figure S6:

Scatterplot showing significantly enriched pathways from DMRs and DPMGs for all 16 WGBS samples, divided according to experimental groups (as shown in Figure 1) identified using KEGG. Corrected P values (Q value) is indicated by the colour and the numbers of genes represented in the pathway is indicated by the size of the circle. The Rich factor identifies the ratio of differentially expressed gene numbers annotated in this pathway term relative to all gene numbers annotated in this pathway term

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. 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 in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Bhat, S.A., Parveen, A., Gormley, E. et al. Extensive differential DNA methylation between tuberculosis skin test positive and skin test negative cattle. BMC Genomics 25, 762 (2024). https://doi.org/10.1186/s12864-024-10574-x

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12864-024-10574-x

Keywords

BMC Genomics

ISSN: 1471-2164

Contact us