key: cord-0254930-h12ed8qv
authors: Dai, Manman; Feng, Min; Li, Ziwei; Chen, Weisan; Liao, Ming
title: Chicken peripheral blood lymphocyte response to ALV-J infection assessed by single-cell RNA sequencing
date: 2021-01-13
journal: bioRxiv
DOI: 10.1101/2021.01.12.426350
sha: 6de077f578cf738aedf9c7632786a7896cdc4188
doc_id: 254930
cord_uid: h12ed8qv

Chicken peripheral blood lymphocytes (PBLs) exhibit wide-ranging cell types, but current understanding of their subclasses, immune cell classification, and function is limited and incomplete. Previously, we found that viremia caused by avian leukosis virus subgroup J (ALV‐J) was eliminated by 21 days post infection (DPI), accompanied by increased CD8+ T cell ratio in PBLs and low antibody levels. Here we performed single-cell RNA sequencing (scRNA-seq) of PBLs in ALV-J infected and control chickens at 21 DPI to determine chicken PBL subsets and their specific molecular and cellular characteristics, before and after viral infection. Eight cell clusters and their potential marker genes were identified in chicken PBLs. T cell populations (clusters 6 and 7) had the strongest response to ALV-J infection at 21 DPI, based on detection of the largest number of differentially expressed genes (DEGs). T cell populations of clusters 6 and 7 could be further divided into four subsets: activated CD4+ T cells (cluster A0), Th1-like cells (cluster A2), Th2-like cells (cluster A1), and cytotoxic CD8+ T cells. Hallmark genes for each T cell subset response to viral infection were initially identified. Furthermore, pseudotime analysis results suggested that chicken CD4+ T cells could potentially differentiate into Th1-like and Th2-like cells. Moreover, ALV-J infection probably induced CD4+ T cell differentiation into Th1-like cells in which the most immune related DEGs were detected. With respect to the control group, ALV-J infection also had an obvious impact on PBL cell composition. B cells showed inconspicuous response and their numbers decreased in PBLs of the ALV-J infected chickens at 21 DPI. Percentages of cytotoxic Th1-like cells and CD8+ T cells were increased in the T cell population of PBLs from ALV-J infected chicken, which were potentially key mitigating factors against ALV-J infection. More importantly, our results provided a rich resource of gene expression profiles of chicken PBL subsets for a systems-level understanding of their function in homeostatic condition as well as in response to viral infection.

Adaptive immunity is known to play a vital protective role against avian viral 60 infections. However, studies so far have focused on the avian innate immune response, 61 while research on avian T cell or B cell immunity is still in its infancy [1] . For example, 2 (SARS-CoV-2), and influenza infection [7] [8] [9] . In chickens, the cell lineage 82 characteristics in some tissues have been identified using scRNA-seq, such as the 83 developing chicken limb [10], chicken skeletal muscle [11] , and embryonic chicken 84 gonad [12] . However, no reported study using scRNA-seq to determine chicken 85 immune cell subsets or lineages. Moreover, to our knowledge, scRNA-seq technology 86 has not even been applied to study chicken PBL responses to any viral infection. vaccines or drug treatments [13] . A potential vaccine for ALV-J has been reported to 90 induce significantly increased CD4 + and CD8 + T cell percentage as well as IL-4 and 91 IFN-γ levels in immunized chickens [14] . Unfortunately, there have been few studies 92 exploring the specific T cell functions against ALV. A full understanding of the ALV-93 specific cellular immune response in chickens is likely the premise for developing 94 effective vaccines. In our previous study, we found that ALV-J viremia was eliminated Neighbor Embedding (t-SNE) [19] or Uniform Manifold Approximation and Projection 157 (UMAP) [20] in Seurat were used to visualize and explore these datasets.

(3) Differentially expressed gene (DEG) (up-regulation) analysis per cluster 159 We used the Wilcoxon rank sum test [21] to identify differential expression for a single 160 cluster, compared to all other cells. We identified differentially expressed genes (DEGs) Therefore, we defined cluster 7 as Th2-like γδT cell (CD3 + KK34 + DRD4 + ; Fig. 2C and 229 D). In addition, we found that a few cytotoxic CD8 + T cells (CD3 + CD8 + GNLY + ) mixed 230 in clusters 6 and 7 (Fig. 2D) . Therefore, we planned to regroup clusters 6 and 7 for a 231 more detailed display of data in the following analysis. Additionally, cluster 8 should 232 be B cells according to the expression of the known gene markers, BCL11A and Bu-1 233 (ENSGALG00000015461) (Fig. 2C) . Interestingly, the top five genes in cluster 4 were 234 9 mainly interferon stimulating genes (ISGs), including RSAD2, TRIM25, OASL, and 235 IFIT5, which implied that cluster 4 also contained a type of important antiviral immune 236 cell (Fig. 2C) . Therefore, we defined cluster 4 as ISG expressing cells in PBL, which 237 needs to be further verified in future studies. Unfortunately, we were unable to define 238 clusters 0, 1, 2, 3, and 5 based on the few known chicken cell-type markers and the top 239 five genes expressed by cells in these clusters. (Fig. 1C, Fig. 3B , and 250 file S5). Next, an immune related DEG analysis was performed after regrouping 251 clusters 6 and 7 in the study outlined below.

Meanwhile, some DEGs in other clusters were analyzed. We found that an 253 important immune gene, interferon regulator 7 (IRF7), revealed up-regulation in 254 clusters 0, 1, 2, and 3 in response to ALV-J infection at 21 DPI ( Fig. 3B and file S5) . 255 We also found that the anti-apoptotic gene BCL2L10 [31] was up-regulated in clusters 256 0, 1, 3, 4, and 5 of PBLs from ALV-J infected chicken ( Fig. 3B and file S5) . Besides,

Pearson's correlation analysis indicated that clusters 0, 1, 2, 3, 4, and 5 were strongly 258 correlated between each other, but displayed very low correlation with clusters 6, 7, and The above-mentioned T cell populations of clusters 6 and 7, both in PBLs from the 265 ALV-J infected or control chicken, were further divided into four distinct clusters 266 (named as clusters A0 to A3) and visualized using UMAP (Fig. 4A-C) . Compared to 267 those in the control PBLs, the percentages of clusters A0 and A1 in the PBLs from 268 ALV-J infected chicken was decreased. Conversely, the percentage of clusters A2 and 269 A3 was significantly increased in PBLs from ALV-J infected chicken (Fig. 4C) . Next, 270 the top five genes expressed in clusters A0 to A3 were picked as potential marker genes 271 and they are shown in a dot plot (Fig. 4D, Fig. S2 , and File S6). Judging according to 272 the classical marker genes and our above analysis, we considered cluster A0 as 273 activated CD4 + T cells (CD3 + CD4 + IL7R + CD28 + ), cluster A3 as cytotoxic CD8 + T cells 274 (CD3 + CD8 + IL2RB + GNLY + ), and cluster A1 as Th2 like γδT cells (CD3 + KK34 + DRD4 + ) 275 ( Fig. 4D and E) . are mainly located at the late stage ( Fig. 5A and B) . These results suggest that clusters 283 A2 and A1 are probably differentiated from A0 and imply that cluster A2 may represent 293 Our previous results show that clusters 6 and 7 play a potentially important role in 294 the antiviral response during ALV-J infection based on their largely immune-related 295 DEG expression (Fig. 3A and file S5) . Meanwhile, the T cell populations in clusters 6 296 and 7 could be grouped into four distinct sub clusters (Fig. 4) . To further investigate 297 each T cell population response, we firstly calculated the DEGs between corresponding 298 T cell population of the PBLs from ALV-J infected and control chickens at 21 DPI using 299 Seurat. Strikingly, the highest numbers of total DEGs and the DEGs enriched in the GO 300 terms "response to stimulus" and "cell proliferation" were predominant in cluster A2 301 (Th1-like T cells; Fig. 6A and file S7), suggesting that cluster A2 represents likely the 302 vital effectors among PBLs of the ALV-J infected host.

Next, the specific DEGs of the four sub clusters enrichment in " response to 304 stimulus" (GO:0050896) were exhibited using a volcano plot (Fig. S4) . Moreover, a 305 total of 33 important immune-related genes were screened from clusters A0 to A3 and 306 their expression levels were presented as a heat map and a dot plot (Fig. 6B -C and file 307 S7). We also observed that most of these immune-related genes were highly expressed Finally, we conducted an interaction network analysis of the 33 candidate DEGs 322 based on the STRING database (Fig. 7) . The results implied that the 13 DEGs marked 323 as red were likely the more important hub genes. Specifically, 4 hub genes including 324 ITGA2, IL1RAP, NOX1 and CDK17 were up-regulated in cluster A2 (Th1-like cells, in T cell apoptosis in PBLs from ALV-J infected chicken (Fig. 1C, Fig. 3B and file S5) .

On the other hand, the percentage of clusters 0, 1, and 3 were obviously increased in 358 PBLs from ALV-J infected chicken compared to the control PBLs, which may be 359 associated with the up-regulation of the anti-apoptotic gene, BCL2L10 [31] . Competing financial interests 451 The authors declare that they have no conflicts of financial interest. 

Progress on chicken T cell immunity to viruses

Th1/Th2 polarization by viral 461 and helminth infection in birds

Systematic Identification of Host Immune Key Factors 463 Influencing Viral Infection in PBL of ALV

Reticuloendotheliosis Virus Inhibits the Immune 465

Response Acting on Lymphocytes from Peripheral Blood of Chicken

Comparative analysis of key immune protection 468 factors in H9N2 avian influenza viruses infected and immunized specific pathogen-free 469 chicken

Immune cell type 'fingerprints' at the basis of outcome diversity

A single-cell atlas of the 473 peripheral immune response in patients with severe COVID-19

Single-Cell Sequencing of 476

Peripheral Mononuclear Cells Reveals Distinct Immune Response Landscapes of COVID-477 19 and Influenza Patients

Predicting bacterial 479 infection outcomes using single cell RNA-sequencing analysis of human immune cells

A single-cell transcriptomic atlas of the 482 developing chicken limb

Identification of diverse cell populations in 484 skeletal muscles and biomarkers for intramuscular fat of chicken by single-cell RNA 485 sequencing

Differentiation Provided by Single-Cell Transcriptomics in the Chicken Embryo

Immunity to Avian Leukosis Virus: Where Are We Now and What 490 Should We Do?

Evaluation of a multi-epitope subunit vaccine against 492 avian leukosis virus subgroup J in chickens

Integrating single-cell transcriptomic 496 data across different conditions, technologies, and species

Identification of cell types from single-cell transcriptomes using a novel 498 clustering method

Dissection of AML Reveals Progenitor-like Cells that Correlate with Prognosis

for improved visualization of single-cell RNA-seq data

Dimensionality reduction for 505 visualizing single-cell data using UMAP

Multilineage 507 communication regulates human liver bud development from pluripotency

The Gene Ontology Resource: 20 years and still GOing strong. 510

The dynamics and 512 regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells

MAST: a flexible statistical 515 18 framework for assessing transcriptional changes and characterizing heterogeneity in 516 single-cell RNA sequencing data

STRING v10: protein-518 protein interaction networks, integrated over the tree of life

Cytoscape: a software 521 environment for integrated models of biomolecular interaction networks

Identification of a 524 novel cytokine-like transcript differentially expressed in avian gamma delta T cells

Innervation Augments T Helper 2-Type Allergic Inflammation in the Postnatal Lung

ITPR2 as a susceptibility 530 gene in sporadic arnyotrophic lateral sclerosis: a genorne-wide association study

Marek's disease: 533 Genetic regulation of gallid herpesvirus 2 infection and latency

Transgenic 4-1BBL-engineered vaccine 536 stimulates potent Gag-specific therapeutic and long-term immunity via increased priming 537 of CD44(+)CD62L(high) IL-7R(+) CTLs with up-and downregulation of anti-and pro-538 apoptosis genes

H-ras and N-ras are 540 dispensable for T-cell development and activation but critical for protective Th1 immunity

Overexpression of miR-210 promotes the 543 potential of cardiac stem cells against hypoxia

2020) T cell subset profile and inflammatory 546 cytokine properties in the gut-associated lymphoid tissues of chickens during infectious 547 bursal disease virus (IBDV) infection

Systematic identification of chicken type I, II and 549 III interferon-stimulated genes

IRF7 Is Involved in Both STING 551 and MAVS Mediating IFN-beta Signaling in IRF3-Lacking Chickens

CD4(+) Th1 cells promote CD8(+) Tc1 cell 554 survival, memory response, tumor localization and therapy by targeted delivery of 555 interleukin 2 via acquired pMHC I complexes

Clonal anergy of CD117(+)chB6(+) B cell 557 progenitors induced by avian leukosis virus subgroup

A high-throughput screen for genes 560 essential for PRRSV infection using a piggyBac-based system

Suppresses Influenza A Virus-Induced Lung Inflammation and Oxidative Stress

Gastroenteritis Virus (TGEV)-induced Mitochondrial Damage Via Targeting RB1

Upregulating Interleukin-1 Receptor Accessory Protein (IL1RAP), and Activating p38

Mitochondria-localised ZNFX1 functions as a 569 dsRNA sensor to initiate antiviral responses through MAVS

miR-146a-5p promotes replication of 571 infectious bronchitis virus by targeting IRAK2 and TNFRSF18

CD4+ regulatory T cells 574 require CTLA-4 for the maintenance of systemic tolerance

Foxo transcription 576 factors control regulatory T cell development and function

TGF-beta receptor 578 maintains CD4 T helper cell identity during chronic viral infections

Association of Marek's Disease 581 induced immunosuppression with activation of a novel regulatory T cells in chickens

Regulatory T cell properties of chicken CD4+CD25+ 584 cells

Single-cell profiling of cell populations in chicken peripheral blood 602 lymphocytes (PBLs) collected from avian leukosis virus subgroup-J (ALV-J)-603 infected and control chickens at 21 days post infection (DPI). (A) Uniform Manifold 604

Approximation and Projection (UMAP) displaying all cell clusters and their percentage 605 in the control PBLs. (B) UMAP displaying all cell clusters and their percentage in PBLs 606 from ALV-J infected chicken. (C) Cell distribution and percentage of each cluster in 607 PBLs from ALV-J

Analysis of cell types of each cluster in chicken PBLs

Stochastic Neighbor Embedding (t-SNE) projection representing the eight clusters of 614 cells identified in the chicken PBL pools (unified set of control and ALV-J infection 615 samples). (B) The statistics of genes involved in the GO terms "immune system 616 process

Expression levels of characteristic marker genes (CD3D, CD8A, GNLY, IL2RB, and 621 BLB2) in PBL clusters

Histogram showing all up-626 regulated (red) and down-regulated (green) DEGs, and the number of DEGs involved 627 in the GO terms "response to stimulus" and "cell proliferation" in PBLs from ALV-J 628 infected chicken compared to control PBLs within each cluster. (B) Dot plot 629 representing selected DEGs (IRF7, ITPR2, JARID2, and BCL2L10) expressed in eight 630 clusters within which cells from the ALV-J infected sample were compared with the 631 control sample. The intensity represents the expression level

UMAP displaying the cell clusters and their percentage of 639 the regrouped T cell populations (clusters 6 and 7) in the control PBLs. (B) UMAP 640 displaying the cell clusters and their percentage of the regrouped T cell population 641 (clusters 6 and 7) in PBLs from ALV-J infected chicken. (C) The regrouped T cell 642 distribution and percentage of each cluster in PBLs from ALV-J

Expression levels of characteristic marker genes (CD3D, CD28, and CD8A

Pseudotime analysis of clusters A0, A1, and A2. (A) Mapping of clusters

Pseudotime 651 trajectory calculated from all cells of clusters A0-A2 in the control and ALV-J infected 652 samples. Darker colored dots represent a shorter pseudotime and earlier differentiation 653 period. (C) Mapping of cells in control and ALV-J infected samples to the pseudotime 654 trajectory. (D) The cell states of pseudotime trajectory partitioned from all cells of 655 clusters A0-A2 in the control and ALV-J infected samples. (E) Cell percentage of each 656 state between the ALV-J

Histogram showing all the up-664 regulated (red) and down-regulated (green) DEGs, and the number of DEGs involved 665 in GO terms "response to stimulus" and "cell proliferation" in ALV-J stimulated cells 666 compared to control cells in clusters A0-A3. (B) Heatmap showing the normalized 667 expression (Z-score) of all immune-related DEGs in various cells of clusters A0-A3 668 within which cells from the ALV-J infected sample are compared with the control 669 sample. (C) Dot plot representing DEGs expressing in clusters A0-A3 within which 670 cells from the ALV-J

The interaction network analysis of the candidate 33 DEGs based on the 675 STRING database

Supplementary Fig.1 (Fig.S1). t-SNE displaying top 5 DEGs

Supplementary Fig.2 (Fig.S2). t-SNE displaying top 5 DEGs

Heat map of branch-dependent differential gene

Volcano plot of identified DEGs involved in the GO 689 term "response to stimulus" (GO:0050896) between ALV-J infected and control 690 samples in clusters

Supplementary file 1 (file S1). scRNA-seq data statistics

Up-regulated differentially expressed genes (DEGs) 695 in each cluster of PBLs

Up-regulated DEGs involved in "immune system 698 process

Top 5 DEGs of each chicken PBL cluster

DEGs found in each cluster between the ALV-J 704 infected and control samples

Top 5 DEGs in clusters A0-A3

DEGs in clusters A0-A3 between the ALV-J infected 709 and control samples