Accurate and robust classification of Mycobacterium bovis-infected cattle using peripheral blood RNA-seq data

Accurate and robust classification of Mycobacterium bovis-infected cattle using peripheral blood RNA-seq data | bioRxiv /* */ /* */ <!-- <!-- /*! * yepnope1.5.4 * (c) WTFPL, GPLv2 */ (function(a,b,c){function d(a){return"[object Function]"==o.call(a)}function e(a){return"string"==typeof a}function f(){}function g(a){return!a||"loaded"==a||"complete"==a||"uninitialized"==a}function h(){var a=p.shift();q=1,a?a.t?m(function(){("c"==a.t?B.injectCss:B.injectJs)(a.s,0,a.a,a.x,a.e,1)},0):(a(),h()):q=0}function i(a,c,d,e,f,i,j){function k(b){if(!o&&g(l.readyState)&&(u.r=o=1,!q&&h(),l.onload=l.onreadystatechange=null,b)){"img"!=a&&m(function(){t.removeChild(l)},50);for(var d in y[c])y[c].hasOwnProperty(d)&&y[c][d].onload()}}var j=j||B.errorTimeout,l=b.createElement(a),o=0,r=0,u={t:d,s:c,e:f,a:i,x:j};1===y[c]&&(r=1,y[c]=[]),"object"==a?l.data=c:(l.src=c,l.type=a),l.width=l.height="0",l.onerror=l.onload=l.onreadystatechange=function(){k.call(this,r)},p.splice(e,0,u),"img"!=a&&(r||2===y[c]?(t.insertBefore(l,s?null:n),m(k,j)):y[c].push(l))}function j(a,b,c,d,f){return q=0,b=b||"j",e(a)?i("c"==b?v:u,a,b,this.i++,c,d,f):(p.splice(this.i++,0,a),1==p.length&&h()),this}function k(){var a=B;return a.loader={load:j,i:0},a}var l=b.documentElement,m=a.setTimeout,n=b.getElementsByTagName("script")[0],o={}.toString,p=[],q=0,r="MozAppearance"in l.style,s=r&&!!b.createRange().compareNode,t=s?l:n.parentNode,l=a.opera&&"[object Opera]"==o.call(a.opera),l=!!b.attachEvent&&!l,u=r?"object":l?"script":"img",v=l?"script":u,w=Array.isArray||function(a){return"[object Array]"==o.call(a)},x=[],y={},z={timeout:function(a,b){return b.length&&(a.timeout=b[0]),a}},A,B;B=function(a){function b(a){var a=a.split("!"),b=x.length,c=a.pop(),d=a.length,c={url:c,origUrl:c,prefixes:a},e,f,g;for(f=0;f<d;f++)g=a[f].split("="),(e=z[g.shift()])&&(c=e(c,g));for(f=0;f<b;f++)c=x[f](c);return c}function g(a,e,f,g,h){var i=b(a),j=i.autoCallback;i.url.split(".").pop().split("?").shift(),i.bypass||(e&&(e=d(e)?e:e[a]||e[g]||e[a.split("/").pop().split("?")[0]]),i.instead?i.instead(a,e,f,g,h):(y[i.url]?i.noexec=!0:y[i.url]=1,f.load(i.url,i.forceCSS||!i.forceJS&&"css"==i.url.split(".").pop().split("?").shift()?"c":c,i.noexec,i.attrs,i.timeout),(d(e)||d(j))&&f.load(function(){k(),e&&e(i.origUrl,h,g),j&&j(i.origUrl,h,g),y[i.url]=2})))}function h(a,b){function c(a,c){if(a){if(e(a))c||(j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}),g(a,j,b,0,h);else if(Object(a)===a)for(n in m=function(){var b=0,c;for(c in a)a.hasOwnProperty(c)&&b++;return b}(),a)a.hasOwnProperty(n)&&(!c&&!--m&&(d(j)?j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}:j[n]=function(a){return function(){var b=[].slice.call(arguments);a&&a.apply(this,b),l()}}(k[n])),g(a[n],j,b,n,h))}else!c&&l()}var h=!!a.test,i=a.load||a.both,j=a.callback||f,k=j,l=a.complete||f,m,n;c(h?a.yep:a.nope,!!i),i&&c(i)}var i,j,l=this.yepnope.loader;if(e(a))g(a,0,l,0);else if(w(a))for(i=0;i (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0];var j=d.createElement(s);var dl=l!='dataLayer'?'&l='+l:'';j.src='//www.googletagmanager.com/gtm.js?id='+i+dl;j.type='text/javascript';j.async=true;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-M677548'); Skip to main content Home About Submit ALERTS / RSS Search for this keyword Advanced Search New Results Accurate and robust classification of Mycobacterium bovis -infected cattle using peripheral blood RNA-seq data View ORCID Profile John F. O’Grady , Adriana Ivich , View ORCID Profile Gillian P. McHugo , View ORCID Profile Adnan Khan , View ORCID Profile Thomas J. Hall , Sarah L. F. O’Donnell , View ORCID Profile Carolina N. Correia , View ORCID Profile John A. Browne , View ORCID Profile Valentina Riggio , View ORCID Profile James G. D. Prendergast , View ORCID Profile Emily L. Clark , View ORCID Profile Hubert Pausch , View ORCID Profile Kieran G. Meade , Isobel C. Gormley , View ORCID Profile Eamonn Gormley , View ORCID Profile Stephen V. Gordon , View ORCID Profile Casey S. Greene , View ORCID Profile David E. MacHugh doi: https://doi.org/10.1101/2025.09.01.673540 John F. O’Grady 1 UCD School of Agriculture and Food Science, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for John F. O’Grady Adriana Ivich 2 Department of Biomedical Informatics, University of Colorado Anschutz Medical Campus , Aurora, CO, USA Find this author on Google Scholar Find this author on PubMed Search for this author on this site Gillian P. McHugo 1 UCD School of Agriculture and Food Science, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Gillian P. McHugo Adnan Khan 1 UCD School of Agriculture and Food Science, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Adnan Khan Thomas J. Hall 1 UCD School of Agriculture and Food Science, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Thomas J. Hall Sarah L. F. O’Donnell 1 UCD School of Agriculture and Food Science, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site Carolina N. Correia 1 UCD School of Agriculture and Food Science, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Carolina N. Correia John A. Browne 1 UCD School of Agriculture and Food Science, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for John A. Browne Valentina Riggio 3 The Roslin Institute and Royal (Dick) School of Veterinary Studies, University of Edinburgh , Midlothian, EH25 9RG, UK 4 Centre for Tropical Livestock Genetics and Health (CTLGH), Roslin Institute, University of Edinburgh , Easter Bush Campus, EH25 9RG, UK Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Valentina Riggio James G. D. Prendergast 3 The Roslin Institute and Royal (Dick) School of Veterinary Studies, University of Edinburgh , Midlothian, EH25 9RG, UK 4 Centre for Tropical Livestock Genetics and Health (CTLGH), Roslin Institute, University of Edinburgh , Easter Bush Campus, EH25 9RG, UK Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for James G. D. Prendergast Emily L. Clark 3 The Roslin Institute and Royal (Dick) School of Veterinary Studies, University of Edinburgh , Midlothian, EH25 9RG, UK 4 Centre for Tropical Livestock Genetics and Health (CTLGH), Roslin Institute, University of Edinburgh , Easter Bush Campus, EH25 9RG, UK Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Emily L. Clark Hubert Pausch 5 Animal Genomics, ETH Zurich, Universitaetstrasse 2, 8006, Zurich, Switzerland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Hubert Pausch Kieran G. Meade 1 UCD School of Agriculture and Food Science, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland 6 UCD Conway Institute of Biomolecular and Biomedical Research, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland 7 UCD One Health Centre, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Kieran G. Meade Isobel C. Gormley 8 UCD School of Mathematics and Statistics, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site Eamonn Gormley 7 UCD One Health Centre, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland 9 UCD School of Veterinary Medicine, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Eamonn Gormley Stephen V. Gordon 6 UCD Conway Institute of Biomolecular and Biomedical Research, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland 7 UCD One Health Centre, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland 9 UCD School of Veterinary Medicine, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Stephen V. Gordon Casey S. Greene 2 Department of Biomedical Informatics, University of Colorado Anschutz Medical Campus , Aurora, CO, USA Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Casey S. Greene David E. MacHugh 1 UCD School of Agriculture and Food Science, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland 6 UCD Conway Institute of Biomolecular and Biomedical Research, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland 7 UCD One Health Centre, University College Dublin , Belfield, Dublin, D04 V1W8, Ireland Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for David E. MacHugh For correspondence: david.machugh{at}ucd.ie Abstract Full Text Info/History Metrics Supplementary material Data/Code Preview PDF Abstract The zoonotic bacterium, Mycobacterium bovis , causes bovine tuberculosis (bTB) and is closely related to Mycobacterium tuberculosis , the primary cause of human tuberculosis (hTB). Bovine TB remains recalcitrant to eradication in endemic countries where current diagnostics fail to identify all infected animals. While blood-based RNA biomarkers identified through machine learning have shown accurate discrimination of hTB-positive and hTB-negative individuals, similar approaches have not been explored for bTB. Here, we use RNA-seq and machine learning to investigate the utility of peripheral blood mRNA as a host-response biomarker for bTB using data from Ireland, the UK and the US. We identify a 30-gene signature and a 273-gene elastic net classifier that differentiate bTB-positive from bTB-negative cattle, achieving area under the curve (AUC) values of 0.986/0.900 for the former and 0.968/0.938 for the latter in training and testing, respectively. These two classifiers produced high sensitivity and specificity values (≥ 0.853 for both metrics) in the testing set. Additionally, we show that they robustly distinguish bTB+ animals from those infected with other bacterial or viral pathogens (AUC ≥ 0.819). These RNA-based classifiers accurately diagnose bTB and differentiate bTB from other diseases, representing a promising tool for augmenting current diagnostics to advance bTB eradication efforts in endemic regions. Introduction Mycobacterium bovis causes bovine tuberculosis (bTB), a damaging and endemic disease principally associated with cattle that is conservatively estimated to cost the global agricultural industry $3 billion USD annually ( Waters et al . 2012 ). Moreover, M. bovis possesses the capacity to infect humans, causing zoonotic tuberculosis (zTB), a disease that is most prevalent in the Global South and due to direct aerosol transmission of the bacillus to animal keepers, or by consumption of contaminated animal by-products ( Cosivi et al . 1998 ; Michel et al . 2015 ; Olea-Popelka et al . 2017 ; Sabuz et al . 2025 ). Although variable estimates have been reported for the proportion of human tuberculosis (hTB) cases attributed to M. bovis (i.e., zTB), including as high as 28% in Mexico ( Portillo-Gómez & Sosa-Iglesias 2011 ), the most recent estimates from the World Health Organization (WHO) attributed 140,000 of newly diagnosed active hTB cases (1.4% of 10.0 million) and 11,400 hTB deaths (8.1% of 1.2 million) to zTB ( World Health Organization 2020 ). In some industrialised countries such as Ireland, England, and Wales, bTB represents a persistent economic and social problem to the agricultural industry, inflicting combined financial costs in excess of €200 million per year ( Bennett 2017 ; More 2024 ), and causing significant distress to herdowners ( Clarke et al . 2022 ). The high economic cost associated with bTB disease management in these countries stems from state-sponsored bTB eradication programmes. These initiatives, which commenced in the mid-20 th century, rely on compulsory skin testing of all herds to identify infected animals that are then removed, slaughtered, and compensated for ( Atkins & Robinson 2013 ; Ryan et al . 2023 ). Despite initial success in curtailing the herd incidence of bTB disease ( Allen et al . 2018 ), in England and Wales, the bTB herd incidence rate has remained consistently above 6% for the past 15 years ( Department for Environment Food & Rural Affairs 2025 ). Similarly, in Ireland, the herd incidence rate fell to a record low of 3.37% in 2015 ( More 2024 ); however, since then, the 12-month rolling bTB herd incidence rate has almost doubled, rising to 6.40% by June 2025 ( Department of Agriculture Food and the Marine 2025 ). The increase in herd incidence and the recalcitrant nature of bTB to broader eradication efforts are multifaceted and underpinned by numerous risk factors, including, but not limited to: stocking density ( Humblet et al . 2009 ); the presence of wild life reservoirs of M. bovis (e.g., the European badger – Meles meles ) ( Griffin et al . 2005 ); intrapopulation genomic variation ( Ring et al . 2019 ); age and sex ( Broughan et al . 2016 ); and the imperfect properties of current bTB diagnostic tests for correctly classifying M. bovis -infected animals (i.e., low sensitivity), which enables the disease to remain within cattle herds ( Hayton et al . 2025 ). There are currently two diagnostic tests used for national bTB eradication programmes: 1) The in vivo field-based single intradermal comparative tuberculin test (SICTT) used in isolation or in conjunction with; 2) an ancillary in vitro enzyme-linked immunosorbent assay-based interferon-γ (IFN-γ) release assay (IGRA) test to increase the sensitivity of the diagnosis by detecting infected animals as early as 14 days post-infection (dpi) ( Buddle et al . 1995 ; Gormley et al . 2006 ). Exploiting data from a range of studies, detailed meta-analyses have yielded low sensitivity estimates of the SICTT under field conditions, ranging from 0.5 (95% posterior credible interval (Crl): γ = 95%, 0.26 ≤ μ ≤ 0.78) based on data from several different countries ( Nunez-Garcia et al . 2018 ), 0.568 (95% Crl: γ = 95%, 0.53 ≤ μ ≤ 0.63) in Ireland ( Clegg et al . 2011 ), to 0.66 (95% confidence interval (CI); 0.52 − 0.8) in Great Britain ( Conlan et al . 2012 ). In contrast, the IGRA test is considered more sensitive than the SICTT, with estimates ranging from 0.63 to 0.88 in Ireland ( Gormley et al . 2006 ; Clegg et al . 2011 ), 0.67 (95% Crl: γ = 95%, 0.49 ≤ μ ≤ 0.82) based on a statistical meta-analysis of data from different countries ( Nunez-Garcia et al . 2018 ), and 0.93 (95% Crl: γ = 95%, 0.90 ≤ μ ≤ 0.96) in Northern Ireland ( Lahuerta-Marin et al . 2018 ). Conversely, the SICTT is considered more specific than the IGRA, expected to misclassify one bTB− animal as bTB+ for every 5000 tests conducted and as such, is selected as the primary population screening test for bTB ( Goodchild et al . 2015 ). Previous work in humans has shown that peripheral blood (PB) transcriptional biomarkers identified using microarray or RNA sequencing (RNA-seq) technologies can discriminate between active pulmonary hTB cases and asymptomatic latent tuberculosis infection (LTBI) patients or other infectious diseases, and also individuals with active hTB from healthy controls with a high degree of accuracy ( Sweeney et al . 2016 ; Tabone et al . 2021 ; Chang et al . 2024 ). In cattle, PB immune responses have been shown to reflect those at the site of infection ( Rhodes et al . 2000 ), and recent work has shown that there is substantial overlap of differentially expressed genes (DEGs) in bovine alveolar macrophages (AMs) experimentally infected with M. bovis and PB RNA-seq data from cattle naturally infected with M. bovis ( Hall et al . 2021 ; O’Grady et al . 2025 ). However, few studies have attempted to characterise PB transcriptional biomarkers capable of discriminating between bTB− and bTB+ cattle with appropriate discovery (i.e., training) and validation (i.e., testing) cohorts. For example, using microarray data, a set of 15 DEGs derived from peripheral blood mononuclear cells (PBMCs) was proposed as being capable of indicating M. bovis infection status ( Meade et al . 2007 ), which failed to replicate in an analysis of peripheral blood leukocytes (PBLs) from cattle naturally infected with M. bovis ( Killick et al . 2011 ). A recent time series study of cattle experimentally infected with M. bovis revealed a set of 17 DEGs exhibiting increased expression that were consistently differentially expressed (DE) in PB at +1 week post-infection (wpi), +2 wpi, +6 wpi, +10 wpi, and +12 wpi, which showed striking similarity to PBL RNA-seq and microarray results in an independent cohort of animals ( Killick et al . 2011 ; McLoughlin et al . 2014 ; McLoughlin et al . 2021 ). However, only three of these genes were observed to be DE in an analysis of substantially large numbers of bTB− and bTB+ cattle, and only one gene ( FOSB ) displayed increased expression ( O’Grady et al . 2025 ). Given the zoonotic threat of M. bovis , the global economic cost of bTB, the increasing prevalence of the disease in Ireland and the persistence of infection in England and Wales ( George et al . 2025 ), there is an urgent need to develop new diagnostics that can accurately discriminate between bTB− and bTB+ cattle with a high degree of specificity and sensitivity. Therefore, the objectives of this study were to leverage PB transcriptome data from bTB− cattle and cattle experimentally or naturally infected with M. bovis to train and evaluate a range of machine learning (ML) models in terms of their ability to discriminate between bTB− and bTB+ cattle, and also cattle infected with other viral and bacterial pathogens, including bovine herpes virus 1 (BoHV-1), bovine respiratory syncytial virus (BRSV) and M. avium ssp. paratuberculosis ( Fig. 1 ). Download figure Open in new tab Fig. 1: Experimental and computational workflow. Data resources for the project included peripheral blood (PB) and peripheral blood leukocyte (PBL) RNA-seq data from two studies, comprised of control non-infected animals (bTB−) and cattle naturally infected with M. bovis (bTB+) and two studies consisting of bTB− cattle and bTB+ animals experimentally infected with M. bovis across an experimental time course. Data analysis procedures included: (1) an intra-dataset differential expression (DE) and functional enrichment analysis; (2) development and tuning of machine-learning models in the training set via five-fold cross-validation using differentially expressed genes identified in the training set; (3) evaluation of hyperparameter-tuned models in the testing set and assessing performance by calculating the area under the receiver operating characteristic curve (AUROC), sensitivity and specificity, respectively; and (4) evaluation of hyperparameter-tuned models in external datasets comprised of non-infected animals and cattle infected with M. avium spp. paratuberculosis (MAP), bovine herpes virus (BoHV−1), or bovine respiratory syncytial virus (BRSV), respectively. Results Intra-dataset differential expression analysis highlights substantial heterogeneity present among bTB+ animals We performed a case-control differential gene expression analysis (DEA) to examine PB and PBL transcriptional changes caused by M. bovis infection and bTB disease. To do this, we analysed PB and PBL RNA-seq data from a total of n = 230 animals either naturally or experimentally infected with M. bovis (bTB+) and control non-infected cattle (bTB−) that were enrolled from four different studies conducted in Ireland, the UK, and the US ( Table 1 , Fig. 1 ). View this table: View inline View popup Download powerpoint Table 1: Overview of datasets used in this manuscript. We first analysed each dataset individually to characterise differentially expressed genes (DEGs) that were consistent between studies. A principal component analysis (PCA) of the top 750 most variable genes within each dataset showed variable separation of bTB+ and bTB− animals across datasets, likely reflecting differences in disease stage (e.g., acute versus chronic), experimental infection dose, age, or genetic variation due to breed composition differences. With respect to the naturally infected cattle from the MCL14-BTB dataset, we observed a distinct separation between bTB+ and bTB− animals, with PC1 accounting for 49% of the total variance explained in the PBL transcriptomes ( Fig. S1a ). Conversely, in the OGR25-BTB dataset, we observed indistinct separation of the bTB+ and bTB− cattle, with PC1 accounting for 12% of the total variance ( Fig. S1b ). For the PB MCL21-BTB time series study, after experimental infection with M. bovis , we observed weak separation of bTB+ and bTB− cattle at +2 wpi; however, a clear division between these two cohorts was evident at +10 wpi, which diminished markedly by +12 wpi ( Fig. S2 ). Similarly, for the PBL WIA20-BTB, the distinction between bTB+ and bTB− cattle was clear at +4 wpi but was substantially reduced at +10 wpi ( Fig. S3 ). Within each dataset, we subsequently conducted a DEA between the bTB+ and bTB− cattle, incorporating the top two PCs derived from the RNA-seq-called SNP variants within each dataset to account for functional genetic variation due to breed/population differences among the animals ( Fig. S4 ). We observed that for variants called across all datasets, and which were present in the O’Grady et al . (2025) imputed and filtered WGS dataset, broad population structure was captured: the coordinates of eigenvector 1 and 2 corresponding to PC1 and PC2, respectively, for the n = 123 animals described by O’Grady et al . (2025) displayed a significant and positive Spearman correlation ( ρ ) with the eigenvector coordinates for the same set of animals derived from pruned high-density SNP-array data ( Fig. S5a and Fig. S5b ) (eigenvector 1; ρ spearman = 0.97, P = 3.22 × 10 −78 ; eigenvector 2; ρ spearman = 0.68, P = 5.48 × 10 −18 ). In total, 13,897 genes were consistently expressed across all datasets, and we identified 4,533 DEGs in the MCL14-BTB dataset, with 2,223 and 2,310 genes displaying increased and decreased expression, respectively, in the bTB+ group compared to the bTB− cohort. For the OGR25-BTB dataset, we characterised 3,139 DEGs and 1,969 displayed increased expression, whereas 1,170 displayed decreased expression in the bTB+ set compared to the bTB− group. For the MCL21-BTB time series experiment, we identified 2,769 DEGs across all time points (2,402 unique DEGs). Of these, 19, 55, 226, 1,095, and 48 displayed increased expression in the bTB+ group compared to the bTB− cohort (−1 wpi) at +1 wpi, +2 wpi, +6 wpi, +10 wpi, and +12 wpi, respectively. Conversely, 3, 2, 93, 1,212, and 16 displayed decreased expression in the bTB+ group compared to the bTB− animals at +1 wpi, +2 wpi, +6 wpi, +10 wpi, and +12 wpi, respectively. For the WIA20-BTB dataset, across all time points, we identified 3,124 DEGs (2,434 unique DEGs). Of these, 1,246 and 623 displayed increased and decreased expression, respectively, in the bTB+ group at +4 wpi versus the bTB− group. At 10 wpi in the bTB+ group, 627 and 628 DEGs displayed increased and decreased expression, respectively, compared to the bTB− group ( Fig. 2a , Table S1 ). Download figure Open in new tab Fig. 2: Intra-dataset differential expression and functional enrichment analysis. a Upset plot illustrating the number of significantly ( P adj. < 0.05) differentially expressed genes (DEGs) identified in the four bovine tuberculosis (bTB) datasets and the overlap of these DEGs among the studies. b Heatmap showing the log 2 fold-change (LFC) in effect size estimates between control non-infected cattle and animals infected with M. bovis for 50 significant DEGs between bTB− and bTB+ animals identified in all naturally infected datasets and in at least one time point in each of the time course datasets, respectively. Red colours indicate increased expression, and blue colours indicate decreased expression in the bTB+ group relative to the bTB− group. The * sign indicates that the LFC value was significantly greater than or less than 0 in the corresponding comparison between bTB+ and bTB− cattle. c Bi-directional jitter plots for significantly impact pathways (SIPs) perturbed by DEGs exhibiting increased and separately, decreased expression in each dataset, respectively for four databases: (1) gene ontology (GO) biological processes (GO:BP); (2) GO cellular component (GO:CC); (3) Kyoto Encyclopaedia of Genes and Genomes (KEGG); and (4) Reactome (REAC). Dashed dotted lines indicates the −log 10 P adj. threshold of 0.05 for characterising a SIP. The total number of input DEGs and identified SIPs are also detailed for each analysis of increased and decreased DEGs in each dataset, respectively. The overlap of DEGs identified among the studies is shown in Fig. 2a , and in total, we identified 50 DEGs in the naturally infected M. bovis datasets that were also significantly DE in at least one time point of the time series experiments ( Fig. 2a , Fig. 2b ). Analysis of the log 2 fold-change (LFC) values of the DEGs for each contrast revealed a heterogeneous pattern ( Fig. 2b ); For example, the PLAUR gene was significantly DE across all nine contrasts and had a positive LFC value in the MCL14-BTB and MCL21-BTB datasets but a negative LFC value in the OGR25-BTB and WIA20-BTB datasets, respectively. Conflicting results were also observed for a range of other immune related genes including CXCL16 , IFI6 , RIGI , and CXCR4 . However, albeit not significant across all comparisons, six genes ( SESN2, VPS54 , RB1 , CSF1 , ABTB2 , and RGS1 ) possessed a positive LFC value in each contrast whereas one gene ( ADRA2A ) displayed negative expression in bTB+ animals for each contrast. We used g:Profiler to evaluate the biological pathways perturbed by DEGs separately exhibiting increased and decreased expression within each dataset. We observed that these DEGs significantly impacted ( P adj. < 0.05) multiple immune response pathways, illustrating the broad spectrum of the host response to M. bovis infection, which ranged from a weak/diminished profile to a robust innate immune response ( Fig. 2c , Table S2 ). In the MCL14-BTB dataset, many of the decreased DEGs were in significantly impacted pathways (SIPs) such as Inflammatory response , NOD-like receptor signalling pathway , and Immune system , likely reflecting an exhausted and anergic immune response to chronic M. bovis infection. For the OGR25-BTB dataset, the results were more heterogeneous, likely reflecting the diverse bTB disease states within this study cohort. For example, SIPs perturbed by the increased DEGs included classical biological processes perturbed during M. bovis infection ( Hall et al . 2021 ) such as RIG-I-like receptor signaling , ISG15 antiviral mechanism , and positive regulation of T-helper 2 cell cytokine production . Conversely, the decreased DEGs impacted the Cellular responses to stress and FCgamma receptor dependent phagocytosis signalling pathways. For the time series experiments, many of the increased DEGs were involved in biological processes reflecting an acute and robust immune response to M. bovis infection. These pathways included Response to type II interferon , Cytosolic DNA-sensing pathway , IL-17 signaling pathway , TLR4 cascade , Leukocyte activation , and NF−κB signaling ( Table S2 ). Consolidation of bTB datasets identifies transcriptional biomarkers that accurately discriminate between bTB+ and bTB− animals We consolidated the four bTB datasets and split all samples ( n = 230) into a training set (70%; n = 162; n = 69 (bTB−), n = 93 (bTB+)) and a testing set (30%; n = 68; n = 34 (bTB−), n = 34 (bTB+)). The animal IDs and associated metadata of samples split into the training and testing sets are described in Tables S3 and S4 . Consolidation of the four bTB datasets revealed the presence of a strong batch effect, whereby samples clustered in a study-specific manner, a technical artifact that was preserved in both the training and testing set ( Fig. S6 ). Within the training set, a DEA was conducted between bTB+ and bTB− cattle accounting for age, study/batch, and the top two genotype PCs inferred from SNPs called from the RNA-seq data, which were intersected and merged with the filtered imputed WGS data from O’Grady et al . (2025) ( Table S3 , Fig. S7 ). From this analysis, a total of 1,635 DEGs were identified: 949 genes exhibited increased expression and 686 displayed decreased expression in the bTB+ group relative to the bTB− cohort. Of these, 730 increased and 581 decreased DEGs, respectively, were retained for the machine learning analysis as they possessed a mean median-of-ratios expression value > 100 across samples in the training set ( Fig. 3a , Table S5 ). Significantly impacted pathways ( P adj. < 0.05) enriched for the retained increased DEGs included: Innate immune response , Regulation of tumor necrosis factor production , the Interferon-alpha production GO:BP terms, the NOD-like receptor signaling pathway KEGG term, and the Interleukin-21 signaling Reactome term. For the retained decreased DEGs, SIPs enriched for these genes included: the Regulation of macroautophagy and Regulation of interleukin-1-mediated signaling pathway GO:BP terms, and the Mitotic Spindle Checkpoint Reactome term ( Fig. 3b , Table S6 ). Download figure Open in new tab Fig. 3: Machine learning (ML) model training and evaluation. a Volcano plot illustrating significantly differentially expressed genes (DEGs) identified in the training set for the bTB+ ( n = 93) vs bTB− ( n = 69) contrast with thresholds determined by FDR- P adj. 0. Genes are coloured to indicate increased (red) or decreased (blue) expression in the bTB+ group relative to the bTB− group, respectively, and triangular datapoints indicate the genes retained for subsequent ML analysis with the criterion of mean median-of-ratios normalised expression count >100 in the training dataset. b Bi-directional jitter plot for significantly impacted pathways (SIPs) perturbed by DEGs exhibiting increased and separately, decreased expression in the training set for four databases: (1) gene ontology (GO) biological processes (GO:BP); (2) GO cellular component (GO:CC); (3) Kyoto Encyclopaedia of Genes and Genomes (KEGG); and (4) Reactome (REAC). Dashed dotted lines indicate the −log 10 P adj. threshold of 0.05 for characterising a SIP. The total number of input DEGs and identified SIPs are also detailed for each analysis of increased and decreased genes in each dataset, respectively. c Area under the receiver operating characteristic curve (AUROC) values for eight hyperparameter-tuned ML models in each of the five folds of the training set. d An AUROC depicting the performance of each of the hyperparameter-tuned models in the testing set. e Boxplot showing the distribution of AUROC estimates of animals in the testing dataset depending on whether the dataset comprised of animals experimentally or naturally infected with M. bovis . f Heatmap showing the predicted class of animals in the testing set when selecting an optimum threshold to maximise the sensitivity (> 0.90), where possible, of each classifier. We trained and performed hyperparameter tuning of eight models (unpenalised logistic regression [GLM], logistic regression with either a lasso [LASSO], ridge [RIDGE] or elastic net [ENET] penalisation, respectively, random forest [RF] and random forest specifying the extra trees algorithm [RF-ET], naïve Bayes [NB], and a multi-layered perceptron [MLP]) in the training set using five-fold cross-validation. The selected hyperparameters for each model are described in Table S7 . The performance of each tuned model across the five folds is illustrated in Fig. 3c and described in Table S8 . The best performing model was the ENET model, which achieved an average AUROC of 0.968 (range: 0.898–0.993). Other models that performed almost as well included the LASSO (0.958, [0.906–0.981]), the MLP (0.955, [0.891–0.987]), the RIDGE (0.949, [0.895–0.974]), the RF-ET (0.926, [0.869–0.949]), and the RF (0.917, [0.859–0.955]). The NB classifier produced a substantially lower average AUROC value of 0.766 (0.679–0.850). The poorest-performing model in the training set was the GLM model, which gave an average AUROC value across the five folds of 0.545 (0.469– 0.696) ( Fig. 3c , Table S8 ). Evaluation of the hyperparameter-tuned models in the testing set illustrated that all the models generalised well, achieving comparable performance to that observed in the training set ( Fig. 3d ). The best-performing models were the ENET and the LASSO models, which both achieved an AUROC of 0.938 (95% confidence interval (CI): 0.87–1.00). These models were followed closely by the MLP (0.933, 95% CI: 0.87–1.00), RF-ET (0.921, 95% CI: 0.85–0.99), RIDGE (0.910, 95% CI: 0.83–0.99), and RF (0.898, 95% CI: 0.82–0.98) models, respectively. The NB model, similar to observations in the training set, gave substantially lower performance in the testing set, obtaining an AUROC of 0.826 (95% CI: 0.74–0.92). Again, consistent with findings from the training set analysis, the worst-performing model was the GLM model, which performed essentially like random guessing (AUROC = 0.501, 95% CI: 0.43–0.57). Taking the AUROC value as the performance metric, we observed that for all models except the NB and GLM classifiers, the ability to classify animals correctly as bTB+ or bTB−, respectively, in the experimentally infected datasets tended to be higher than for the naturally infected datasets; however, this difference was not statistically significant (Wilcoxon paired rank-sum test; P = 0.64; Fig. 3e ). Current bTB diagnostics suffer from a lack of sensitivity. We therefore exploited the ROC curves depicted in Fig. 3d and for each model, selected an optimum probability threshold ( θ ) for classifying an animal as bTB+ that yielded a sensitivity value > 90% where possible, while also balancing specificity estimates ( Table S9 ). At these selected θ values, the models that obtained the joint highest sensitivity estimates were the GLM and ENET classifiers, each achieving a sensitivity value of (0.941); however, the GLM model misclassified a high proportion of bTB− animals as bTB+ (31 out of 34 bTB− animals), yielding a specificity estimate of 0.088 ( Fig. 3f , Table S10 ). Overall, the models that achieved the optimum balance between sensitivity and specificity were the ENET classifier (sensitivity = 0.941, specificity = 0.853), followed closely by the LASSO, MLP, and RF-ET models, respectively (sensitivity = 0.912, specificity = 0.882). Similar sensitivity estimates were observed for the RF and RIDGE models, with reduced specificity estimates (0.853 and 0.824, respectively). The NB model obtained comparable sensitivity estimates to the LASSO, MLP, RF-ET, RF and RIDGE models (0.912) but had substantially reduced specificity estimates (0.706). The ENET model incorrectly classified two out of 34 bTB+ animals (5.88%) as bTB−, with both animals being derived from the OGR25-BTB dataset. All models, apart from the GLM classifier, correctly classified the animals experimentally infected with M. bovis as bTB+ from the MCL21-BTB time series dataset, which included 3 bTB+ animals sampled at +1 wpi and +2 wpi ( Fig. 3f , Table S10 ). This is an important result because the current IGRA bTB diagnostic test has a minimum +2 wpi chronological detection limit, because an adaptive immune response to mycobacterial antigens is required for the assay to function correctly ( Buddle et al . 1995 ). A greedy forward search algorithm identifies robust discriminatory gene sets in the training set with poor generalisability in the testing set Previous studies have shown the utility of employing a greedy forward search (GFS) algorithm for identifying a small number of discriminatory genes between non-infected individuals and individuals diagnosed with a range of diseases, including sepsis ( Sweeney et al . 2015 ), dengue fever ( Robinson et al . 2019 ), and hTB using both microarray ( Sweeney et al . 2016 ) and RNA-seq data ( Chang et al . 2024 ), respectively. We therefore implemented the GFS algorithm on the entire set of retained DEGs (1,311 genes; Fig. 3a ). This analysis identified a set of 13 genes ( SATB1 , ALAS1 , NR4A2 , SMG8 , ELP2 , DCAF5 , ZNF75D , PRICKLE3 , MYB , BCDIN3D , TEX2 , MAD2L2 , and KATNBL1 ; 13-gene set, Fig. 4a ), which achieved an average AUROC across the five folds of 0.996 (range: 0.990–1.000) ( Fig. 4a , Fig. 4b , Table S11 ). We subsequently reimplemented the GFS strategy on the set of retained DEGs without the genes identified in the first initialisation and identified a set of 17 genes ( PTK2 , TNFRSF13B , UBE2G1 , SENP3 , UQCC2 , TAP2 , MON1A , NFE2L3 , GTSE1 , TRIQK , ADPGK , GDPG1 , CCR5 , ZNF628 , NABP1 , GEMIN5 , and GDAP2 ; 17-gene set, Fig. 4a ) which achieved near perfect discrimination across the five folds in the testing set, with an average AUROC across the five folds of 0.999 (range: 0.993–1.000) ( Fig 4a , Fig. 4b , Table S12 ). We also evaluated a combination of the 13-gene and 17-gene sets (30-gene set), respectively, which yielded an average AUROC across the five folds in the training set of 0.986 (range: 0.975– 1.000) ( Table S13 ). Download figure Open in new tab Fig. 4: Greedy-forward search strategy evaluation in training and testing sets. a Line plot showing the average area under the receiver operating characteristic curve (AUROC) across the five folds in the training set for each combination of genes identified in the first pass or second pass of the greedy-forward search algorithm, respectively. b Line plot showing the 95% confidence interval of AUROC estimates for the 13-gene set, 17-gene set and combined 30-gene set in each of the five cross-validation folds in the training set. c Boxplots showing the distribution of the scaled infection z -score derived from each of the three gene sets for all bTB− (blue; circle) and bTB+ (red; triangle) animals in the training set, separated by study of origin. d Boxplots showing the distribution of the scaled infection z -score derived from each of the three gene sets for all bTB− (blue) and bTB+ (red) animals in the testing set, separated by study of origin. e An AUROC depicting the performance of each of the three gene sets in the testing set. f Heatmap showing the predicted class of animals in the testing set when selecting an optimum threshold to maximise the sensitivity (> 0.85) for each of the three gene sets. Evaluation of the infection z -scores inferred from each of the gene sets across all four cohorts in the training dataset revealed significant differences between bTB+ and bTB− cattle, with all intra-dataset-condition pairs showing a significant difference in the inferred infection z- score (Wilcoxon rank-sum test, P adj. < 0.05) ( Fig. 4c , Table S14) . For each gene set, we next inferred the infection z- score of bTB+ and bTB− animals in the testing set and observed significant differences in the inferred bTB scores for these conditions in the MCL21-BTB and OGR25-BTB datasets (Wilcoxon rank-sum test, P adj. < 0.05); however, this difference was not statistically significant in the MCL14-BTB dataset ( P adj. = 0.20) owing to the small sample size ( n = 2 bTB+, n = 3 bTB−) ( Fig. 4d , Table S15 ). We were unable to compare the bTB score between bTB+ and bTB− animals from the WIA20-BTB cohort in the testing set because the original train/test splitting of the animals yielded no bTB+ animals from this cohort in the testing set. Evaluation of the predictive capabilities of the three gene sets revealed significantly lower AUROC estimates in comparison to the training set with the 13-gene set, the 17-gene set and the 30-gene set achieving AUROC values of 0.856 (95% CI: 0.77−0.95), 0.887 (95% CI: 0.80−0.97), and 0.900 (95% CI: 0.82−0.98), respectively ( Fig. 4e ). Following the same approach used for the analysis of the other ML algorithms ( Fig. 3f ), for each gene set, we selected the threshold that yielded a sensitivity value > 0.85 and that also maximised the specificity ( Table S16 ). We selected a lower sensitivity value in this analysis as selecting a sensitivity value > 0.9 for all three gene sets yielded very low specificity estimates (< 0.65 for all sets; Table S16 ) Using these selected threshold values, the 13-gene set achieved a sensitivity value of 0.882 and a specificity value of 0.559, with a large number of bTB− cattle being misclassified as bTB+ (15 out of 34 bTB− animals). The 17-gene and 30-gene set signatures both achieved balanced sensitivity and specificity estimates of 0.853 for both metrics in both classifiers ( Fig. 4f , Table S17 ). Similar to the other ML models, at these thresholds, the 13-gene and 30-gene sets correctly classified the animals experimentally infected with M. bovis at +1 wpi and +2 wpi, respectively, in the MCL21-BTB dataset as bTB+. Evaluation on comparator pathogen datasets demonstrates the robustness of biosignatures To evaluate the specificity of the bTB ML models and inferred bTB gene sets, we obtained RNA-seq data for three other infectious diseases, including Johne’s disease caused by M. avium spp. paratuberculosis (MAP) ( n = 3 MAP−, n = 11 MAP+) and bovine respiratory disease (BRD) caused by infection with either BoHV-1 ( n = 6 BoHV-1−, n = 12 BoHV-1+) or BRSV ( n = 6 BRSV−, n = 12 BRSV+), respectively ( Table 1 , Fig. 1 ). Within each of these datasets, we first evaluated the ability of the eight ML classifiers to discriminate between control non-infected animals and cattle infected with each of the three pathogens. With regards to the MAP dataset, we observed that all eight ML models performed poorly, with the RIDGE model achieving the highest AUROC value of 0.667 and the RF-ET model achieving the lowest AUROC of 0.455, no better than random allocation; however, wide confidence interval estimates were inferred because of the small number of animals in this dataset ( Fig. 5a ). For animals infected with BoHV-1, we observed that the models displayed improved discriminative ability, achieving a maximum AUROC from the ENET model of 0.931 (95% CI: 0.79–1.00), whereas the NB and GLM models achieved the lowest AUROC of 0.50 (95% CI: 0.5−0.5) ( Fig. 5b ). For cattle experimentally infected with BRSV, we observed high discriminatory ability with the LASSO, RIDGE, RF, RF−ET, RIDGE, and ENET models, all achieving an AUROC ≥ 0.944, whereas as the GLM and NB models achieved an AUROC of 0.667 and 0.500, respectively ( Fig. 5c ). Download figure Open in new tab Fig. 5: Evaluation of machine learning and greedy forward search strategies in external datasets. a An area under the receiver operating characteristic curve (AUROC) showing the ability of the hyperparameter-tuned models to discriminate between control non-infected animals and cattle infected with M. avium spp. paratuberculosis (MAP) from the ALO19-MAP dataset. b An AUROC showing the ability of the hyperparameter-tuned models to discriminate between control non-infected animals and cattle infected with bovine herpes virus (BoHV-1) from the ODO23-BRD dataset. c An AUROC showing the ability of the hyperparameter-tuned models to discriminate between control non-infected animals and cattle infected with bovine respiratory syncytial virus (BRSV) from the JOH21-BRD dataset. d Boxplots showing the distribution of the scaled infection z -score derived from each of the three gene sets for all control non-infected animals and cattle infected with either MAP, BoHV-1, or BRSV, respectively. e Three AUROCs illustrating the performance of the three gene sets in terms of discriminating between control non-infected animals and cattle infected with either MAP, BoHV−1, or BRSV, respectively. f An AUROC curve showing the discriminatory power of all hyperparameter-tuned ML models and gene sets when comparing bTB+ cattle to MAP+, BoHV-1+ and BRSV+ cattle, respectively. For each external disease dataset, we next inferred an infection z -score for all non-infected and infected animals using the 13-gene set biosignature inferred from the first iteration of the GFS algorithm on the bTB training dataset, the 17-gene set signature from the second iteration, and the combined 30-gene set inferred from both iterations ( Fig. 4a , Fig. 5d ). We observed significant differences (Wilcoxon rank-sum test; P adj. < 0.05) between the inferred infection z -scores from all three gene set signatures for control non-infected animals and cattle infected with BoHV-1( Fig. 5d , Table S18 ). Conversely, this difference was not significant between MAP+ and MAP− animals for any of the three gene set signatures ( P adj. ≥ 0.586), or between BRSV+ and BRSV− animals using either the 13- or 30-gene set signatures, respectively ( P adj. ≥ 0.120) ( Table S18 ). However, the difference in the inferred infection z -score between BRSV+ and BRSV− animals was significant based on the 17-gene set signature ( P adj. = 1.9 × 10 −3 ) ( Fig. 5d , Table S18 ). We next inferred the discriminatory ability of these gene sets between control non-infected and infected animals in each external dataset using ROC curves ( Fig. 5e ). We observed that the three gene sets performed poorly at discriminating between MAP+ and MAP− animals, achieving an AUROC ≤ 0.667; however, these estimates were accompanied by wide confidence intervals. For the BoHV-1 dataset, similarly to the observations from the ML models, we observed better discriminatory capacity between BoHV-1− and BoHV-1+ animals with the 30-, 17- and 13-gene set signatures each achieving AUROC estimates of 0.861, 0.931, and 1.00 (perfect discrimination), respectively ( Fig. 5e ). With respect to the BRSV dataset, we observed high discriminatory power for the 17-gene set signature, which achieved an AUROC of 0.972 (95% CI: 0.91–1.00). Conversely the 13- and 30-gene set signatures performed poorly at differentiating between BRSV− and BRSV+ cattle, yielding AUROC estimates of 0.778 and 0.639, respectively ( Fig. 5e ). We lastly evaluated the capacity of all eight ML models and the three inferred GFS gene sets to discriminate between bTB+ animals from the testing set, and animals infected with MAP, BoHV-1, or BRSV, respectively. We observed variable performance among ML models and gene sets in their ability to differentiate between bTB+ animals and animals infected with the other infectious agents ( Fig. 5f ). The best-performing model/gene set was the MLP model. This model yielded an AUROC of 0.868 (95% CI: 0.78–0.96), indicating that it is highly specific in differentiating between animals with bTB and animals with other bacterial and viral infections. This was closely followed by the ENET model (the best performing model for distinguishing between bTB+ and bTB− cattle; Fig. 3c ), which achieved an AUROC of 0.826 (95% CI: 0.72–0.93) ( Fig. 5f ) The 30- and 13-gene set signatures, and the RF-ET, RIDGE, and RF models, performed comparably well, obtaining AUROC estimates ranging from 0.765−0.819. We subsequently observed a drop-off in performance for the LASSO and 17-gene set signature, which achieved AUROC values of 0.711 and 0.584, respectively. The NB and unpenalized GLM algorithms performed extremely poorly, achieving AUROC estimates < 0.5 ( Fig. 5f ). Discussion We present a comprehensive systematic evaluation and integration of PB and PBL RNA-seq data from cattle naturally or experimentally infected with M. bovis (bTB+) and control non-infected animals (bTB−). It is well established that the peripheral immune responses of cattle infected with M. bovis reflect those at the site of infection ( Rhodes et al . 2000 ). In this regard, we have shown that the peripheral transcriptomic response following M. bovis infection can exist on a spectrum, ranging from an anergic/exhausted response in animals chronically infected with M. bovis to an immune signature that is biased towards an innate response in animals shortly after being experimentally infected with M. bovis . Through consolidation of multiple bTB RNA-seq datasets and using eight different ML algorithms and a GFS strategy, we observed that accurate discrimination of bTB+ and bTB− cattle can be achieved with high sensitivity and specificity, and gave results that generalise well to unseen testing data. We also showed that a 30-gene set signature, identified through a greedy-forward search algorithm, and an ENET classifier are robust in differentiating between bTB+ cattle and animals infected with either MAP, BoHV−1 or BRSV. Divergent peripheral blood transcriptional profile reveals a broad spectrum of bTB disease in cattle Analysis of RNA-seq data from bTB+ and bTB− cattle within each bTB dataset identified thousands of DEGs ( Fig. 2a ), indicating that the peripheral transcriptome is substantially perturbed following M. bovis infection—an observation that is consistent with previous work focused on bovine alveolar macrophages (bAMs) experimentally infected with M. bovis and M. tuberculosis ( Nalpas et al . 2015 ; Malone et al . 2018 ; Hall et al . 2021 ). However, we observed that 35.1–47.2% of DEGs identified in each bTB dataset were cohort-specific, suggesting that analysis of varying bTB disease states (e.g., acute versus chronic) can result in different transcriptional patterns being observed. In this regard, the DEGs that displayed decreased expression in bTB+ cattle relative to the bTB− cattle in the MCL14-BTB dataset, were enriched in pathways that are typically associated with the initial host-pathogen interaction, as reflected, for example, in the proinflammatory responses observed in bAMs infected with M. bovis ( Nalpas et al . 2015 ; Malone et al . 2018 ; Hall et al . 2021 ) ( Fig. 2c ). Our results for bTB+ cattle naturally infected with M. bovis are consistent with a previous microarray study that also detected suppression of key innate immune genes in naturally infected cattle ( Meade et al . 2007 ). Conversely, for the WIA20-BTB and MCL21-BTB datasets, which encompassed animals that were recently experimentally infected with M. bovis , the DEGs that displayed increased expression were significantly enriched in SIP sets that were dominated by interferon-mediated signalling pathways, which is consistent with an acute phase of infection that has been reported in PB data for hTB ( Singhania et al . 2018 ; Tabone et al . 2021 ). Our analysis of SIPs enriched by DEGs exhibiting increased and decreased expression in the bTB+ cohort compared to the bTB− cohort revealed conflicting gene expression patterns in the OGR25-BTB dataset. Many of the SIPs enriched for genes showing increased expression following M. bovis infection were classical innate immune response and cytosolic sensing pathways, including the RIG-I-like receptor signalling pathway , the Type 1 interferon-mediated signalling pathway , and the ISG15 antiviral mechanism Reactome term, all of which play key roles in the host response to mycobacterial challenge ( Bogunovic et al . 2012 ; Chai et al . 2020 ; Vance 2025 ). On the other hand, the Chemokine signaling , the FCgamma receptor dependent phagocytosis and Cellular responses to stress pathways were enriched for DEGs exhibiting decreased expression in the bTB+ cohort ( Fig. 2c ; Table S2 ), indicating suppression of these pathways. Collectively, these findings indicate that the animals from the OGR25-BTB dataset are clinically heterogeneous, and this variation is reflected in the PB transcriptome, which is unsurprising given the large number of bTB+ animals ( n = 63) and lack of information regarding the chronology of M. bovis infection. In total, across all four bTB datasets, given the observed transcriptional heterogeneity, we characterised only 50 shared DEGs, the majority of which displayed substantial inter-dataset differences in terms of inferred effect sizes (LFC values) between bTB− and bTB+ cattle, respectively ( Fig. 2b ). These findings show the ineffectiveness of solely assessing single or small numbers of DEGs among different bTB studies as potential biomarkers for discrimination of bTB+ and bTB− cattle. Accurate discrimination of bTB+ and bTB− cattle is possible using peripheral blood RNA-seq data We trained a total of eight ML models using 1,311 DEGs identified in the consolidated training set encompassing animals from all four bTB datasets ( Fig. 3a ). Using a five-fold cross validation framework for hyperparameter tuning, we observed high average AUROC values > 0.90 for the penalised logistic regression models (ENET, LASSO, and RIDGE), and the MLP, RF and RF-ET models, which generalised well and achieved comparable performance in the testing dataset ( Fig. 3c , Fig. 3d ), indicating that the PB transcriptome can accurately discriminate bTB+ and bTB− cattle. We also observed a trend whereby the best performing models yielded higher AUROC estimates in the experimentally infected compared to the naturally infected groups; however, this difference was not significant ( P = 0.64, Fig. 4e ). This observation may indicate that the inferred RNA biosignatures are being driven by an innate immune-biased signature, as the transcriptional profiles of animals recently experimentally infected with M. bovis exhibited perturbed innate immune response pathways ( Fig. 2c ). This is supported by the observation that all experimentally infected animals from the MCL21-BTB dataset could be classified as bTB+ when optimum thresholds were selected to yield sensitivity values ≥ 0.90 in the testing set ( Fig. 3d , Fig. 3f , Table S10 ). Moreover, the best performing model (ENET), which achieved a sensitivity and specificity of 0.941 and 0.853, respectively, used expression information from 273 (20.8%) out of all input DEGs to determine the probability that an animal was bTB+. This biosignature was composed of key innate immune response genes, including TLR4 , MX1 , IFNGR1 , IFITM3 , GBP2, GBP4 , TNFAIP2 , IFI6 , and UBASH3A ; these genes have been shown to have pivotal roles in the host response to mycobacterial infections or in transcriptional biosignatures for hTB using whole blood ( Sánchez et al . 2010 ; Ranjbar et al . 2015 ; Hall et al . 2019 ; Gupta et al . 2020 ; Akter et al . 2022 ; Sun et al . 2024 ). We also generated three infection z -scores derived from three gene sets that were identified with a GFS algorithm ( Sweeney et al . 2015 ; Sweeney et al . 2016 ; Chang et al . 2024 ), which resulted in a 13-gene set signature, a 17-gene set signature, and a combined 30-gene set signature ( Fig. 4a , Fig. 4b ). These three gene sets achieved high discrimination between bTB− and bTB+ cattle in the training set with AUROC values ≥ 0.96 across the five folds ( Fig. 4b , Fig. 4c ). For the training set, these gene sets achieved slightly lower but highly accurate discriminatory power between bTB+ and bTB− cattle, with all combinations achieving an AUROC > 0.85. The most discriminatory gene that contributed to the infection z -score identified in the 13-gene set was SATB1 , which has previously been shown to be upregulated in macrophages upon infection with virulent M. tuberculosis , contributing to enhanced bacterial survival through repression of NADPH oxidase 2 (NOX2 – the catalytic core of NADPH oxidase encoded by CYBB ) with a concomitant downregulation of bactericidal reactive oxygen species (ROS) production ( Subuddhi et al . 2020 ). The PTK2 gene encoding the focal adhesion kinase (FAK) protein was identified as the most discriminatory gene in the 17-gene set ( Fig. 4a ). The concentration of FAK in macrophages infected with M. tuberculosis has been shown to decrease in a time-dependent manner orchestrated by the bacterium through inhibition of ROS production, culminating in increased macrophage necrosis and uncontrolled replication of M. tuberculosis ( Afriyie-Asante et al . 2021 ). Other genes identified by the two GFS iterations included NR4A2 , a transcription factor that regulates macrophage polarity from an M1 phenotype to an M2 phenotype ( Mahajan et al . 2015 ), and which was previously identified as being upregulated in monocytes derived from non-infected household contacts heavily exposed to M. tuberculosis in comparison to monocytes isolated from individuals latently infected with M. tuberculosis ( Simmons et al . 2022 ); CCR5 , which was upregulated in macrophages infected with M. tuberculosis resulting in enhanced IL-10 production ( Das et al . 2014 ); and UBASH3A described previously. These results indicate that ML models and the GFS approach have identified genes expressed in PB, which can effectively discriminate between bTB+ and bTB− cattle, and that have mechanistic roles in mycobacterial infection biology and host-pathogen interaction. Discrimination of bTB+ cattle from animals infected with other infectious agents highlights shared and divergent peripheral responses Analysis of ML models for other infectious diseases revealed variable discrimination between infected and control non-infected animals within each dataset. For example, in the ALO19-MAP dataset, the RIDGE classifier was the best performing model, achieving an AUROC of 0.667 (95% CI: 0.01–1.00) indicating poor to moderate discriminatory capacity ( Fig. 5a ). Although infection with MAP upregulated the expression of key immune-related genes such as CXCL8 and IFI27 , and perturbed common biological pathways observed in the bTB dataset, such as Defence response and Immune response ( Alonso-Hearn et al . 2019 ), these results indicate that MAP induces a different PB transcriptional profile in cattle compared to M. bovis , which resulted in the observed suboptimal performance of ML models in this dataset. Unlike MAP infection, we observed substantially higher discriminatory power between BoHV-1− and BoHV-1+ cattle, respectively, with the ENET model achieving high discriminatory capacity (AUROC = 0.931; 95% CI: 0.79−1.00) followed closely by the RF, RIDGE, LASSO, RF−ET and MLP models, which achieved AUROC values ranging from 0.847−0.917 ( Fig. 5b ). Additionally, we observed even higher discrimination between BRSV− and BRSV+ animals, ranging from perfect discrimination for the LASSO, MLP and RF-based models (AUROC = 1.00) to excellent performance for the RIDGE and ENET models with these two models each achieving AUROC estimates of 0.986 and 0.944, respectively. These findings suggest that infection with either BoHV-1 or BRSV induces transcriptional perturbations in PB of cattle that are similar to those observed after infection with M. bovis . These observations are supported by the biological pathways perturbed by DEGs in the BoHV-1 and BRSV datasets reported by O’Donoghue et al . (2023) and Johnston et al . (2021) , respectively. For the BoHV-1 dataset, SIPs included the NOD - like receptor signalling pathway, the Defence response to virus GO term, the Interferon signalling Ingenuity Pathway Analysis (IPA ® ) term, and the Regulation of viral process GO term ( O’Donoghue et al . 2023 ). For the BRSV dataset, SIPs included the Defence response to virus , Regulation of viral life cycle , and Innate immune response GO terms ( Johnston et al . 2021 ). We also observed somewhat comparable variability across datasets in terms of discriminatory power between non-infected and infected cattle for the three gene sets inferred from the GFS algorithm (13-gene set, 17-gene set, and 30-gene set). The three gene sets achieved poor discrimination in the MAP dataset with AUROC values ranging between 0.545−0.667, extremely high discrimination in the BoHV-1 (AUROC range: 0.861–1.00) and extremely poor to excellent discriminatory power in the BRSV cohort (AUROC range: 0.639–0.972) ( Fig. 5e ). These results indicate that the gene sets inferred using the GFS algorithm closely resemble peripheral transcriptional perturbations of these genes induced by BoHV-1 but are distinct from those observed in MAP- and BRSV-infected cattle compared to healthy controls. Lastly, we evaluated the ability of all ML models and GFS gene sets to discriminate between bTB+ animals and MAP+, BRSV+, and BoHV-1+ cattle as alternative disease-positive controls. We observed varying discriminatory power across all ML models and gene sets with AUROC values ranging from 0.367 for the unpenalized logistic regression model to 0.868 for the MLP model ( Fig. 5f ). The second and third-ranked classifiers were the ENET classifier and the 30-gene set signature, which achieved AUROC estimates of 0.826 and 0.819, respectively. These results are noteworthy as the ENET classifier and the 30-gene set signature both achieved a specificity of 0.853 and sensitivity of ≥ 0.853 for the bTB testing group ( Fig. 3f and Fig. 4f ). The performances of these classifiers observed here indicate that they both have the capacity to effectively differentiate between bTB+ and bTB− cattle and are also robust at discriminating between bTB+ animals and cattle infected with other bacterial and viral pathogens. Transcriptomic biosignatures and classifiers for augmenting bTB control and eradication programmes In current bTB endemic countries, control and eradication of the disease is underpinned primarily by the implementation of test and slaughter programmes, whereby infected animals are identified and removed from the herd ( McCormack 2023 ; More 2024 ). The in vivo SICTT diagnostic is used as a population screening tool for bTB prevalence due to its high specificity (> 0.99 in Britain), and positive predictive value (PPV – the probability that a positive result is truly positive; > 0.91 in Britain) ( Goodchild et al . 2015 ). In addition, the PPV improves as disease prevalence rises because confidence increases that a positive result is a true positive ( Altman & Bland 1994 ; Monaghan et al . 2021 ). The ancillary in vitro IGRA test is employed typically when ≥ 5 SICTT reactor animals in a herd are identified, and is interpreted in parallel with the SICTT to maximise the detection of infected cattle in these high-risk herds ( Wood et al . 1991 ; European Food Safety Authority 2012 ; Madden et al . 2024 ). In this regard, the sensitivity of the SICTT and IGRA combined is estimated to be 93% ( Gormley et al . 2006 ). The SICTT and IGRA tests are based on cell-mediated immunological responses and, as such, are limited to the detection of animals at least +2 wpi in the case of the IGRA test ( Buddle et al . 1995 ), and between +3 and +6 wpi for the SICTT ( de la Rua-Domenech et al . 2006 ). Although the sample size is limited in the testing set, our ENET model and 30-gene set signature correctly classified experimentally infected animals as being bTB+ at +1 and +2 wpi, respectively ( Fig. 3f , Fig. 4f ), indicating that these classifiers have different diagnostic properties compared to the SICTT and IGRA test; these observations suggest that a third diagnostic capable of detecting early disease cases, may increase the ability to identify all bTB+ animals in a herd and contribute significantly to control and eradication programmes. Notwithstanding these encouraging results, it is important to emphasise that to validate these observations, a much larger PB transcriptomics study using cattle experimentally infected with M. bovis at several early post-infection time points will be required. While the bTB+ animals at +1 wpi and +2 wpi were correctly classified as such in the testing set by the ENET and 30-gene set classifiers, the two bTB+ animals that were consistently misclassified as bTB− based on the selected probability threshold values for these two models were T059 and T062, both of which were derived from the OGR25-BTB dataset ( Table S9 and S16 ). Interferon-γ release assay test results demonstrated that T059 yielded a positive result, whereas T062 produced a negative result for this diagnostic ( O’Grady et al . 2025 ). These results indicate that while the ENET and 30-gene set classifiers, respectively, augment current bTB diagnostics through correctly classifying animals infected with M. bovis ≤ +2wpi as bTB+, they illustrate that training alternative PB transcriptomic models using data from SICTT+/IGRA−/culture+ animals or SICTT−/IGRA−/culture+ animals, respectively, is necessary, and represents a parallel approach to strengthen current bTB testing strategies further and to address the persistent issue of residual infection in bTB breakdown herds ( Wiseman et al . 2024 ). Similar to the in vitro IGRA test, the ENET classifier and infection z -score thresholds from the 30-gene set signature can be readily adjusted in high-risk herds where test specificity is less of an issue, because maximisation of sensitivity in these populations can aid in removing all infected animals, which will substantially improve control and eradication efforts ( Gormley et al . 2006 ). However, a larger cohort of naturally infected animals will be required to validate the findings and diagnostic estimates obtained in this study. This may be achieved through additional high-throughput transcriptomic analysis, which has seen a dramatic reduction in costs over recent years. In addition to conventional bulk RNA-seq, targeted transcriptomic methods could also be used ( Kulkarni 2011 ; Li et al . 2015 ). In this regard, a targeted RNA-seq strategy was used to validate RNA-seq-derived biomarkers for Lyme disease in humans ( Servellita et al . 2022 ). Methods Data acquisition and preprocessing For this study, publicly available raw RNA-seq data were obtained from six experiments investigating transcriptional differences in peripheral blood (PB) or peripheral blood leukocyte (PBL) samples between control non-infected cattle and animals infected naturally or experimentally with: 1) M. bovis (bTB+), encompassing three datasets reported by McLoughlin et al . (2014) (MCL2014-BTB), Wiarda et al . (2020) (WIA20-BTB), McLoughlin et al . (2021) (MCL21-BTB); 2) Bovine respiratory syncytial virus (BRSV) generated by Johnston et al . (2021) (JOH21-BRD); 3) Bovine herpes virus (BoHV-1) analysed by O’Donoghue et al . (2023) (ODO23-BRD) and; 4) M. avium ssp. paratuberculosis (MAP) reported by Alonso-Hearn et al . (2019) (ALO19-MAP) ( Table 1 ). These raw datasets were downloaded using fastq-dl v. 2.0.4 ( Petit III et al . 2025 ). Lastly, an additional pre-processed RNA-seq count matrix file derived from n = 63 control (bTB−) and n = 60 (bTB+) cattle from O’Grady et al . (2025) (OGR25-BTB) was downloaded from the Gene Expression Omnibus (GEO) database using the BioProject accession number GSE255724 ( Table 1 ). The quality of the raw RNA-seq data from MCL14-BTB, WIA20-BTB, MCL21-BTB, JOH21-BTB, ODO23-BRD and ALO19-MAP datasets was assessed using fastqc (v. 0.11.9) ( Andrews 2010 ) and multiqc (v. 1.25.1) ( Ewels et al . 2016 ). Illumina adaptors were trimmed using trimmomatic (v. 0.39) ( Bolger et al . 2014 ) and raw reads were further filtered to retain reads ≥ 36 bp in length and remove trailing bases with a phred-scale quality score < 30. The trimmed RNA-seq reads were aligned to the ARS-UCD1.2 bovine reference genome ( Rosen et al . 2020 ) using STAR (v. 2.7.1) ( Dobin et al . 2013 ). The reference genome index was generated using the genomeGenerate command with the Bos taurus ARS-UCD1.2 Ensembl annotation file (v. 110) specified in addition to the -- sjdbOverhang 99 parameter. Read counts for each gene were quantified using featureCounts (v. 2.0.0) ( Liao et al . 2014 ). Variant calling from RNA-seq data and population genomics analyses To call variants directly from the RNA-seq data, we followed the Pig Genotype Tissue Expression (GTEx) Consortium’s pipeline ( https://github.com/FarmGTEx/PigGTEx-Pipeline-v0/blob/master/02_RNA-Seq/03_SNP_calling.smk ) ( Teng et al . 2024 ) and used the Genome Analysis Tool Kit (GATK) (v. 4.3.0.0) ( McKenna et al . 2010 ). Note: imputed whole-genome sequence data (WGS) comprising a total of 3,866,506 filtered autosomal SNPs (minor allele frequency (MAF) ≥ 5%, non-significant deviation from Hardy-Weinberg Equilibrium (HWE) P ≥ 1 × 10 −6 , imputed dosage R 2 ≥ 0.6) was available for the n = 123 animals used in the experiment conducted by O’Grady et al . (2025) ( https://zenodo.org/records/13752056 ). Briefly, for the RNA-seq data derived from the time series experiments ( Wiarda et al . 2020 ; McLoughlin et al . 2021 ), we first merged animal-specific BAM files across time points to increase read coverage for the variant call. Next, for all independent BAM samples across the RNA-seq-only experiments ( n = 39), we applied the GATK Best Practices guidelines for variant discovery ( DePristo et al . 2011 ; Van der Auwera G.A. & O’Connor B.D. 2020 ). We subsequently filtered out low-quality variants using the options: FisherStrand (FS) test > 30.0, QualByDepth (QD) < 2.0, read depth (DP) < 4.0, and restricted our analysis to biallelic autosomal SNPs. We subsequently merged the sample and study-specific VCF files using the BCFtools (v. 1.15.1) ( Danecek et al . 2021 ) merge command. Overall, this analysis yielded a total of 163,589, 989,438 and 1,482,214 raw variants called in the MCL14-BTB, WIA20-BTB and MCL21-BTB datasets, respectively. We next intersected the above RNA-only VCF files with the filtered WGS data from O’Grady et al . (2025) using the BCFtools isec command with options ‘ - n +4’ and ‘-c none’ to retain a total of 31,312 SNPs called/present across all datasets. Next, we merged all datasets and filtered the 31,312 SNPs for those with a MAF > 5% and that did not deviate from HWE ( P < 1 × 10 −6 ) and that were not missing in more than 20% of samples. This procedure resulted 29,740 SNPs being retained. From these, we extracted all the samples from the MCL14BTB, WIA20-BTB, MCL21-BTB and OGR25-BTB datasets, respectively, into study specific VCF files. The filtered VCF files resulting from the workflow detailed above were loaded into R (v. 4.3.2) ( R Core Team 2025 ) with the SNPRelate (v. 1.34.1) toolset ( Zheng et al . 2012 ) using the snpgdsVCF2GDS and snpgdsOpen functions. A principal component analysis (PCA) was conducted on each SNP dataset using the snpgdsPCA function, and the eigenvectors corresponding to principal component 1 (PC1) and PC2 were plotted using ggplot2 (v. 3.4.4) ( Wickham 2016 ). As a validation step, we also compared the eigenvector coordinates of PC1 and PC2 derived from the 29,740 SNPs in the OGR25-BTB dataset to the eigenvector coordinates for PC1 and PC2 from the same set of animals derived from 34,272 pruned genome-wide array SNPs initially reported by O’Grady et al . (2025) using the Spearman correlation ( ρ spearman ). Intra-bTB dataset differential expression analysis A DEA was conducted separately between the reactor (bTB+) and control (bTB−) animal groups for each bTB dataset using DESeq2 (v. 1.40.2) ( Love et al . 2014 ). For the generalised linear model, a design matrix was specified, which included the following covariates: age in months, where applicable; RNA-seq sequencing batch, where appropriate; and population genetic structure in the form of PC1 and PC2 from the PCA of the RNA-seq SNPs for each dataset with infection status selected as the variable of interest. The PC1 and PC2 covariates were included because the crossbred/multibreed nature of the animals in each data set study population should be incorporated in the DEA model design. Genes with raw expression counts ≥ 6 in at least 20% of samples were retained prior to the DEA. For the DEA, the null hypothesis was that the LFC between the bTB+ and bTB− groups, for the expression of a particular gene is exactly 0. To account for potential heteroscedasticity of LFCs, we implemented the approximate posterior estimation for generalised linear model coefficients (APEGLM) method ( Zhu et al . 2019 ) using the lfcShrink function. Genes with a Benjamini-Hochberg (BH) false discovery rate (FDR) adjusted P -value ( Benjamini & Hochberg 1995 ) (FDR- P adj. ) 0 were considered significantly DEGs. For the OGR25-BTB and MCL2014-BTB datasets, a straightforward comparison between bTB− and bTB+ cattle was conducted. For the MCL21-BTB dataset, five comparisons were performed to compare the bTB+ animals at each time point (+1 wpi, +2 wpi, +6 wpi, +10 wpi and +12 wpi) with the bTB− animals pre-infection (−1 wpi). For the WIA20-BTB dataset, the bTB+ animals at each time point (+4 wpi, +10 wpi) were compared to the combined n = 13 bTB+ and bTB− animals at 0 wpi; the “bTB+” animals at 0 wpi were sampled prior to experimental infection with M. bovis and are therefore designated as bTB−. The DE results were plotted as volcano plots using ggplot2 (v. 3.4.4) ( Wickham 2016 ) and the intersections of DEGs among studies were visualised using the UpSetR R package (v. 1.4.0) ( Conway et al . 2017 ) Functional enrichment analysis of differentially expressed genes A functional enrichment analysis was conducted for DEGs inferred from each bTB dataset using the g:GOSt tool in g:Profiler (v. 0.2.2) ( Raudvere et al . 2019 ). The organism selected was Bos taurus and genes exhibiting increased, and separately decreased expression, were used to generate ordered input query lists based on the DEA FDR- P adj. values. Note: For the analysis of time series data, we aggregated results from each time point, filtered for DEGs that had a consistent direction of effect (either LFC > 0 or LFC < 0) across time points and ordered the retained genes based on the DEA FDR- P adj. values. In addition, we followed best practice recommendations to account for tissue-specific sampling biases in gene set overrepresentation and functional enrichment analyses ( Timmons et al . 2015 ). Consequently, for analysis of DEGs, the background set consisted of all expressed genes that were tested for differential expression across the bTB datasets. For analyses of our query DEG sets, we selected the gene ontology biological process (GO:BP) and cellular component (GO:CC) categories from the GO knowledgebase ( Ashburner et al . 2000 ) in addition to the Kyoto Encyclopaedia of Genes and Genomes (KEGG) ( Kanehisa et al . 2025 ) and Reactome ( Milacic et al . 2024 ) resources. To identify significantly enriched/overrepresented pathways, a BH-FDR multiple testing correction was applied (FDR- P adj. < 0.05). Generation of training and testing sets We divided the complete bTB dataset ( n = 230) into a training set (70%; n = 162: n = 69 (bTB−) and n = 93 (bTB+)) and a testing set (30%; n = 68: n = 34 (bTB−) and n = 34 (bTB+)), ensuring that the repeated measures of RNA-seq data derived from animals involved in the time series experiments were included in only one of the two sets to avoid overlap. For the training of ML models and evaluation of hyperparameters in the training set, we further divided the training set into five folds: Fold 1 ( n = 35: n = 14 (bTB−) and n = 21 (bTB+)); Fold 2 ( n = 36: n = 17 (bTB−) and n = 19 (bTB+)); Fold 3 ( n = 32: n = 16 (bTB−) and n = 16 (bTB+)); Fold 4 ( n = 34: n = 10 (bTB−) and n = 24 (bTB+)); and Fold 5 ( n = 25: n = 12 (bTB−) and n = 13 (bTB+)). This procedure was implemented in a similar fashion to the broad train/test splitting procedure, ensuring that repeated measures of RNA-seq data from animals involved in the time series experiments were not present across folds. Differential expression analysis on bTB training data Within the entire training dataset, a DEA was conducted between the reactor (bTB+) and control (bTB−) animal groups using DESeq2 (v. 1.40.2) accounting for age in months binned into four categories (1–6 months, 6–12 months, 12–24 months, and > 24 months); study/batch (five levels: two for the OGR25-BTB dataset and one each for the MCL14-BTB, WIA20-BTB, and the MCL21-BTB datasets, respectively); and PC1 and PC2 derived from the PCA of the pruned RNA-seq called SNPs across all samples in the training set as covariates, with condition (bTB+ or bTB−) set as the variable of interest. For this DEA, the null hypothesis was that the LFC between the bTB+ and bTB− groups, for the expression of a particular gene, is exactly 0. Similar to the intra-dataset DE analysis, to account for potential heteroscedasticity of LFCs, we implemented the approximate posterior estimation for generalised linear model coefficients (APEGLM) method using the lfcShrink function. Genes with an FDR- P adj. 0 were considered to be significant DEGs. For the identification of features for inclusion in the development of ML classifiers, we post-hoc filtered DEGs and retained those with a mean median-of-ratios normalised expression value ( Love et al . 2014 ) > 100 ( Chang et al . 2024 ). Training of integrated classifiers on bTB training data Machine learning (ML) model training and hyperparameter grid search were performed using R (v. 4.3.2) ( R Core Team 2025 ). All ML classification models were generated using the variance-stabilised transformed (VST) counts for the post-hoc filtered DEGs obtained using the varianceStabilizingTransformation function from the DESeq2 R package as input in the caret R package (v. 6.0-90) ( Kuhn 2008 ) with results evaluated using the pROC package (v. 1.18.4) ( Robin et al . 2011 ) and visualised using ggplot2 (v. 3.4.4). The classification models and hyperparameters evaluated in the training set were: generalised logistic regression (GLM; no hyperparameters); generalised logistic regression with ridge feature selection (RIDGE; hyperparameters: α = 0, λ = [10 -5 –10 0 , for 10,000 values]); generalised logistic regression with least absolute shrinkage and selection operator feature selection (LASSO; hyperparameters: α = 1, λ = [range = 10 -5 –10 0 , for 10,000 values]); generalised logistic regression with elastic net feature selection (ENET; hyperparameters: α = [range = 0.05−0.95, in steps of 0.05], λ = [range = 10 -5 –10 0 , for 10,000 values]); random forest (RF; hyperparameters: mtry = [range = 2−35, in steps of 5], min.node.size = [range = 2−10, in steps of 1], num.trees = [501, 751, 1001, 1501, 3001, 5001], splitrule = “gini”]); random forest with the ExtraTrees algorithm (RF-ET; hyperparameters: mtry = [range = 5−100, in steps of 5], min.node.size = [range = 2−10, in steps of 1], num.trees = [501, 751, 1001, 1501, 3001, 5001], splitrule = “extratress”]); naïve Bayes (NB; hyperparameters: laplace = [in steps of 0 – 3 in steps of 0.1], adjust = [10 -1 −10 1 , 0.1], usekernel = TRUE); and a multi-layered perceptron (MLP; hyperparameters: size = [2−10, 1], decay = [range = -10 -9 −10 1 , in steps of 10 n where n = [-9−1]). The confidence intervals of AUROC estimates were computed using the ci function from the pROC R package which performed stratified bootstrap resampling 2000 times ( Carpenter & Bithell 2000 ; Robin et al . 2011 ). Greedy forward search strategy and generation of a bTB infection z -score Previous studies have shown across a wide range of diseases, including for TB ( Sweeney et al . 2016 ; Chang et al . 2024 ), sepsis ( Sweeney et al . 2018 ), and dengue fever ( Robinson et al . 2019 ) that a greedy forward search (GFS) algorithm can be used to generate a disease-specific “infection z -score” ( Sweeney et al . 2015 ), which is derived from a panel of DEGs with high discriminatory power between cases and controls. The infection Z score is calculated as the mean of the VST normalised expression levels for positive DEGs (LFC > 0) minus the mean of the VST normalised expression levels for all negative DEGs (LFC < 0) ( Sweeney et al . 2015 ; Chang et al . 2024 ). In this study, the positive and negative DEGs, respectively, are derived from the comparison of bTB+ and bTB− cattle in the training set. The GFS algorithm starts with a single gene exhibiting the best discriminatory power represented as the average AUC (area under the curve) across the five folds used in the training set, and then at each subsequent step adds the gene with the best possible increase in the average AUC to the set of genes, until no further additions can increase the average AUC. Since there was a substantial amount of transcriptional heterogeneity present across bTB datasets, we wanted to generate a robust diagnostic signature set rather than aiming for extreme parsimony. Therefore, we ran the forward search a second time, such that once the first gene set had been identified, those genes were removed from the remaining pool and the forward search was run again. We subsequently generated a third combined list of 30 genes derived from the sets of genes identified in the first and second iterations, respectively. Assessment of classifier performance on bTB testing data In order to normalise the testing data, we first froze the VST estimates of the training data using the dispersionFunction function from DESeq2 and subsequently applied these frozen VST estimates to the testing data again using the dispersionFunction function. We then VST-normalised the testing data with the updated dispersion estimates using the varianceStabilizingTransformation function from DESeq2 setting the blind parameter to FALSE. Within the testing set, we then extracted the VST counts of DEGs used in the training set and selected the hyperparameter-tuned models and gene sets from the GFS strategy for evaluation. Evaluation of classifier performance on other infectious diseases We evaluated the performance of the hyperparameter-tuned models and inferred GFS gene sets using publicly available RNA-seq data for three other diseases, including Johne’s disease caused by M. avium spp. paratuberculosis (MAP) ( n = 3 MAP−, n = 11 MAP+) ( Alonso-Hearn et al . 2019 ) and bovine respiratory disease (BRD) caused by infection with either BoHV-1 ( n = 6 BoHV-1−, n = 12 BoHV-1+) ( O’Donoghue et al . 2023 ), or BRSV ( n = 6 BRSV−, n = 12 BRSV+) ( Johnston et al . 2021 ), respectively. We first evaluated the performance of each classifier in an intra-dataset manner, whereby we compared the discriminatory power of each classifier between control non-infected animals and cattle infected with each of the respective agents, as this may highlight common and divergent host peripheral transcriptional changes between M. bovis , MAP, BoHV-1, or BRSV. For each dataset, we followed the same procedure used for normalising the bTB testing dataset, which included applying the frozen VST estimates from the bTB training dataset to each external dataset via the dispersionFunction function from DESeq2, followed by obtaining the VST normalised counts using the varianceStabilizingTransformation function from DESeq2, setting the blind parameter to FALSE. We then extracted the VST counts of DEGs used in the bTB training set and selected the hyperparameter-tuned models and gene sets from the bTB training set GFS strategy for evaluation. Finally, we evaluated the specificity of each classifier by determining the discriminatory power between bTB+ animals ( n = 34) in the testing set and MAP+ ( n = 11), BRSV+ ( n = 12) and BoHV-1+ ( n = 12) animals. To do this we extracted and integrated all bTB+ animals from the testing set and all MAP+, BRSV+ and BoHV-1+ animals from the three external datasets and set the TB_status variable of all bTB+ animals to Infected and all other animals (MAP+, BRSV+ and BoHV-1+) to Control. We then normalised the expression matrix following the procedure implemented for the testing set and the external dataset evaluations. Declarations Ethics approval and consent to participate Not applicable. Consent for publication Not applicable. Availability of data and materials These raw and pre-processed RNA-seq datasets used in this study are available under the following accession numbers: 1) OGRA25-BTB ( GSE255724 ) ( O’Grady et al . 2025 ); 2) MCL14-BTB ( PRJNA257841 ) ( McLoughlin et al . 2014 ); 3) MCL21-BTB ( PRJEB27764 (bTB−) and PRJEB44470 (bTB+)) ( McLoughlin et al . 2021 ); 4) WIA20-BTB ( PRJNA600004 ) ( Wiarda et al . 2020 ); 5) ALO19-MAP ( PRJNA565369 ) ( Alonso-Hearn et al . 2019 ); 6) JOH21-BRD ( GSE152959 ) ( Johnston et al . 2021 ); and 7) ODO23-BRD ( GSE199108 ) ( O’Donoghue et al . 2023 ). Imputed and filtered WGS data from O’Grady et al . (2025) is available at Zenodo (( https://zenodo.org/records/13752056 ). Raw genotype information derived from the RNA-seq variant call for bTB studies described here is available at Zenodo ( https://doi.org/10.5281/zenodo.16887104 ). The computer code and scripts used in this study, in addition to the trained models are available at the following GitHub link: https://github.com/jfogrady1/ML4Tb Competing interests The authors declare that they have no competing interests. Funding J.F.O.’G. was supported by Research Ireland through the Research Ireland Centre for Research Training in Genomics Data Science (grant no. 18/CRT/6214). This study was also supported by Science Foundation Ireland (SFI) Investigator Programme Awards to D.E.M. and S.V.G. (grant nos. SFI/08/IN.1/B2038 and SFI/15/IA/3154); the University College Dublin— University of Edinburgh Strategic Partnership in One Health awarded to D.E.M., S.V.G., J.G.D.P., and E.L.C.; the Department of Agriculture, Food and the Marine (DAFM) project awards to D.E.M. (TARGET-TB: grant no. 17/RD/US-ROI/52; TB-ORNOT-TB: grant no. 2023RP902); the European Network on Livestock Phenomics (EU-LI-PHE) COST Action (CA22112); and a Fulbright-Teagasc Ireland-USA grant awarded to J.F.O.’G.. Note: Since the 1st of August 2024, Science Foundation Ireland (SFI) has been part of Taighde Éireann— Research Ireland ( www.researchireland.ie ). Authors’ contributions J.F.O.’G., K.G.M., E.G., I.C.G., S.V.G., C.S.G., and D.E.M. conceived and designed the study. D.E.M., C.S.G., I.C.G., S.V.G., J.G.D.P., E.L.C, H.P., and J.F.O.’G. acquired funding for the study. J.F.O.’G. performed all bioinformatics analyses and A.I., C.N.C., J.G.D.P., S.L.F.O.’D., V.R., H.P., J.A.B., G.P.M., A.K., and T.J.H., provided important scientific input. J.F.O.’G. wrote the first draft of the manuscript, and prepared all figures, tables, and supplementary information with input from D.E.M. and C.S.G. All authors read and approved the final manuscript. Acknowledgements We thank Michael McDonald for assistance with animal handling and blood collection. Funder Information Declared Taighde Éireann - Research Ireland , 18/CRT/6214 Science Foundation Ireland, https://ror.org/0271asj38 , SFI/08/IN.1/B2038 , SFI/15/IA/3154 Department of Agriculture Food and the Marine, https://ror.org/008gjgb19 , 17/RD/US-ROI/52 , 2023RP902 European Commission, https://ror.org/00k4n6c32 , CA22112 Fulbright Program, https://ror.org/014geb094 Footnotes https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE255724 https://zenodo.org/records/13752056 https://doi.org/10.5281/zenodo.16887104 References ↵ Afriyie-Asante A. , Dabla A. , Dagenais A. , Berton S. , Smyth R. & Sun J . ( 2021 ) Mycobacterium tuberculosis Exploits Focal Adhesion Kinase to Induce Necrotic Cell Death and Inhibit Reactive Oxygen Species Production . Frontiers in Immunology 12 , 742370 . OpenUrl PubMed ↵ Akter S. , Chauhan K.S. , Dunlap M.D. , Choreño-Parra J.A. , Lu L. , Esaulova E. , Zúñiga J. , Artyomov M.N. , Kaushal D. & Khader S.A. ( 2022 ) Mycobacterium tuberculosis infection drives a type I IFN signature in lung lymphocytes . Cell Reports 39 , 110983 . OpenUrl PubMed ↵ Allen A.R. , Skuce R.A. & Byrne A.W . ( 2018 ) Bovine Tuberculosis in Britain and Ireland - A Perfect Storm? the Confluence of Potential Ecological and Epidemiological Impediments to Controlling a Chronic Infectious Disease . Frontiers in Veterinary Science 5 , 109 . OpenUrl PubMed ↵ Alonso-Hearn M. , Canive M. , Blanco-Vazquez C. , Torremocha R. , Balseiro A. , Amado J. , Varela-Martinez E. , Ramos R. , Jugo B.M. & Casais R. ( 2019 ) RNA-Seq analysis of ileocecal valve and peripheral blood from Holstein cattle infected with Mycobacterium avium subsp. paratuberculosis revealed dysregulation of the CXCL8/IL8 signaling pathway . Scientific Reports 9 , 14845 . OpenUrl PubMed ↵ Altman D.G. & Bland J.M . ( 1994 ) Diagnostic tests 2: predictive values . BMJ 309 , 102 . OpenUrl FREE Full Text ↵ Andrews S. ( 2010 ) FastQC: a quality control tool for high throughput sequence data . URL https://www.bioinformatics.babraham.ac.uk/projects/fastqc/ . ↵ Ashburner M. , Ball C.A. , Blake J.A. , Botstein D. , Butler H. , Cherry J.M. , Davis A.P. , Dolinski K. , Dwight S.S. , Eppig J.T. , Harris M.A. , Hill D.P. , Issel-Tarver L. , Kasarskis A. , Lewis S. , Matese J.C. , Richardson J.E. , Ringwald M. , Rubin G.M. & Sherlock G . ( 2000 ) Gene ontology: tool for the unification of biology. The Gene Ontology Consortium . Nature Genetics 25 , 25 – 9 . OpenUrl CrossRef PubMed Web of Science ↵ Atkins P.J. & Robinson P.A . ( 2013 ) Bovine tuberculosis and badgers in Britain: relevance of the past . Epidemiology and Infection 141 , 1437 – 44 . OpenUrl CrossRef ↵ Benjamini Y. & Hochberg Y . ( 1995 ) Controlling the false discovery rate: a practical and powerful approach to multiple testing . Journal of the Royal Statistical Society Series B, Methodological 57 , 289 – 300 . OpenUrl CrossRef PubMed Web of Science ↵ Bennett R.M . ( 2017 ) The political economy of bovine tuberculosis in Great Britain . Revue Scientifique et Technique, Office International des Epizooties 36 , 105 – 14 . OpenUrl ↵ Bogunovic D. , Byun M. , Durfee L.A. , Abhyankar A. , Sanal O. , Mansouri D. , Salem S. , Radovanovic I. , Grant A.V. , Adimi P. , Mansouri N. , Okada S. , Bryant V.L. , Kong X.F. , Kreins A. , Velez M.M. , Boisson B. , Khalilzadeh S. , Ozcelik U. , Darazam I.A. , Schoggins J.W. , Rice C.M. , Al-Muhsen S. , Behr M. , Vogt G. , Puel A. , Bustamante J. , Gros P. , Huibregtse J.M. , Abel L. , Boisson-Dupuis S. & Casanova J.L . ( 2012 ) Mycobacterial disease and impaired IFN-γ immunity in humans with inherited ISG15 deficiency . Science 337 , 1684 – 8 . OpenUrl Abstract / FREE Full Text ↵ Bolger A.M. , Lohse M. & Usadel B . ( 2014 ) Trimmomatic: a flexible trimmer for Illumina sequence data . Bioinformatics 30 , 2114 – 20 . OpenUrl CrossRef PubMed Web of Science ↵ Broughan J.M. , Judge J. , Ely E. , Delahay R.J. , Wilson G. , Clifton-Hadley R.S. , Goodchild A.V. , Bishop H. , Parry J.E. & Downs S.H . ( 2016 ) A review of risk factors for bovine tuberculosis infection in cattle in the UK and Ireland . Epidemiology and Infection 144 , 2899 – 926 . OpenUrl CrossRef ↵ Buddle B.M. , de Lisle G.W. , Pfeffer A. & Aldwell F.E. ( 1995 ) Immunological responses and protection against Mycobacterium bovis in calves vaccinated with a low dose of BCG . Vaccine 13 , 1123 – 30 . OpenUrl CrossRef PubMed Web of Science ↵ Carpenter J. & Bithell J . ( 2000 ) Bootstrap confidence intervals: when, which, what? A practical guide for medical statisticians . Statistics in Medicine 19 , 1141 – 64 . OpenUrl CrossRef PubMed Web of Science ↵ Chai Q. , Wang L. , Liu C.H. & Ge B . ( 2020 ) New insights into the evasion of host innate immunity by Mycobacterium tuberculosis . Cellular and Molecular Immunology 17 , 901 – 13 . OpenUrl ↵ Chang A. , Loy C.J. , Eweis-LaBolle D. , Lenz J.S. , Steadman A. , Andgrama A. , Nhung N.V. , Yu C. , Worodria W. , Denkinger C.M. , Nahid P. , Cattamanchi A. & De Vlaminck I. ( 2024 ) Circulating cell-free RNA in blood as a host response biomarker for detection of tuberculosis . Nature Communications 15 , 4949 . OpenUrl PubMed ↵ Clarke A. , Byrne A.W. , Maher J. , Ryan E. , Farrell F. , McSweeney C. & Barrett D . ( 2022 ) Engaging with farmers to explore correlates of bovine tuberculosis risk in an internationally important heritage landscape: the Burren, in the west of Ireland . Frontiers in Veterinary Science 9 , 791661 . OpenUrl PubMed ↵ Clegg T.A. , Duignan A. , Whelan C. , Gormley E. , Good M. , Clarke J. , Toft N. & More S.J . ( 2011 ) Using latent class analysis to estimate the test characteristics of the γ-interferon test, the single intradermal comparative tuberculin test and a multiplex immunoassay under Irish conditions . Veterinary Microbiology 151 , 68 – 76 . OpenUrl CrossRef PubMed ↵ Conlan A.J. , McKinley T.J. , Karolemeas K. , Pollock E.B. , Goodchild A.V. , Mitchell A.P. , Birch C.P. , Clifton-Hadley R.S. & Wood J.L . ( 2012 ) Estimating the hidden burden of bovine tuberculosis in Great Britain . PLoS Computational Biology 8 , e1002730 . OpenUrl ↵ Conway J.R. , Lex A. & Gehlenborg N . ( 2017 ) UpSetR: an R package for the visualization of intersecting sets and their properties . Bioinformatics 33 , 2938 – 40 . OpenUrl CrossRef PubMed ↵ Cosivi O. , Grange J.M. , Daborn C.J. , Raviglione M.C. , Fujikura T. , Cousins D. , Robinson R.A. , Huchzermeyer H.F. , de Kantor I. & Meslin F.X. ( 1998 ) Zoonotic tuberculosis due to Mycobacterium bovis in developing countries . Emerging Infectious Diseases 4 , 59 – 70 . OpenUrl CrossRef PubMed Web of Science ↵ Danecek P. , Bonfield J.K. , Liddle J. , Marshall J. , Ohan V. , Pollard M.O. , Whitwham A. , Keane T. , McCarthy S.A. , Davies R.M. & Li H . ( 2021 ) Twelve years of SAMtools and BCFtools . Gigascience 10 . ↵ Das S. , Banerjee S. , Majumder S. , Chowdhury B.P. , Goswami A. , Halder K. , Chakraborty U. , Pal N.K. & Majumdar S . ( 2014 ) Immune subversion by Mycobacterium tuberculosis through CCR5 mediated signaling: involvement of IL-10 . PLoS ONE 9 , e92477 . OpenUrl CrossRef PubMed ↵ de la Rua-Domenech R. , Goodchild A.T. , Vordermeier H.M. , Hewinson R.G. , Christiansen K.H. & Clifton-Hadley R.S. ( 2006 ) Ante mortem diagnosis of tuberculosis in cattle: a review of the tuberculin tests, γ-interferon assay and other ancillary diagnostic techniques . Research in Veterinary Science 81 , 190 – 210 . OpenUrl CrossRef PubMed Web of Science ↵ Department for Environment Food & Rural Affairs ( 2025 ) Tuberculosis (TB) in cattle in Great Britian . ↵ Department of Agriculture Food and the Marine ( 2025 ) Quarter 2 2025 Bovine TB Disease Trends and Expenditure . (ed. by F.a.t.M. Department of Agriculture ). Goverment of Ireland , Dublin, Ireland . ↵ DePristo M.A. , Banks E. , Poplin R. , Garimella K.V. , Maguire J.R. , Hartl C. , Philippakis A.A. , del Angel G. , Rivas M.A. , Hanna M. , McKenna A. , Fennell T.J. , Kernytsky A.M. , Sivachenko A.Y. , Cibulskis K. , Gabriel S.B. , Altshuler D. & Daly M.J. ( 2011 ) A framework for variation discovery and genotyping using next-generation DNA sequencing data . Nature Genetics 43 , 491 – 8 . OpenUrl CrossRef PubMed Web of Science ↵ Dobin A. , Davis C.A. , Schlesinger F. , Drenkow J. , Zaleski C. , Jha S. , Batut P. , Chaisson M. & Gingeras T.R . ( 2013 ) STAR: ultrafast universal RNA-seq aligner . Bioinformatics 29 , 15 – 21 . OpenUrl CrossRef PubMed Web of Science ↵ European Food Safety Authority ( 2012 ) Scientific Opinion on the use of a gamma interferon test for the diagnosis of bovine tuberculosis . EFSA Journal 10 , 2975 . OpenUrl ↵ Ewels P. , Magnusson M. , Lundin S. & Kaller M . ( 2016 ) MultiQC: summarize analysis results for multiple tools and samples in a single report . Bioinformatics 32 , 3047 – 8 . OpenUrl CrossRef PubMed ↵ George B.M. , Duncan D. , Waller E. , Marriott E. , Payne M.C. , Withenshaw S. , Brouwer A. , Harris K.A. , May H.E. & Avigad R . ( 2025 ) Bovine TB infection status in cattle in Great Britain in 2023 . Veterinary Record 196 , e5696 . OpenUrl PubMed ↵ Goodchild A.V. , Downs S.H. , Upton P. , Wood J.L. & de la Rua-Domenech R. ( 2015 ) Specificity of the comparative skin test for bovine tuberculosis in Great Britain . Veterinary Record 177 , 258 . OpenUrl Abstract / FREE Full Text ↵ Gormley E. , Doyle M.B. , Fitzsimons T. , McGill K. & Collins J.D . ( 2006 ) Diagnosis of Mycobacterium bovis infection in cattle by use of the gamma-interferon (Bovigam) assay . Veterinary Microbiology 112 , 171 – 9 . OpenUrl CrossRef PubMed Web of Science ↵ Griffin J.M. , Williams D.H. , Kelly G.E. , Clegg T.A. , O’Boyle I. , Collins J.D. & More S.J . ( 2005 ) The impact of badger removal on the control of tuberculosis in cattle herds in Ireland . Preventive Veterinary Medicine 67 , 237 – 66 . OpenUrl CrossRef PubMed Web of Science ↵ Gupta R.K. , Turner C.T. , Venturini C. , Esmail H. , Rangaka M.X. , Copas A. , Lipman M. , Abubakar I. & Noursadeghi M . ( 2020 ) Concise whole blood transcriptional signatures for incipient tuberculosis: a systematic review and patient-level pooled meta-analysis . Lancet Respir Med 8 , 395 – 406 . OpenUrl PubMed ↵ Hall T.J. , Mullen M.P. , McHugo G.P. , Killick K.E. , Ring S.C. , Berry D.P. , Correia C.N. , Browne J.A. , Gordon S.V. & MacHugh D.E . ( 2021 ) Integrative genomics of the mammalian alveolar macrophage response to intracellular mycobacteria . BMC Genomics 22 , 343 . OpenUrl CrossRef PubMed ↵ Hall T.J. , Vernimmen D. , Browne J.A. , Mullen M.P. , Gordon S.V. , MacHugh D.E. & O’Doherty A.M . ( 2019 ) Alveolar macrophage chromatin is modified to orchestrate host response to Mycobacterium bovis infection . Frontiers in Genetics 10 , 1386 . OpenUrl PubMed ↵ Hayton A. , O’Brien A. , Adler A. , Cutler K. , Clarke J. , Shaw D.J. , Watt N.J. & Harkiss G.D . ( 2025 ) Diagnostic performance of the Enferplex bovine tuberculosis antibody test using bulk tank milk samples from dairy cattle . Journal of Dairy Science 108 , 5296 – 312 . OpenUrl PubMed ↵ Humblet M.F. , Boschiroli M.L. & Saegerman C . ( 2009 ) Classification of worldwide bovine tuberculosis risk factors in cattle: a stratified approach . Veterinary Research 40 , 50 . OpenUrl CrossRef PubMed ↵ Johnston D. , Earley B. , McCabe M.S. , Kim J. , Taylor J.F. , Lemon K. , Duffy C. , McMenamy M. , Cosby S.L. & Waters S.M . ( 2021 ) Messenger RNA biomarkers of bovine respiratory syncytial virus infection in the whole blood of dairy calves . Scientific Reports 11 , 9392 . OpenUrl PubMed ↵ Kanehisa M. , Furumichi M. , Sato Y. , Matsuura Y. & Ishiguro-Watanabe M . ( 2025 ) KEGG: biological systems database as a model of the real world . Nucleic Acids Research 53 , D672 – D7 . OpenUrl CrossRef PubMed ↵ Killick K.E. , Browne J.A. , Park S.D. , Magee D.A. , Martin I. , Meade K.G. , Gordon S.V. , Gormley E. , O’Farrelly C. , Hokamp K. & MacHugh D.E . ( 2011 ) Genome-wide transcriptional profiling of peripheral blood leukocytes from cattle infected with Mycobacterium bovis reveals suppression of host immune genes . BMC Genomics 12 , 611 . OpenUrl CrossRef PubMed ↵ Kuhn M . ( 2008 ) Building predictive models in R using the caret package . Journal of Statistical Software 28 , 1 – 26 . OpenUrl CrossRef PubMed ↵ Kulkarni M.M . ( 2011 ) Digital multiplexed gene expression analysis using the NanoString nCounter system . Curr Protoc Mol Biol Chapter 25 , Unit25B 10. ↵ Lahuerta-Marin A. , Milne M.G. , McNair J. , Skuce R.A. , McBride S.H. , Menzies F.D. , McDowell S.J.W. , Byrne A.W. , Handel I.G. & de C.B.B.M . ( 2018 ) Bayesian latent class estimation of sensitivity and specificity parameters of diagnostic tests for bovine tuberculosis in chronically infected herds in Northern Ireland . Veterinary Journal 238 , 15 – 21 . OpenUrl PubMed ↵ Li W. , Turner A. , Aggarwal P. , Matter A. , Storvick E. , Arnett D.K. & Broeckel U . ( 2015 ) Comprehensive evaluation of AmpliSeq transcriptome, a novel targeted whole transcriptome RNA sequencing methodology for global gene expression analysis . BMC Genomics 16 , 1069 . OpenUrl CrossRef PubMed ↵ Liao Y. , Smyth G.K. & Shi W . ( 2014 ) featureCounts: an efficient general purpose program for assigning sequence reads to genomic features . Bioinformatics 30 , 923 – 30 . OpenUrl CrossRef PubMed Web of Science ↵ Love M.I. , Huber W. & Anders S . ( 2014 ) Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2 . Genome Biology 15 , 550 . OpenUrl CrossRef PubMed ↵ Madden J.M. , O’Donovan J. , Casey-Bryars M. , Sweeney J. , Messam L.L. , McAloon C.G. , More S.J. , Kenny K. , Ryan E. & Gormley E . ( 2024 ) The impact of changing the cut-off threshold of the interferon-gamma (IFN-γ) assay for diagnosing bovine tuberculosis in Ireland . Preventive Veterinary Medicine 224 , 106129 . OpenUrl PubMed ↵ Mahajan S. , Saini A. , Chandra V. , Nanduri R. , Kalra R. , Bhagyaraj E. , Khatri N. & Gupta P . ( 2015 ) Nuclear receptor Nr4a2 promotes alternative polarization of macrophages and confers protection in sepsis . Journal of Biological Chemistry 290 , 18304 – 14 . OpenUrl Abstract / FREE Full Text ↵ Malone K.M. , Rue-Albrecht K. , Magee D.A. , Conlon K. , Schubert O.T. , Nalpas N.C. , Browne J.A. , Smyth A. , Gormley E. , Aebersold R. , MacHugh D.E. & Gordon S.V . ( 2018 ) Comparative ’omics analyses differentiate Mycobacterium tuberculosis and Mycobacterium bovis and reveal distinct macrophage responses to infection with the human and bovine tubercle bacilli . Microbial Genomics 4 . ↵ McCormack J . ( 2023 ) Horizon scanning: what next for bovine TB control in England? Irish Veterinary Journal 76 , 18 . OpenUrl PubMed ↵ McKenna A. , Hanna M. , Banks E. , Sivachenko A. , Cibulskis K. , Kernytsky A. , Garimella K. , Altshuler D. , Gabriel S. , Daly M. & DePristo M.A . ( 2010 ) The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data . Genome Research 20 , 1297 – 303 . OpenUrl Abstract / FREE Full Text ↵ McLoughlin K.E. , Correia C.N. , Browne J.A. , Magee D.A. , Nalpas N.C. , Rue-Albrecht K. , Whelan A.O. , Villarreal-Ramos B. , Vordermeier H.M. , Gormley E. , Gordon S.V. & MacHugh D.E . ( 2021 ) RNA-seq transcriptome analysis of peripheral blood From cattle infected with Mycobacterium bovis across an experimental time course . Frontiers in Veterinary Science 8 , 662002 . OpenUrl PubMed ↵ McLoughlin K.E. , Nalpas N.C. , Rue-Albrecht K. , Browne J.A. , Magee D.A. , Killick K.E. , Park S.D. , Hokamp K. , Meade K.G. , O’Farrelly C. , Gormley E. , Gordon S.V. & MacHugh D.E . ( 2014 ) RNA-seq transcriptional profiling of peripheral blood leukocytes from cattle infected with Mycobacterium bovis . Frontiers in Immunology 5 , 396 . OpenUrl PubMed ↵ Meade K.G. , Gormley E. , Doyle M.B. , Fitzsimons T. , O’Farrelly C. , Costello E. , Keane J. , Zhao Y. & MacHugh D.E . ( 2007 ) Innate gene repression associated with Mycobacterium bovis infection in cattle: toward a gene signature of disease . BMC Genomics 8 , 400 . OpenUrl CrossRef PubMed ↵ Michel A.L. , Geoghegan C. , Hlokwe T. , Raseleka K. , Getz W.M. & Marcotty T . ( 2015 ) Longevity of Mycobacterium bovis in raw and traditional souring milk as a function of storage temperature and dose . PLoS ONE 10 , e0129926 . OpenUrl PubMed ↵ Milacic M. , Beavers D. , Conley P. , Gong C. , Gillespie M. , Griss J. , Haw R. , Jassal B. , Matthews L. , May B. , Petryszak R. , Ragueneau E. , Rothfels K. , Sevilla C. , Shamovsky V. , Stephan R. , Tiwari K. , Varusai T. , Weiser J. , Wright A. , Wu G. , Stein L. , Hermjakob H. & D’Eustachio P . ( 2024 ) The Reactome Pathway Knowledgebase 2024 . Nucleic Acids Research 52 , D672 – D8 . OpenUrl CrossRef PubMed ↵ Monaghan T.F. , Rahman S.N. , Agudelo C.W. , Wein A.J. , Lazar J.M. , Everaert K. & Dmochowski R.R . ( 2021 ) Foundational statistical principles in medical research: sensitivity, specificity, positive predictive value, and negative predictive value . Medicina (Kaunas ) 57 . ↵ More S.J . ( 2024 ) Eradication of bovine tuberculosis in Ireland: is it a case of now or never? Irish Veterinary Journal 77 , 22 . OpenUrl PubMed ↵ Nalpas N.C. , Magee D.A. , Conlon K.M. , Browne J.A. , Healy C. , McLoughlin K.E. , Rue-Albrecht K. , McGettigan P.A. , Killick K.E. , Gormley E. , Gordon S.V. & MacHugh D.E . ( 2015 ) RNA sequencing provides exquisite insight into the manipulation of the alveolar macrophage by tubercle bacilli . Scientific Reports 5 , 13629 . OpenUrl PubMed ↵ Nunez-Garcia J. , Downs S.H. , Parry J.E. , Abernethy D.A. , Broughan J.M. , Cameron A.R. , Cook A.J. , de la Rua-Domenech R. , Goodchild A.V. , Gunn J. , More S.J. , Rhodes S. , Rolfe S. , Sharp M. , Upton P.A. , Vordermeier H.M. , Watson E. , Welsh M. , Whelan A.O. , Woolliams J.A. , Clifton-Hadley R.S. & Greiner M. ( 2018 ) Meta-analyses of the sensitivity and specificity of ante-mortem and post-mortem diagnostic tests for bovine tuberculosis in the UK and Ireland . Preventive Veterinary Medicine 153 , 94 – 107 . OpenUrl CrossRef PubMed ↵ O’Donoghue S. , Earley B. , Johnston D. , McCabe M.S. , Kim J.W. , Taylor J.F. , Duffy C. , Lemon K. , McMenamy M. , Cosby S.L. , Morris D.W. & Waters S.M . ( 2023 ) Whole blood transcriptome analysis in dairy calves experimentally challenged with bovine herpesvirus 1 (BoHV-1) and comparison to a bovine respiratory syncytial virus (BRSV) challenge . Frontiers in Genetics 14 , 1092877 . OpenUrl PubMed ↵ O’Grady J.F. , McHugo G.P. , Ward J.A. , Hall T.J. , Faherty O’Donnell S.L. , Correia C.N. , Browne J.A. , McDonald M. , Gormley E. , Riggio V. , Prendergast J.G.D. , Clark E.L. , Pausch H. , Meade K.G. , Gormley I.C. , Gordon S.V. & MacHugh D.E. ( 2025 ) Integrative genomics sheds light on the immunogenetics of tuberculosis in cattle . Communications Biology 8 , 479 . OpenUrl PubMed ↵ Olea-Popelka F. , Muwonge A. , Perera A. , Dean A.S. , Mumford E. , Erlacher-Vindel E. , Forcella S. , Silk B.J. , Ditiu L. , El Idrissi A. , Raviglione M. , Cosivi O. , LoBue P. & Fujiwara P.I . ( 2017 ) Zoonotic tuberculosis in human beings caused by Mycobacterium bovis —a call for action . Lancet Infectious Diseases 17 , e21 – e5 . OpenUrl CrossRef PubMed ↵ Petit III R.A. , Hall M.B. , Tonkin-Hill G. , Zhu J. & Read T.D. ( 2025 ) fastq-dl: efficiently download FASTQ files from SRA or ENA repositories . URL https://github.com/rpetit3/fastq-dlER . ↵ Portillo-Gómez L. & Sosa-Iglesias E.G . ( 2011 ) Molecular identification of Mycobacterium bovis and the importance of zoonotic tuberculosis in Mexican patients . Int J Tuberc Lung Dis 15 , 1409 – 14 . OpenUrl CrossRef PubMed ↵ R Core Team ( 2025 ) R: A language and environment for statistical computing . R Foundation for Statistical Computing , Vienna, Austria . ↵ Ranjbar S. , Haridas V. , Jasenosky L.D. , Falvo J.V. & Goldfeld A.E . ( 2015 ) A role for IFITM proteins in restriction of Mycobacterium tuberculosis infection . Cell Reports 13 , 874 – 83 . OpenUrl PubMed ↵ Raudvere U. , Kolberg L. , Kuzmin I. , Arak T. , Adler P. , Peterson H. & Vilo J . ( 2019 ) g:Profiler: a web server for functional enrichment analysis and conversions of gene lists (2019 update) . Nucleic Acids Research 47 , W191 – W8 . OpenUrl CrossRef PubMed ↵ Rhodes S.G. , Buddle B.M. , Hewinson R.G. & Vordermeier H.M . ( 2000 ) Bovine tuberculosis: immune responses in the peripheral blood and at the site of active disease . Immunology 99 , 195 – 202 . OpenUrl CrossRef PubMed Web of Science ↵ Ring S.C. , Purfield D.C. , Good M. , Breslin P. , Ryan E. , Blom A. , Evans R.D. , Doherty M.L. , Bradley D.G. & Berry D.P . ( 2019 ) Variance components for bovine tuberculosis infection and multi-breed genome-wide association analysis using imputed whole genome sequence data . PLoS ONE 14 , e0212067 . OpenUrl CrossRef PubMed ↵ Robin X. , Turck N. , Hainard A. , Tiberti N. , Lisacek F. , Sanchez J.C. & Müller M . ( 2011 ) pROC: an open-source package for R and S+ to analyze and compare ROC curves . BMC Bioinformatics 12 , 77 . OpenUrl CrossRef PubMed ↵ Robinson M. , Sweeney T.E. , Barouch-Bentov R. , Sahoo M.K. , Kalesinskas L. , Vallania F. , Sanz A.M. , Ortiz-Lasso E. , Albornoz L.L. , Rosso F. , Montoya J.G. , Pinsky B.A. , Khatri P. & Einav S . ( 2019 ) A 20-gene set predictive of progression to severe dengue . Cell Reports 26 , 1104 – 11 e4. OpenUrl PubMed ↵ Rosen B.D. , Bickhart D.M. , Schnabel R.D. , Koren S. , Elsik C.G. , Tseng E. , Rowan T.N. , Low W.Y. , Zimin A. , Couldrey C. , Hall R. , Li W. , Rhie A. , Ghurye J. , McKay S.D. , Thibaud-Nissen F. , Hoffman J. , Murdoch B.M. , Snelling W.M. , McDaneld T.G. , Hammond J.A. , Schwartz J.C. , Nandolo W. , Hagen D.E. , Dreischer C. , Schultheiss S.J. , Schroeder S.G. , Phillippy A.M. , Cole J.B. , Van Tassell C.P. , Liu G. , Smith T.P.L. & Medrano J.F. ( 2020 ) De novo assembly of the cattle reference genome with single-molecule sequencing . Gigascience 9 . ↵ Ryan E. , Breslin P. , O’Keeffe J. , Byrne A.W. , Wrigley K. & Barrett D . ( 2023 ) The Irish bTB eradication programme: combining stakeholder engagement and research-driven policy to tackle bovine tuberculosis . Irish Veterinary Journal 76 , 32 . OpenUrl PubMed ↵ Sabuz S.H. , Jahan I. , Debnath B.K. , Islam M.M. & Islam M.S . ( 2025 ) Zoonotic Tuberculosis and Dairy Products: Comprehensive Meta-Analysis of Prevalence and Public Health Implications . Vet Med Sci 11 , e70556 . OpenUrl ↵ Sánchez D. , Rojas M. , Hernández I. , Radzioch D. , García L.F. & Barrera L.F . ( 2010 ) Role of TLR2- and TLR4-mediated signaling in Mycobacterium tuberculosis -induced macrophage death . Cellular Immunology 260 , 128 – 36 . OpenUrl CrossRef PubMed ↵ Servellita V. , Bouquet J. , Rebman A. , Yang T. , Samayoa E. , Miller S. , Stone M. , Lanteri M. , Busch M. , Tang P. , Morshed M. , Soloski M.J. , Aucott J. & Chiu C.Y . ( 2022 ) A diagnostic classifier for gene expression-based identification of early Lyme disease . Commun Med (Lond ) 2 , 92 . OpenUrl PubMed ↵ Simmons J.D. , Dill-McFarland K.A. , Stein C.M. , Van P.T. , Chihota V. , Ntshiqa T. , Maenetje P. , Peterson G.J. , Benchek P. , Nsereko M. , Velen K. , Fielding K.L. , Grant A.D. , Gottardo R. , Mayanja-Kizza H. , Wallis R.S. , Churchyard G. , Boom W.H. & Hawn T.R. ( 2022 ) Monocyte transcriptional responses to Mycobacterium tuberculosis associate with resistance to tuberculin skin test and interferon gamma release assay conversion . mSphere 7 , e0015922 . OpenUrl PubMed ↵ Singhania A. , Verma R. , Graham C.M. , Lee J. , Tran T. , Richardson M. , Lecine P. , Leissner P. , Berry M.P.R. , Wilkinson R.J. , Kaiser K. , Rodrigue M. , Woltmann G. , Haldar P. & O’Garra A . ( 2018 ) A modular transcriptional signature identifies phenotypic heterogeneity of human tuberculosis infection . Nature Communications 9 , 2308 . OpenUrl PubMed ↵ Subuddhi A. , Kumar M. , Majumder D. , Sarkar A. , Ghosh Z. , Vasudevan M. , Kundu M. & Basu J . ( 2020 ) Unraveling the role of H3K4 trimethylation and lncRNA HOTAIR in SATB1 and DUSP4-dependent survival of virulent Mycobacterium tuberculosis in macrophages . Tuberculosis (Edinb ) 120 , 101897 . OpenUrl PubMed ↵ Sun M. , Phan J.M. , Kieswetter N.S. , Huang H. , Yu K.K.Q. , Smith M.T. , Liu Y.E. , Wang C. , Gupta S. , Obermoser G. , Maecker H.T. , Krishnan A. , Suresh S. , Gupta N. , Rieck M. , Acs P. , Ghanizada M. , Chiou S.H. , Khatri P. , Boom W.H. , Hawn T.R. , Stein C.M. , Mayanja-Kizza H. , Davis M.M. & Seshadri C . ( 2024 ) Specific CD4 + T cell phenotypes associate with bacterial control in people who ‘resist’ infection with Mycobacterium tuberculosis . Nature Immunology 25 , 1411 – 21 . OpenUrl CrossRef PubMed ↵ Sweeney T.E. , Azad T.D. , Donato M. , Haynes W.A. , Perumal T.M. , Henao R. , Bermejo-Martin J.F. , Almansa R. , Tamayo E. , Howrylak J.A. , Choi A. , Parnell G.P. , Tang B. , Nichols M. , Woods C.W. , Ginsburg G.S. , Kingsmore S.F. , Omberg L. , Mangravite L.M. , Wong H.R. , Tsalik E.L. , Langley R.J. & Khatri P . ( 2018 ) Unsupervised analysis of transcriptomics in bacterial sepsis across multiple datasets reveals three robust clusters . Critical Care Medicine 46 , 915 – 25 . OpenUrl CrossRef PubMed ↵ Sweeney T.E. , Braviak L. , Tato C.M. & Khatri P . ( 2016 ) Genome-wide expression for diagnosis of pulmonary tuberculosis: a multicohort analysis . Lancet Respir Med 4 , 213 – 24 . OpenUrl PubMed ↵ Sweeney T.E. , Shidham A. , Wong H.R. & Khatri P . ( 2015 ) A comprehensive time-course-based multicohort analysis of sepsis and sterile inflammation reveals a robust diagnostic gene set . Science Translational Medicine 7 , 287r a71 . OpenUrl ↵ Tabone O. , Verma R. , Singhania A. , Chakravarty P. , Branchett W.J. , Graham C.M. , Lee J. , Trang T. , Reynier F. , Leissner P. , Kaiser K. , Rodrigue M. , Woltmann G. , Haldar P. & O’Garra A . ( 2021 ) Blood transcriptomics reveal the evolution and resolution of the immune response in tuberculosis . Journal of Experimental Medicine 218 . ↵ Teng J. , Gao Y. , Yin H. , Bai Z. , Liu S. , Zeng H. , Pig GTEx Consortium , Bai L. , Cai Z. , Zhao B. , Li X. , Xu Z. , Lin Q. , Pan Z. , Yang W. , Yu X. , Guan D. , Hou Y. , Keel B.N. , Rohrer G.A. , Lindholm-Perry A.K. , Oliver W.T. , Ballester M. , Crespo-Piazuelo D. , Quintanilla R. , Canela-Xandri O. , Rawlik K. , Xia C. , Yao Y. , Zhao Q. , Yao W. , Yang L. , Li H. , Zhang H. , Liao W. , Chen T. , Karlskov-Mortensen P. , Fredholm M. , Amills M. , Clop A. , Giuffra E. , Wu J. , Cai X. , Diao S. , Pan X. , Wei C. , Li J. , Cheng H. , Wang S. , Su G. , Sahana G. , Lund M.S. , Dekkers J.C.M. , Kramer L. , Tuggle C.K. , Corbett R. , Groenen M.A.M. , Madsen O. , Godia M. , Rocha D. , Charles M. , Li C.J. , Pausch H. , Hu X. , Frantz L. , Luo Y. , Lin L. , Zhou Z. , Zhang Z. , Chen Z. , Cui L. , Xiang R. , Shen X. , Li P. , Huang R. , Tang G. , Li M. , Zhao Y. , Yi G. , Tang Z. , Jiang J. , Zhao F. , Yuan X. , Liu X. , Chen Y. , Xu X. , Zhao S. , Zhao P. , Haley C. , Zhou H. , Wang Q. , Pan Y. , Ding X. , Ma L. , Li J. , Navarro P. , Zhang Q. , Li B. , Tenesa A. , Li K. , Liu G.E. , Zhang Z. & Fang L . ( 2024 ) A compendium of genetic regulatory effects across pig tissues . Nature Genetics 56 , 112 – 23 . OpenUrl CrossRef PubMed ↵ Timmons J.A. , Szkop K.J. & Gallagher I.J . ( 2015 ) Multiple sources of bias confound functional enrichment analysis of global -omics data . Genome Biology 16 , 186 . OpenUrl CrossRef PubMed ↵ Van der Auwera G.A. & O’Connor B.D . ( 2020 ) Genomics in the Cloud: Using Docker, GATK, and WDL in Terra (1st Edition) . O’Reilly Media . ↵ Vance R.E . ( 2025 ) Tuberculosis as an unconventional interferonopathy . Current Opinion in Immunology 92 , 102508 . OpenUrl CrossRef PubMed ↵ Waters W.R. , Palmer M.V. , Buddle B.M. & Vordermeier H.M . ( 2012 ) Bovine tuberculosis vaccine research: historical perspectives and recent advances . Vaccine 30 , 2611 – 22 . OpenUrl CrossRef PubMed Web of Science ↵ Wiarda J.E. , Boggiatto P.M. , Bayles D.O. , Waters W.R. , Thacker T.C. & Palmer M.V . ( 2020 ) Severity of bovine tuberculosis is associated with innate immune-biased transcriptional signatures of whole blood in early weeks after experimental Mycobacterium bovis infection . PLoS ONE 15 , e0239938 . OpenUrl CrossRef PubMed ↵ Wickham H . ( 2016 ) ggplot2: Elegant Graphics for Data Analysis . Springer-Verlag , New York . ↵ Wiseman J. , Cassidy J.P. & Gormley E . ( 2024 ) The problem that residual Mycobacterium bovis infection poses for the eradication of bovine tuberculosis . Veterinary Journal 308 , 106266 . OpenUrl PubMed ↵ Wood P.R. , Corner L.A. , Rothel J.S. , Baldock C. , Jones S.L. , Cousins D.B. , McCormick B.S. , Francis B.R. , Creeper J. & Tweddle N.E . ( 1991 ) Field comparison of the interferon-gamma assay and the intradermal tuberculin test for the diagnosis of bovine tuberculosis . Australian Veterinary Journal 68 , 286 – 90 . OpenUrl PubMed Web of Science ↵ World Health Organization ( 2020 ) Global Tuberculosis Report 2020 . World Health Organization , Geneva, Switzerland . ↵ Zheng X. , Levine D. , Shen J. , Gogarten S.M. , Laurie C. & Weir B.S . ( 2012 ) A high-performance computing toolset for relatedness and principal component analysis of SNP data . Bioinformatics 28 , 3326 – 8 . OpenUrl CrossRef PubMed Web of Science ↵ Zhu A. , Ibrahim J.G. & Love M.I . ( 2019 ) Heavy-tailed prior distributions for sequence count data: removing the noise and preserving large differences . Bioinformatics 35 , 2084 – 92 . OpenUrl CrossRef PubMed View the discussion thread. Back to top Previous Next Posted September 04, 2025. Download PDF Supplementary Material Data/Code Email Thank you for your interest in spreading the word about bioRxiv. NOTE: Your email address is requested solely to identify you as the sender of this article. Your Email * Your Name * Send To * Enter multiple addresses on separate lines or separate them with commas. You are going to email the following Accurate and robust classification of Mycobacterium bovis-infected cattle using peripheral blood RNA-seq data Message Subject (Your Name) has forwarded a page to you from bioRxiv Message Body (Your Name) thought you would like to see this page from the bioRxiv website. Your Personal Message CAPTCHA This question is for testing whether or not you are a human visitor and to prevent automated spam submissions. Share Accurate and robust classification of Mycobacterium bovis -infected cattle using peripheral blood RNA-seq data John F. O’Grady , Adriana Ivich , Gillian P. McHugo , Adnan Khan , Thomas J. Hall , Sarah L. F. O’Donnell , Carolina N. Correia , John A. Browne , Valentina Riggio , James G. D. Prendergast , Emily L. Clark , Hubert Pausch , Kieran G. Meade , Isobel C. Gormley , Eamonn Gormley , Stephen V. Gordon , Casey S. Greene , David E. MacHugh bioRxiv 2025.09.01.673540; doi: https://doi.org/10.1101/2025.09.01.673540 Share This Article: Copy Citation Tools Accurate and robust classification of Mycobacterium bovis -infected cattle using peripheral blood RNA-seq data John F. O’Grady , Adriana Ivich , Gillian P. McHugo , Adnan Khan , Thomas J. Hall , Sarah L. F. O’Donnell , Carolina N. Correia , John A. Browne , Valentina Riggio , James G. D. Prendergast , Emily L. Clark , Hubert Pausch , Kieran G. Meade , Isobel C. Gormley , Eamonn Gormley , Stephen V. Gordon , Casey S. Greene , David E. MacHugh bioRxiv 2025.09.01.673540; doi: https://doi.org/10.1101/2025.09.01.673540 Citation Manager Formats BibTeX Bookends EasyBib EndNote (tagged) EndNote 8 (xml) Medlars Mendeley Papers RefWorks Tagged Ref Manager RIS Zotero Tweet Widget Facebook Like Google Plus One Subject Area Genomics Subject Areas All Articles Animal Behavior and Cognition (7635) Biochemistry (17690) Bioengineering (13892) Bioinformatics (41936) Biophysics (21451) Cancer Biology (18588) Cell Biology (25499) Clinical Trials (138) Developmental Biology (13378) Ecology (19899) Epidemiology (2067) Evolutionary Biology (24320) Genetics (15609) Genomics (22506) Immunology (17736) Microbiology (40394) Molecular Biology (17181) Neuroscience (88603) Paleontology (666) Pathology (2832) Pharmacology and Toxicology (4824) Physiology (7641) Plant Biology (15152) Scientific Communication and Education (2045) Synthetic Biology (4294) Systems Biology (9825) Zoology (2271)