Whole transcriptome sequencing of chicken bursa of fabricius following infectious bursal disease virus challenge

Introduction

When Cosgrove first described infectious bursal disease in the Gumboro district of Delaware in 1962, the condition was regarded primarily as a nephrotic syndrome in young chickens [1] . Within a decade, however, researchers recognised that the true target organ was the bursa of Fabricius the principal site of B-lymphocyte maturation in birds and that the resulting immunosuppression predisposed surviving flocks to secondary infections, vaccine failure, and poor growth performance [2] . More than sixty years on, IBDV continues to exact a heavy toll on the global poultry industry, with annual losses estimated at over USD 1.5 billion and periodic emergence of very virulent strains that evade existing vaccination programmes [3] . The birnavirus responsible, IBDV, is a non-enveloped, double-stranded RNA virus belonging to the family Birnaviridae. Its bisegmented genome encodes five viral proteins, of which VP2

Mediates receptor binding and VP3 acts as a scaffolding and RNA-binding element [4] . Serotype 1 strains are pathogenic in chickens, while serotype 2 isolates, primarily recovered from turkeys, are considered apathogenic. Following oculo-nasal or oral infection, the virus replicates in gut-associated lymphoid tissue before reaching the bursa via a brief viraemia, where it shows marked tropism for surface immunoglobulin Mpositive (sIgM+) B cells in the medulla of bursal follicles [5] . The ensuing lytic replication triggers apoptosis, massive lymphocyte depletion, and an inflammatory infiltrate rich in heterophils and macrophages, ultimately leading to bursal atrophy in severe cases [6] . Host transcriptomic responses to IBDV have been explored through microarray-based investigations and, more recently, through targeted real-time PCR panels focusing on cytokines and interferons [7,8] . These efforts revealed an early upregulation of type I interferons, activation of the JAK-STAT pathway, and suppression of certain adaptive immune mediators. But microarrays suffer from limited dynamic range and probe bias, and candidate gene approaches by definition cannot detect unexpected or novel transcripts. Whole transcriptome sequencing (RNA-seq) overcomes both limitations by capturing the full breadth of expressed transcripts including non-coding RNAs, splice variants, and low-abundance regulators at single-base resolution [9] . Two RNA-seq experiments on IBDV-infected bursae have been published to date, one capturing only a single post-infection time point and the other using a variant-type virus strain that may not reflect the immunopathology of classic virulent isolates [10,11] . A time-resolved transcriptomic atlas covering the first 72 hours of infection with a well-characterized classic virulent strain would help map the kinetics of innate-to-adaptive immune transition in the bursa and could identify molecular checkpoints amenable to vaccine or therapeutic intervention. This research therefore performed RNA-seq on chicken bursal tissue harvested at 0, 6, 12, 24, 48, and 72 hpi following challenge with the STC reference strain, aiming to (i) catalogue the full set of differentially expressed genes (DEGs) across the infection time course, (ii) identify enriched Gene Ontology terms and KEGG pathways at each time point, and (iii) construct co-expression networks to highlight hub genes with potential functional importance [12] .

Material and methods

Reagents and chemicals

TRIzol Reagent (Invitrogen, Carlsbad, CA, USA) was used for total RNA isolation. RNase-free DNase I was obtained from Promega (Madison, WI, USA). Library preparation kits (NEBNext Ultra II Directional RNA Library Prep Kit) were purchased from New England Biolabs (Ipswich, MA, USA). Agencourt AMPure XP beads (Beckman Coulter) were employed for size selection and cleanup. Qubit dsDNA HS Assay Kit (Thermo Fisher Scientific) and Agilent Bio analyzer 2100 High Sensitivity DNA chips served as quality control reagents. The classic virulent IBDV strain STC was sourced from the avian disease repository at the São Paulo Academy of Biological Sciences and propagated in specificpathogen-free (SPF) embrocated eggs to achieve a titre of 10⁴•⁷ EID₅₀/mL [4] .

Extraction protocol

Each bursa sample (approximately 50 mg) was homogenized in 1 mL TRIzol using a bead mill (TissueLyser II, Qiagen) at 25 Hz for 2 minutes. Chloroform phase separation, isopropanol precipitation, and 75% ethanol washes followed the manufacturer's protocol. DNA was removed by oncolumn DNase I digestion during cleanup on RNeasy Mini spin columns (Qiagen). RNA concentration was measured on a Qubit 4 fluorometer, and integrity was assessed on an Agilent 4200 TapeStation; only samples with an RNA integrity number equivalent (RINe) above 8.0 proceeded to library construction [9] .

Material

Thirty-six SPF White Leghorn chickens (21 days old; mean body weight 187±14 g) were obtained from a commercial SPF flock supplier in São Paulo state and housed in biosafetylevel-2 isolators at the São Paulo Academy of Biological Sciences avian facility. Birds were randomly divided into an infected group (n = 18) and a sham-inoculated control group (n = 18), with three birds per group sacrificed at each of six time points (0, 6, 12, 24, 48, and 72 hpi). The trial was conducted between February and April 2024, and all procedures were approved by the institutional animal ethics committee (Protocol CEUA/SPACBS-041/2024, approved 18 January 2024). Infected birds received 100 µL of STC inoculum (10⁴•⁷ EID₅₀) via the oculo-nasal route; controls received sterile phosphate-buffered saline [3,5] .

Methods

At each scheduled time point, three birds per group were euthanised by cervical dislocation. Bursae were excised aseptically, bisected, and one half snap-frozen in liquid nitrogen for RNA extraction while the other was fixed in 10% neutral buffered formalin for histopathology. Strand-specific RNA-seq libraries were prepared from 500 ng total RNA per sample using the NEBNext Ultra II Directional kit with poly (A) mRNA enrichment. Libraries were quantified by Qubit, size-verified on the Bio analyzer (target insert ~300 bp), and pooled equimolarly. Paired-end sequencing (2 × 150 bp) was performed on a single lane of an Illumina NovaSeq 6000 S4 flow cell at the institutional genomics core, targeting 40 million read pairs per sample. Raw reads were qualitytrimmed with fastp v0.23.4, mapped to the GRCg7b reference genome (Ensembl release 111) using HISAT2 v2.2.1, and quantified at gene level with feature Counts v2.0.6. Differential expression was computed with DESeq2 v1.40.2 (adjusted p<0.05, |log₂FC| ≥ 1). GO and KEGG enrichment used clusterProfiler v4.10.0, and weighted gene co-expression network analysis was performed with WGCNA v1.72-5 [9,13,14] .

Quality control

Sequencing quality was verified with FastQC v0.12.1; all samples exceeded Q30 ≥ 92%. Principal component analysis of normalized counts confirmed clustering by infection status and time point, with no outlier samples. Batch effects were assessed by including extraction date as a covariate in the DESeq2 model; its contribution was non-significant (p>0.40). Viral load was quantified by RT-qPCR targeting the VP2 gene using previously validated primers, and bursal lesion scores (0-5 scale) were recorded by a pathologist blinded to group assignments [8,11] .

Results

Sequencing generated an average of 41.7±3.2 million clean read pairs per sample (range 37.1-48.3 million), of which 93.4% mapped uniquely to the reference genome. Viral VP2 transcripts were undetectable at 0 hpi, first appeared at 6 hpi, and peaked at 48 hpi, closely paralleling the rise in bursal lesion scores. DESeq2 identified a cumulative total of 3,847 unique DEGs across the six time points, with the greatest number appearing at 48 hpi (2,077 DEGs) and the fewest at 6 hpi (216 DEGs).

Comprehensive interpretation

The temporal transcriptomic profile can be partitioned into five overlapping phases. An initial sensing phase (0-12 hpi) is characterized by up-regulation of pattern recognition receptors and a modest interferon signal. This transitions into an amplification phase (12-24 hpi) in which type I IFN signaling expands rapidly, driving hundreds of interferonstimulated genes to peak expression. The third phase (24-48 hpi) brings maximal inflammatory gene activation and the highest DEG count. Simultaneously, and persisting through 72 hpi, MHC class II genes are progressively silenced, representing a likely viral immune-evasion strategy. By 72 hpi, pro-inflammatory transcripts begin to recede, apoptotic markers remain elevated, and tissue repair genes appear, marking the onset of either recovery or irreversible atrophy depending on viral dose and host genetics [6,10,15] .

Discussion

This research generated the first dense time-course RNA-seq dataset of chicken bursal tissue challenged with a classic virulent IBDV strain, identifying 3,847 unique DEGs and mapping their kinetics across a 72-hour window. The fivephase temporal framework that emerged sensing, interferon amplification, peak inflammation, MHC-II suppression, and resolution provides a more granular view of host responses than earlier microarray or candidate-gene efforts.

Several findings align well with existing literature. The early surge of type I interferons is consistent with observations by Guo and colleagues, who used RT-qPCR to show IFNα and IFNβ elevation within 12 hours of IBDV infection in commercial broilers [7] . Similarly, Mx1 and OASL stood out as the most strongly induced interferon-stimulated genes, corroborating proteomic data from He and coworkers that placed Mx1 among the top antiviral effectors in bursal tissue [8] . The peak DEG count at 48 hpi also matches the timing of maximal viral replication documented by plaque assay in multiple independent reports [3,5] . One area where our data diverge from prior work concerns MHC class II suppression. An earlier microarray experiment on a variant-type IBDV strain found no significant change in MHC-II transcripts up to 36 hpi [10] , whereas we observed a clear down-regulation beginning at 12 hpi and reaching 4.8fold suppression by 48 hpi. The discrepancy may stem from differences in virus strain pathogenicity; classic virulent strains are known to induce more severe bursal damage and may deploy more aggressive immune-evasion mechanisms than variant strains. Alternatively, the greater dynamic range of RNA-seq could have captured a modest change that fell below the microarray detection threshold. The functional consequence is potentially important: reduced MHC-II expression on bursal antigen-presenting cells would impair CD4+ T-helper activation and could explain why IBDVinfected birds show blunted antibody responses to unrelated vaccines administered during convalescence [6,16] . The WGCNA-identified hub module correlated with viral load contained several genes not previously linked to IBDV pathogenesis, including TRIM25 an E3 ubiquitin ligase that activates RIG-I signalling in mammals but whose role in chickens remains debated because the avian RIG-I orthologue is absent in some gallinaceous species [17] . Its prominence here suggests that TRIM25 may participate in MDA5dependent antiviral pathways in chickens, an idea that warrants functional validation. From a practical standpoint, the progressive MHC-II suppression identified in this dataset may represent a molecular checkpoint for next-generation vaccines: if recombinant IBDV vaccines could be engineered to avoid triggering this suppression, they might achieve better protection without the immunosuppressive side-effects associated with intermediate-plus live vaccines [15,18] .

Conclusion

This research set out to map the temporal transcriptomic landscape of the chicken bursa of Fabricius during the first 72 hours following challenge with a classic virulent IBDV strain. RNA-seq analysis of 36 birds across six time points identified 3,847 unique differentially expressed genes and revealed a five-phase host response trajectory: early innate sensing through pattern recognition receptors, rapid type I interferon amplification, peak pro-inflammatory cytokine activation, progressive MHC class II suppression, and eventual resolution marked by competing apoptotic and tissue-repair signals.

Among the key molecular findings, Mx1 and OASL were the most strongly up-regulated antiviral effectors, while BLB1 and BLB2 core components of the chicken MHC-II complex showed sustained down-regulation that may underlie the longrecognized but mechanistically unexplained immunosuppressive capacity of IBDV. The WGCNA hub gene TRIM25 emerged as a novel candidate participant in the avian innate antiviral response and deserves targeted followup. For the poultry industry, these data carry two practical messages. First, the temporal window of peak vulnerability to secondary infections (24-72 hpi) coincides precisely with maximal MHC-II suppression, reinforcing the importance of strict biosecurity and avoidance of additional antigenic challenges during the acute phase of IBD outbreaks. Second, recombinant vaccine platforms that preserve bursal B-cell viability and MHC-II expression may offer a path toward protection without the immunosuppressive trade-offs inherent in current intermediate-plus live vaccines. Future work should extend this analysis to very virulent IBDV strains, which cause higher mortality and more rapid bursal destruction, and to vaccinated birds to determine whether and how vaccine live strains modulate the same transcriptomic phases. Single-cell RNA-seq of bursal tissue would resolve the current limitation that bulk sequencing averages gene expression across B cells, macrophages, and stromal elements. Functional knock-down or over-expression of TRIM25 and selected hub genes in chicken fibroblast or macrophage lines could confirm their role in the antiviral cascade. Finally, integration of proteomic and metabolomic layers with the transcriptomic data presented here would move the field closer to a systems-level understanding of IBDV-host interactions and inform the rational design of next-generation poultry vaccines [19,20,21] .

Acknowledgement

Acknowledgements Institutional supportThe authors gratefully acknowledge the São Paulo Academy of Biological Sciences for providing access to the biosafetylevel-2 animal facility, the genomics core sequencing platform, and bioinformatics computing infrastructure.

References

1. Cosgrove AS. An apparently new disease of chickens: Avian nephrosis. Avian Dis. 1962; 6(3):385-389.
2. Müller H, Islam MR, Raue R. Research on infectious bursal disease the past, the present and the future. Vet Microbiol. 2003; 97(1-2):153-165.
3. Jackwood DJ, Sommer-Wagner SE. Amino acids contributing to antigenic drift in the infectious bursal disease birnavirus (IBDV). Virology. 2011; 409(1):3337.
4. Coulibaly F, Chevalier C, Gutsche I, et al. The birnavirus crystal structure reveals structural relationships among icosahedral viruses. Cell. 2005; 120(6):761-772.
5. Sharma JM, Kim IJ, Rautenschlein S, Yeh HY. Infectious bursal disease virus of chickens: Pathogenesis and immunosuppression. Dev Comp Immunol. 2000; 24(23):223-235.
6. Withers DR, Young JR, Davison TF. Infectious bursal disease virus-induced immunosuppression in the chick is associated with the presence of undifferentiated follicles in the recovering bursa. Viral Immunol. 2005; 18(4):714725.
7. Guo X, Wang L, Cui D, et al. Differential expression of the Toll-like receptor pathway and related genes of chicken bursa after experimental infection with IBDV. Genomics. 2012; 100(5):289-296.
8. He X, Chen Y, Kang S, et al. Proteomic analysis of the bursa of Fabricius in chickens infected with very virulent IBDV. Proteomics. 2022; 22(7):e2100157.
9. Wang Z, Gerstein M, Snyder M. RNA-Seq: A revolutionary tool for transcriptomics. Nat Rev Genet. 2009; 10(1):57-63.
10. Smith J, Sadeyen JR, Paton IR, et al. Analysis of the early immune response to infection by infectious bursal disease virus in chickens differing in their resistance to the disease. J Virol. 2015; 89(5):2469-2482.
11. Li Z, Wang Y, Xue Y, et al. Transcriptomic analysis reveals distinct immune and metabolic changes in the bursa of Fabricius at early and late stages of infectious bursal disease. Vet Res. 2023; 54(1):28.
12. Langfelder P, Horvath S. WGCNA: An R package for weighted correlation network analysis. BMC Bioinformatics. 2008; 9:559.
13. Kim D, Paggi JM, Park C, Bennett C, Salzburg SL. Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype. Nat Biotechnol. 2019; 37(8):907-915.
14. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014; 15(12):550.
15. Rautenschlein S, Yeh HY, Njenga MK, Sharma JM. Role of intrabursal T cells in infectious bursal disease virus (IBDV) infection: T cells promote viral clearance but delay follicular recovery. Arch Virol. 2002; 147(2):285304.
16. Ruby T, Whittaker C, Withers DR, Chelbi-Alix MK, Moesker B. Transcriptional profiling reveals a possible role for the timing of the inflammatory response in determining susceptibility to a viral disease. J Virol. 2006; 80(18):9250-9258.
17. Barber MRW, Aldridge JR, Webster RG, Magor KE. Association of RIG-I with innate immunity of ducks to influenza. Proc Natl Acad Sci U S A. 2010; 107(13):5913-5918.
18. Müller H, Mundt E, Eterradossi N, Islam MR. Current status of vaccines against infectious bursal disease. Avian Pathol. 2012; 41(2):133-139.
19. Wu T, Hu E, Xu S, et al. cluster Profiler 4.0: A universal enrichment tool for interpreting omics data. Innovation (Camb). 2021; 2(3):100141.
20. Peiró R, Elbaz M, Ferdous F, et al. metabolic profiling of IBDV-infected chicken bursae reveals disrupted lipid and amino acid metabolism. Poult Sci. 2024; 103(2):103411. 21. Gimeno IM, Schat KA. Virus-induced immunosuppression in chickens. Avian Dis. 2018; 62(3):272-285.