https://orcid.org/0000-0002-7745-2844
Figures
Abstract
The transcription factor Rora has been shown to be important for the development of ILC2 and the regulation of ILC3, macrophages and Treg cells. Here we investigate the role of Rora across CD4+ T cells in general, but with an emphasis on Th2 cells, both in vitro as well as in the context of several in vivo type 2 infection models. We dissect the function of Rora using overexpression and a CD4-conditional Rora-knockout mouse, as well as a RORA-reporter mouse. We establish the importance of Rora in CD4+ T cells for controlling lung inflammation induced by Nippostrongylus brasiliensis infection, and have measured the effect on downstream genes using RNA-seq. Using a systematic stimulation screen of CD4+ T cells, coupled with RNA-seq, we identify upstream regulators of Rora, most importantly IL-33 and CCL7. Our data suggest that Rora is a negative regulator of the immune system, possibly through several downstream pathways, and is under control of the local microenvironment.
Citation: Haim-Vilmovsky L, Henriksson J, Walker JA, Miao Z, Natan E, Kar G, et al. (2021) Mapping Rora expression in resting and activated CD4+ T cells. PLoS ONE 16(5): e0251233. https://doi.org/10.1371/journal.pone.0251233
Editor: Scott N. Mueller, The University of Melbourne, AUSTRALIA
Received: April 1, 2021; Accepted: April 22, 2021; Published: May 18, 2021
Copyright: © 2021 Haim-Vilmovsky et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: The sequencing data has been deposited at ArrayExpress (E-MTAB-7694, E-MTAB-8000, E-MTAB-8001, E-MTAB-8003), and the image data at Zenodo (https://doi.org/10.5281/zenodo.2554078). The R code used for the analysis is available on GitHub (https://github.com/mahogny/liorarora).
Funding: This study was supported by EMBO in the form of a grant awarded to LHV (ALTF 698-2012), Directorate-General for Research and Innovation in the form of a grant awarded to LHV (FP7-PEOPLE-2010-IEF, ThPLAST 274192), an EMBL Interdisciplinary Postdoctoral fellowship, supported by H2020 Marie Skłodowska Curie Actions and awarded to LHV, the Swedish Research Council in the form of a grant awarded to JH (#2016-06598), the European Research Council in the form of a grant awarded to SAT (ThDEFINE), FUV - Fondazione Umberto Veronesi in the form of funds awarded VP, UK Medical Research Council in the form of a grant awarded to JW and ANJM (U105178805), Wellcome Trust in the form of a grant awarded to JW and ANJM (100963/Z/13/Z) and a Single Cell Gene Expression Atlas grant awarded to ZM (108437/Z/15/Z), and the Open Targets Grant awarded to ZM (OTAR2067). Wellcome Sanger Institute core facilities are supported by a grant (WT206194). Microsoft Research Cambridge, a subsidiary of Microsoft, provided support in the form of a salary for JF. The specific roles of these authors are articulated in the ‘author contributions’ section. The funders had no role in the study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have read the journal’s policy and have the following competing interests: JF was an employee of Microsoft Research Cambridge at the time the study was conducted, but is no longer an employee of the company. EN is an employee of Aleph Labs, but was not employed by the company at the time the study was conducted. GK is an employee of AstraZeneca, but was not employed by the company at the time the study was conducted. This does not alter our adherence to PLOS ONE policies on sharing data and materials. There are no patents, products in development or marketed products associated with this research to declare.
Introduction
Type 2 response is driven by a broad range of stimuli including helminths, allergens, specific bacterial and viral infections, and endogenous host molecules [1]. This reaction is characterized by a cellular repertoire which include the innate lymphoid cells type 2 (ILC2) [2, 3] and adaptive T helper 2 (Th2) which are known to express transcription factor GATA-binding protein 3 (GATA-3). Those cells secrete distinct cytokines including interleukin-4 (IL-4), IL-5, and IL-13, driving eosinophil recruitment and immunoglobulin production.
Retinoic acid receptor-related orphan receptor alpha (Rora), is a nuclear receptor that functions as a ligand-dependent transcription factor [4–6], and is known to regulate functions in immunity in addition to other roles in development, circadian rhythm and metabolism [7]. Rora was found to play a critical role in the development of ILC2 [8, 9] and for cytokine production in ILC3 [10]. In the liver, it controls inflammation by promoting the macrophage M2 polarization [11]. Together with Rorc, it influences the development of Th17 [12]. Recently it was shown to have an important role in skin Tregs, where it is required for proper regulation of the immune response in atopic dermatitis. In this study, the deletion of Rora in Tregs resulted in exaggerated type 2 allergic skin inflammation, increased expression of Il-10, Ahr and Il-4, while Tnfrsf25 expression was reduced [13].
Given this previous literature, we hypothesised that Rora might be involved in regulating an inflammatory response via CD4+ T cells more widely. To investigate this we assessed the function of Rora in four different infection models involving, but not limited to, type 2 immunity (Fig 1).
The immune response to the infection models considered in this paper, over time, and which T helper cells are involved according to literature. The times denote approximate peak of cell counts according to previous studies, not the time points we have sampled.
By employing bulk RNA sequencing (RNA-seq), single cell RNA-seq (scRNA-seq) and mass spectrometry (MS) analysis of T helper cells, we have investigated the role played by Rora in the context of in vitro cell culture and in vivo infection models. Our data suggest that Rora is expressed by activated CD4+ T cells, and its expression correlates with effector function, such as the expression of Th lineage-defining transcription factors and the production of cytokines.
Materials and methods
Mice strains
Rora reporter mice.
Rora T2A Teal reporter mice [14] were generated briefly as follows: A reporter cassette, encoding a short Gly-Ser-Gly linker peptide, FLAG epitope tag, T2A self-cleaving peptide and Teal fluorescent protein, followed by a loxP-flanked Neomycin cassette, was inserted (directly upstream of the) Rora stop codon. Successful targeting of JM8 ES cells was confirmed by Southern blot and the Neomycin cassette was removed from the resultant mice by inter-crossing with a Cre recombinase strain.
Overexpression of Rora and analysis
Rora was cloned into a retroviral vector (M6-mCherry-Cter) as follows. cDNA was made from RNA extracted from mouse spleens, following the protocol of SuperScript II reverse transcriptase (Thermo Fisher Scientific #18064014). Primers were designed against Rora isoform ENSMUST00000113624 (fwd: ATGTATTTTGTGATCGCAGCGATGAAAGCTCAAATTGAA, rev: TTACCCATCGATTTGCATGGCTGGCTCAAATT). After PCR and gel purification, another round of PCR was done to add compatible overhangs (fwd: CGATaagcttATGTATTTTGTGATCGCAGC, rev: GaatgcggccgcCCCATCGATTTGCATGG). The product and backbone were digested with HindIII-HF and NotI-HF (NEB), purified, and ligated using T4 ligase (NEB). Several clones were picked and sequence verified by sanger sequencing (M6-Rora2-mCherry). Two other constructs were also constructed analogously for ENSMUST00000034766 (M6-Rora1-mCherry and M6-Rora1-GFP).
The virus production, activation and transfection was done as previously, with the addition of IL4 [15]. The transfection was done into 3 biological replicates (JaxJ, 8 weeks old), each with the 3 constructs and the no-insert M6 as control. No puromycin was used for selection. Cells were harvested on day 4, FACS sorted, and libraries prepared using the Sanger institute standard pipeline. Sequencing was done on a HiSeq 2500.
The differentially expressed genes were called by DESeq2 [16] (normalized counts ~ treatment + mouse). Each of the constructs were analyzed. Of these, one had Rora as the most DE gene and was retained for further comparison (M6-Rora2-mCherry). Counts and condition matrices are provided in S2 File.
Mice infection and cell isolation
Nippostrongylus Brasiliensis infection.
C57BL/6 female mice were subcutaneously injected with 100ul (300 live third stage Nippostrongylus Brasiliensis (N. Brasiliensis) larvae per dose), over two sites. MedLNs, mesLNs, spleens, lungs, or siLP were taken from infected mice, 3, 5, 7, 10 or 30 days after infection as well as from uninfected mice. At each time point, cells were isolated from medLNs and mesLNs, by smashing the tissue through 70μm cell strainers. Lungs were incubated in collagenase D (0.72mg/ml, Amersham, Bucks, UK) for 30min, smashed through 70μm cell strainers, and suspended in RBC lysis buffer (eBioscience Ltd). Spleens were smashed through 70μm cell strainers, and suspended in RBC lysis buffer. siLP cells were isolated from colonic tissue as previously described [17]. In brief, siLP cells were released by digestion of the tissue with RPMI/4-(2-hydroxyethyl)-1-piperazineethanesulfonic acid (HEPES) supplemented with 60 μg/ml DNaseI (Sigma), and 400 ng/ml of Liberase (Roche Applied Science, Burgess Hill, UK). Isolated cells were stained with conjugated antibodies for sorting (S3 File). For RT-qPCR and bulk RNA-seq experiments, single cells were sorted into the lysis buffer and were stored in -80°C.
Mice were weighed before the infection and daily during the experiment or until weights were deemed to be stable. Mice were monitored daily for clinical signs, including observations of the animal’s extremities, abdomen, coat condition, behavior and weight. If clinical signs were observed, monitoring was carried out more frequently. Mice showing signs of ill health that were likely to exceed a moderate severity, were killed by a Schedule 1 method. Humane endpoints were determined utilising a tailored score sheet that incorporated observations of an animal’s extremities, abdomen, coat condition, behavior and weight. No adverse events were experienced.
Schistosoma mansoni infection.
The complete life cycle of the parasite S. mansoni is maintained at the Wellcome Sanger Institute. Cercariae, the mammalian infective stage, were harvested by exposure of infected Biomphalaria glabrata snails to light for two hours in aquarium-conditioned water. Female C57BL/6 mice were infected with 300-mix sex S. mansoni cercariae via IP as described elsewhere [18]. Mice were checked twice a day for any sign of health deterioration. After 6 weeks post infection, mice were euthanize with an overdose of anaesthetic, adult worms were perfused from the mouse circulatory system as described previously [19], and mouse spleen and mesLNs were recovered surgically and placed in PBS for further processing. Cells were isolated from tissues and were stained for sorting. Activated T cells, CD3+CD4+CD44+SELL-, were sorted into a C1 chip and processed for RNA-seq.
The work using animals was carried out in accordance with the UK Animals (Scientific Procedures) Act under appropriate licenses issued by the Animals in Science Regulation Unit (ASRU) at the Home Office.
Lung immune challenge models.
Mice were anaesthetised by isofluorane inhalation followed by intranasal administration of Ragweed pollen (100 μg per dose; Greer Laboratories, Lenoir, NC), papain (10 μg per dose; Sigma) or 2W1S peptide (10 μg per dose; Designer Bioscience), dissolved in 40 μl PBS at the indicated time points. The PE-conjugated I-A(b) 2W1S tetramer was obtained from the NIH Tetramer core.
On the first day of dosing, more than one health check was conducted in order to monitor any reactions the animal may have had to the substance administered. When it was clear there were no reactions, the mice were health-monitored once daily. The health and welfare of animals were monitored with reference to the limiting clinical sign parameters listed in the “LMB Animal Usage Guidelines". In addition to this, the severity status of the procedure was determined by observations of the animal’s extremities, abdomen, coat condition, behavior and weight. Animals were observed after immunisation and any animals displaying signs of ill-health were killed by a Schedule 1 method. Humane endpoints were determined utilising a tailored score sheet that incorporated observations of an animal’s extremities, abdomen, coat condition, behavior and weight. No adverse events were experienced.
Cell sorting.
Cells were sorted with a BD Influx Cell Sorter or SY3200 Synergy cell sorter (iCyt) and analysed on Fortessa or LSRII BD cell analyzers using FlowJo. Surface and intracellular stainings were carried out according to the eBioscience protocols. A full list of the antibodies used is provided in S3 File.
Histology.
Lung lobes were fixed for 48 hours in 10% neutral-buffered formalin (Sigma), washed in PBS and transported to AML Laboratories (Florida) for embedding, sectioning and staining with haematoxylin and eosin.
Microscopy and image analysis.
Images were captured using a VHX keyence microscope, using objective #2D1510064 at magnification “200” (2.1495um/px). The 2D stitching feature was used to capture the entirety of each lung section. Two sections were imaged from each mouse, with focus automatically calibrated for each section. All images were acquired in one seating, with the same settings for brightness.
Image analysis was done in R using the EBimage package [20]. In brief, pixels corresponding to background, cytoplasm and nuclei were manually annotated in one image. A classifier was set up to separate pixels into the 3 classes based on manual inspection of the RGB space.
To estimate emphysema, the distribution of inner holes was calculated by first detecting non-background area. A binary dilation operation was applied to filter out 1–2 pixel holes caused by improper tissue detection. Segmentation was then performed to find disjunct regions. The area of each region was calculated.
From the histogram H of hole sizes, we estimate the unevenness (indicating emphysema, compensating for overall compression of slides due to fixation) as S = H99%/H80%. Other ranges were tested and gave a similar outcome. Veins, surrounding space, and other artifacts are filtered out above the 1% top percentile, while the lower percentile reflects tissue shrinkage. Alternatively, the score can be interpreted as a measure of hole size unevenness/variance. To test for emphysema in the Rora KO we fit a linear model S ~ replicate + genotype, with genotype = 1 for Cre and +/Flox, 2 for Cre and Flox/Flox, and 0 otherwise (p = 0.018). A simpler t-test of +/+ vs Flox/Flox, not using any of the +/Flox samples, and ignoring replicate batch effects, gives p = 0.17. This is not significant but in the direction of the more inclusive linear model.
In vivo bulk RNA-seq.
RNA was extracted with SPRI beads, followed by Smart-seq2. Read counts were estimated with Kallisto [21]. DE genes were estimated by DESeq2 [16] after counts had been rounded to the nearest integer count. Simple control vs treatment linear models were used. Counts and condition matrices are provided in S2 File.
Quantitative single-cell gene expression analysis
Single cell RT-qPCR measurement and analysis.
Single-cell gene expression analysis was performed using BioMark 96.96 Dynamic Array platform (Fluidigm, San Francisco, CA) and TaqMan Gene Expression Assays (Applied Biosystems, Carlsbad, CA). Single cells were sorted into 5μl of CellsDirect reaction mix and immediately stored in -80C. Control wells containing no cells were included. On thawing, a mix containing 2.5μL gene specific 0.2x TaqMan gene expression assays (Applied Biosystems), 1.2 μL CellsDirect RT/Taq mix, and 0.3 μL TE buffer were added to each well. RT–PCR pre-amplification cycling conditions were: 50°C, 15min; 95°C, 2min; 22x(95°C, 15s; 60°C, 4min). Samples were diluted 1:5 in TE buffer and 6% were mixed with TaqMan Universal PCR Master Mix (Applied Biosystems). The sample mix and TaqMan assays were loaded separately into the wells of 96.96 Gene expression Dynamic Arrays (Fluidigm) in presence of appropriate loading reagents. The arrays were read in a Biomark analysis system (Fluidigm). ΔCt values were calculated in reference to the average of Atp5a1, Hprt1 and Ubc.
The expression of each gene was fit by a bivariate normal distribution (R mixtools package [22]) and a cut-off set at the average position between the lower and upper gaussian midpoints. Genes above this were considered expressed if above this cut-off. qPCR levels and metadata are provided in S1 File.
S. mansoni C1 single cell RNA-seq and analysis.
For S. mansoni, three small (5–10 μm) C1 Single-Cell Auto Prep IFC chips (Fluidigm) were primed and 5000 cells were sorted directly into the chip. To allow estimation of technical variability, 1 μl of a 1:4000 dilution of ERCC (External RNA Controls Consortium) spike-in mix (Ambion, Life Technologies) was added to the lysis reagent. Cell capture sites were visually inspected one by one using a microscope. The capture sites that did not contain single cells were noted and were removed from downstream analysis. Reverse transcription and cDNA preamplification were performed using the SMARTer Ultra Low RNA kit (Clontech) and the Advantage 2 PCR kit according to the manufacturer’s instructions on the C1 device. cDNA was harvested and diluted to 0.1–0.3 ng/μl and libraries were prepared in 96-well plates using a Nextera XT DNA Sample Preparation kit (Illumina) according to the protocol supplied by Fluidigm. Libraries were pooled and sequenced on an Illumina HiSeq2500 using paired-end 75-bp reads.
Salmon was used to estimate gene expression counts [23]. Poor quality libraries were eliminated using Scater [24] based on exonic and mitochondrial read counts. For all queries, a gene was considered expressed if the Log10 (1+normalized count) was above 0.5. Gene overlap was tested using Fisher’s method. Counts and condition matrices are provided in S2 File.
N. brasiliensis C1 single cell RNA-seq.
T cells from N. brasiliensis infected mice were captured with C1 analogously to the S. mansoni mice. The cells were sorted for CD3+CD4+SELL-, from mesLN, medLN, lung and spleens. This data is available but was not used due to low quality.
N. brasiliensis Smart-seq2 single-cell RNA-seq and analysis.
For this analysis we included CD3+CD4+SELL- cells from the lungs at 30 days after infection, CD3+CD4+Rorateal+ and CD3+CD4+Rorateal- cells from N. brasiliensis infected Rora+/teal reporter mice. Further we include cells from lungs, spleens, medLN and mesLN from infected and uninfected Rora+/teal mice, 7 days after infection.
Single cell transcriptomes were generated by using the Smartseq-2 protocol [25]. Single cells were sorted into 96 well plate that contained 5μl of Triton-X lysis buffer, 1 μl of 10 μM olido-dT30-VN, 1 μl dNTP mix (25 mM each) and ERCC controls at a final dilution of 1:64m and immediately stored in -80°C. cDNA was submitted for the amplification for 25 cycles on an Alpha Cycler 4 thermal cycler. Amplified cDNA went through two subsequent rounds of cleaning by using Agencourt AMPure XP beads (Beckman Coulter UK Ltd, High Wycombe, UK) at a 1.0x ratio on a liquid handler Zephyr G3 NGS Workstation (PerkinElmer) and eluted in 20ul of RNAse free water. Purified cDNA was subjected to quality control using 1 μl of cDNA on an Agilent 2100 BioAnalyser (Agilent Technologies, Santa Clara, CA, USA) using the Agilent High Sensitivity DNA kit. Samples were normalised to a concentration of 0.3 ng/μl. Nextera libraries were prepared using Nextera XT DNA Sample Preparation kit (Illumina) according to the protocol supplied by Fluidigm. The library preparation has been done in combination with benchtop liquid handlers: Zephyr G3 NGS Workstation (PerkinElmer) and Mantis (Formulatrix). Purified pool of Nextera libraries was subjected to quality control using 1 μl of cDNA on an Agilent 2100 BioAnalyser using the Agilent High Sensitivity DNA kit.
Salmon v0.8.2 [23] was used to map the reads to ensembl cDNA reference GRCm38 to quantify the gene expression. Counts were used for all the analyses. Low quality cells with >15% mitochondrial genes or <1500 total genes were filtered out. Seurat [26] was used for single-cell analysis: Highly variable genes were calculated (x.low.cutoff = 0, x.high.cutoff = 5, y.cutoff = 0.1), PCA performed on these (pcs.compute = 50), and t-SNE on the first 25 PCA components. Finally clustering was done (reduction.type = "pca", dims.use = 1:25, resolution = 0.5). Common marker genes were plotted on top of the clusters and the cluster identities were assigned manually (Gata3, Tbx21, Foxp3, Sell, and Rorc). The TCR sequences for each single T cell were assembled using TraCeR [27] which allowed the reconstruction of the TCRs from scRNA-seq data and their expression abundance, in addition to identification of the size, diversity and lineage relation of clonal subpopulations. Cells for which more than two alpha or beta chains were identified were excluded from further analysis. iNKT cells were detected by their characteristic TCRA gene segments (TRAV11–TRAJ18). Only CD4+ T cells were kept in the final plot; other cells remain annotated in the ArrayExpress submission. Counts and condition matrices are provided in S2 File.
Screening for upstream regulators of Rora
Cell culture.
Naive CD4+ T helper cells were extracted from spleens of JaxJ mice, 8 weeks old males (Stem Cell Technologies #19765, CD8-CD11b-CD11c-CD19-CD24-CD25-CD44-CD45R-CD49b-TCRγ/δ-TER119-), according to the manufacturer’s protocol. Up to 4 spleens were pooled in each biological replicate. The cells were activated in 96-well round bottom plates that were coated with 3ul/ml anti-CD3e (BioLegend #100202) and 5ul/ml anti-CD28 (BioLegend #102102) for 4 hours at 32°C. 24 hours later cytokines were added (See S2 File, condition table). 72 hours after activation the cells were pelleted, resuspended in RLT buffer (Qiagen #79216), and stored at -80°C.
For the repeated screen, culturing was done in the same manner. The cells were then stained with PI and anti-CD4 eFluor 660 (eBioscience 50-0041-82, clone GK1.5) for 30 minutes in IMDM media. The cells were then FACS sorted to retain live CD4+ T cells, on average 5000 cells per well.
During Th0 conditions, it is possible for the cells to become somewhat Th1 polarized; we have looked at markers in the data but can neither confirm nor disprove that this is the case for these experiments.
RNA-sequencing.
For the first screen, RNA was extracted using Agencourt AMPure XP beads (Catalog #A63881) at 1.5:1 ratio into 10ul elution volume. RNA-seq was done following Smart-seq2 [25] but with the following changes: The input was increased to 7ul and Mix 1 consisted of 2ul 10uM dNTP and 2.5ul 10uM oligo-dT. Mix 2 consisted of 5ul SmartScribe buffer, 5ul 5M Betaine, 1.25ul SmartScribe, 1.25ul (100uM) DTT, 0.15ul 1M MgCl2, 0.38ul 100uM TSO, and 0.47ul SUPERase inhibitor. SPRI cleanup was done at ratio 0.8:1 following first-strand synthesis, eluted into 12ul NFW. Mix 3 consisted of 12.5ul KAPA HiFi 2x master mix and 0.5ul 10uM ISPCR. PCR was done with 10 cycles. In the second screen, RNA-seq was done similarly with the following modification: No SPRI-cleanup was done after first-strand synthesis. Instead a larger amount of mix 3 was added directly.
On average 0.2ng cDNA was then used for Nextera XT library preparation, done at 25% of volume recommended by manufacturers protocol. Sequencing was done with 2 lanes of Illumina 50bp PE HiSeq 2500.
Upstream analysis.
Salmon was used to estimate gene expression counts. Poor quality libraries were eliminated using Scater based on exonic and mitochondrial read counts. Genes with average normalized counts below 0.3 were removed. The increase/decrease in Rora expression level was computed by DEseq2 [16] using a linear model: normalized counts ~ treatment + mouseReplicate, where treatment is a factorial variable. The correlation heatmap between treatments was calculated as cor(cor(FC)). The first correlation is used as a dimensionality reduction while the second correlation is used to remove the technical bias toward positive correlation. Clustering was done using R hclust, standard parameters.
The CellTrace FACS data was analyzed with FACSanadu [28] to obtain cells in each division stage. The cell counts were analyzed across replicates and compared to controls using R.
Metabolite LC-MS
Metabolites were extracted using a modification of a previous method [29]. Briefly, roughly 200k cells were resuspended in chloroform/methanol (2:1). After mixing samples, water was added (300 μL) and centrifuged at 13,000 x g for 20 min. The organic layer (lower fraction) was collected and dried under a stream of nitrogen. The organic fraction was reconstituted in chloroform/methanol (1:1; 100 μL) and 10 μL were mixed in 90 μL IPA/acetonitrile/water (2:1:1). The samples were analysed by LC-MS using an LTQ Orbitrap Elite Mass Spectrometer (Thermo Scientific). Sample (5 μL) was injected onto a C18 packed-tip column (75μm x 100mm, Thermo Scientific) at 55°C. The mobile phase A was acetonitrile/water 60:40, 10 mM ammonium formate whilst B was LC–MS-grade acetonitrile/isopropanol 10:90, 10 mmol/L ammonium formate. In negative ion mode, ammonium acetate was used to aid ionisation. In both positive and negative ion mode, the gradient used (S5 File) was run at a flow rate 0.5 mL/min. Data were acquired with a mass range of 100–2000 m/z. Chromatograms were converted to mzML format and annotated using an in-house R script. The annotated metabolites were size factor normalized to enable comparison. PCA was performed on the in vitro and in vivo samples separately (S1 Fig in S1 File). Lungs and spleens from in vivo samples appear to cluster together, suggesting no need to correct for tissue differences. In vitro samples were compared with Limma [30], with all naive cells (negatively selected as CD8-CD11b-CD11c-CD19-CD24-CD25-CD44-CD45R-CD49b-TCRγ/δ-TER119-) compared to all day 6 samples Th0/1/2. For the in vivo samples, we used the simplest model of ~1+isko. Bayesian error correction was performed.
Results and discussion
Rora is expressed in activated CD4+ T helper cells in vivo
We wanted to get an overview of Rora expression in CD4+ T helper cells during worm infection. As most of the mouse models used in this study are type 2 immunity, we focus more specifically in Th2 cells. Mice were infected with Nippostrongylus brasiliensis (N. brasiliensis), a tissue migrating parasitic roundworm of rodents. N. brasiliensis initiates infection at the skin, the site of invasion, and then migrates via the circulatory system to the lungs. Through exploitation of host lung clearance, worms are transported to the intestine, from where they are expelled. In total, this infection is resolved within 7–10 days (Fig 1).
To get an overview of the infection at the molecular level, we performed single cell quantitative reverse transcription PCR (scRT-qPCR, Fluidigm Biomark). We prepared cells from the lungs, small intestine lamina propria (siLP), mediastinal lymph nodes (medLN) and mesenteric lymph nodes (mesLN) of mice, 3, 5 and 7 days after N. brasiliensis infection (Fig 2A). Tissues were also collected from uninfected mice as controls. Both CD3+CD4+SELL+ (naive) and CD3+CD4+SELL- (activated) cells were sorted. We examined 440 individual CD3+CD4+ T cells for the expression of 96 genes including housekeeping genes, different immune cell markers, T helper markers, cytokines, cytokine receptors, chemokine receptors, known and candidate regulators for Th cell differentiation (Fig 2B; S2 File). Consistent with a predominantly type-2 immune response, 70% of the CD3+CD4+SELL- cells expressed Gata3, whilst 26% expressed Foxp3 and only 2% contained Tbx21 or Rorc. Interestingly, in this population, we noted a significant overlap of Gata3 and Rora (Fisher, p = 4*10−4), and also of Rora and Foxp3 (Fisher, p = 2*10−4) (Fig 2B).
(a) Experimental design for single-cell RT-qPCR. (b) Gene expression according to RT-qPCR from CD4+ T cells during N. brasiliensis infection. Gata3hi cells, which we assume are Th2, are also expressing Rora. The category “Cytokines” refers to any of Ifng, Il4, Il13, Il5, Il6, Il10 or Il17a being expressed. Top row denotes spatial location: Lung, gut (siLIP), mesenteric and mediastinal lymph nodes (n = 4 mice). Every column is one individual cell. (c) Flow cytometric analysis of CD4+ T cells from the spleen and lung of Rora+/teal naive mice (black histogram = Rorateal–cells, green histogram = Rorateal+ cells) (n = 23). All error bars are ±sample s.d. (d) Frequency of Rorateal expression across tissues, with siLIP having the highest number (n = 12). (e) Fraction of naive (CD4+CD25-CD44-SELL+) or activated cells expressing Rorateal. RORA is mainly expressed in activated cells in vivo (n = 12). (f) N. brasiliensis Rora expression over time, overlaid with previously published mouse and human in vitro time course data of Rora, showing no or low expression of Rora in activated in vitro T cells. (g) Experimental design of scRNA-seq experiment, with (h) tSNE clustering of CD4+ T cells extracted from N. brasiliensis infected mice and measured by single-cell RNA-seq (1670 cells after filtering). We primarily identify Naive, Th1, Th2 and Treg cells, of which Rora is expressed in all non-naive cells. See methods and the provided code for the precise workflow. (i) Expression distribution of selected genes from the scRNA-seq data.
To investigate Rora activity in activated Th cells, we examined Rora expression using a Rora+/teal reporter mouse, in which the gene encoding Teal fluorescent protein has been inserted at the 3’ end of the Rora coding sequence and expressed from a self-cleaving T2A site [14]. Flow cytometric analysis of reporter gene expression indicated that Rora was expressed in a proportion of CD4+ T cells isolated from medLN, spleen, lung and siLP (Fig 2C). We noted that all Rora expressing cells possessed an activated phenotype, as evidenced by their expression of CD44, ICOS and CD38 and down-regulation of SELL (Fig 2D). Accordingly, Rora+/teal T cells were therefore significantly enriched within the activated T cell subset, comprising 40–60% of this population, but absent from the pool of naive cells (CD4+CD25-CD44-SELL+) (Fig 2E). This correlation of Rora expression with T cell activation is in contrast with RNA-seq analysis made during the first 72h of in vitro culture of Th0 and Th2 cells, as previously reported (Fig 2F) [15]. There, in both mouse and human, Rora is expressed in naive cells, with a sudden increase around 2–4 hours of activation, before a gradual decrease in expression (overall slightly higher expression in Th2 over Th0). This agrees with another in vitro activation dataset where the Rora expression level is unchanged and low, in both Th17 and Th1 cells, following activation [31]. To compare with the in vivo situation, we generated RNA-seq data from CD4+ T cells, isolated from the spleens of N. brasiliensis infected wild-type mice 1, 3 and 10 days after infection, and overlapped with the in vitro data (Fig 2F). Here Rora is seen to increase, reaching a high and fairly sustained level around day 3. This is in agreement with a previous single-cell analysis of Th1/Tfh cells during malaria infection where Rora expression is seen to increase over time in Th1 [32] (S2A, S2B Fig in S1 File). It is also consistent with a study on atopic dermatitis [10], a study on Tregs [33] (S2C Fig in S1 File), as well as on T cells in a melanoma model [34] (S2D Fig in S1 File). Taken together, as Rora is consistently present in vivo but not in vitro, these 6 datasets suggest that Rora is induced not primarily by activation, but possibly by an external cue from the micro-environment.
Because the particular (micro-)environments might influence Rora expression, we examined our scRT-qPCR data across different organs (Fig 2B). Overall, the expression of Rora was relatively constant in CD3+CD4+ cells 3–7 days after infection (Fig 2F). However, in terms of tissue, it was expressed more frequently in the siLP and lungs, where 37% of cells expressed Rora, compared to 20% of the cells in lymph nodes. This pattern was confirmed using flow cytometry, where 10.9±1.9% of CD4+ cells in the lung were Rora+/teal, 10.1±1.1% from the spleen, 27.9±5.8% from the small intestine lamina propria and only 4.0±0.6% from medLN (Fig 2D). This tissue distribution likely reflects the relative frequency of activated Th cells in each location.
Rora has so far been seen in activated cells, partly overlapping with primary markers of Th2 and Treg (Gata3 and Foxp3). To fully confirm the presence of Rora among different T helper cell types, as defined by the full transcriptional program, we used scRNA-seq on lymphoid and non-lymphoid tissue, isolating CD3+CD4+SELL- cells from the lungs 30 days after infection (N. brasiliensis infected wild-type mice 1, 3 and 10 days after infection, and overlapped with the in vitro data (Fig 2G). As reference points we included CD3+CD4+Rorateal+ and CD3+CD4+Rorateal- cells from N. brasiliensis infected Rora+/teal reporter mice. Lungs, spleens, medLN and mesLN were also taken from infected and uninfected Rora+/teal mice, 7 days after infection. The cells were clustered and each cluster annotated using common markers (Sell, Ifng, Gata3, Foxp3) as Naive, Th1, Th2 and Treg cells (Fig 2H, all cells are shown in S2E Fig in S1 File). Cells from the lymph nodes and spleen are mainly Naive or Treg, otherwise the cell origin and T cell fate appear fairly uncorrelated. The distribution of Rora expression level across the different CD4+ T cell types was investigated (Fig 2H and 2I). Again, we note that Rora is expressed among all activated T cells (however, less induced in Th1).
To confirm that this distribution of Rora expression is not specific to the response to N. brasiliensis infection we compared it to another worm model, Schistosoma mansoni (S. mansoni). CD3+CD4+SELL- cells were collected 6 weeks after infection from spleens and mesLN (S2F Fig in S1 File). Single cell RNA-seq was performed and analysis shows that Rora is expressed in 59% of the activated cells (expressed defined as greater than zero counts). Approximately all (98%) cells that expressed Rora also expressed any of the key regulators Gata3, Foxp3 or Tbx21.
To conclude, we have performed FACS, scRT-qPCR and scRNA-seq, and confirmed that Rora is generally expressed in a proportion of activated T cells during N. brasiliensis infection. Rora induction is likely to depend on a local environmental cue not present in the usual CD3/CD28 in vitro activation system. Thus an in vivo system is required to study the physiological role of Rora.
Rora is expressed in activated T helper cells during immune responses
To further corroborate these findings, we also investigated a non-worm immunity model, Ragweed pollen (RWP), a common allergen which causes lung inflammation [35]. The Rora+/teal reporter mice were exposed to intranasal RWP or PBS (Phosphate Buffered Saline, a negative control administration), administered every 2–3 days for 2 or 6 weeks (Fig 3A). As anticipated, RWP exposure resulted in increased CD4+ T cell counts in the lung and mediastinal lymph node (medLN) after 2 weeks as compared to PBS controls, and this was maintained after 6 weeks of continuous RWP administration (Fig 3A). In lung, but not medLN, this increase in CD4+ T cell frequency was accompanied by an elevated proportion of T cells expressing Rora, increasing from 30% to 40% between 2 and 6 weeks of RWP administration (Fig 3B and S3 Fig in S1 File). Consistent with our findings in untreated mice, Rora+/teal T cells elicited by RWP stimulation had an activated phenotype, characterised by a CD44+SELL- surface profile (S3 Fig in S1 File). Notably, the Rora+/teal T cell subset was enriched for cells that had acquired effector function, identified by their expression of the cytokines IFNg, IL-13 and IL-17A, or the Treg-associated transcription factor, FoxP3 (Fig 3C). Here again we note that Rora expressing cells include diverse activated cell phenotypes, that is, Th1, Th2 and Th17.
All error bars are ±sample s.d. (a) Design of ragweed flow cytometry experiment. (b) Flow cytometry analysis of CD4+ T cells from mice challenged with Ragweed pollen protein for 2 or 6 weeks, showing that Rora is expressed in lungs also during allergy. The Tetramer+ cells fraction among all T cells ranged between 0.5–2%. Here the number of Rorateal+/- cells are shown over time (n = 12, one repeat). (c) RORA expression in CD4+ cells increases in the lungs but not in MedLN during the ragweed challenge (n = 12, one repeat). (d) Flow cytometry analysis diagrams showing higher proportions of cytokines-expression and Tregs cells in Rorateal+ sorted cells vs the Rorateal-, after Ragweed pollen protein exposure (n = 2). IL17 and IFNg were also measured (S2B Fig in S1 File) (note that we only did this experiment once). (e) 2W1S/papain challenge experimental design. (f) Flow cytometry analysis when T cells are challenged twice with 2W1S/papain, again confirming the presence of Rora in tetramer+ cells. Note that we only performed this experiment once (n = 3+3). (g) Flow cytometry analysis diagrams showing higher proportions of cytokines-expression and Tregs cells in Rorateal+ sorted cells vs the Rorateal-, after 21WS/Papain challenge (note that we only did this experiment once).
To understand the involvement of Rora in effector cells, we examined its expression in cytokine-expressing cells. From the N. brasiliensis model single-cell RT-qPCR, 83% of the cells that expressed either IL-4, IL-13 or IL-10, also expressed Rora (Fig 2B). The corresponding number in the S. mansoni-infected mice is 71% (S2F Fig in S1 File, scRNA-seq).
To investigate whether Rora was expressed by Th cells responding to secondary immune challenge, we challenged Rora+/teal and control mice with 2W1S peptide intranasally, in conjunction with the protease allergen papain, which promotes a type-2 immune response (Fig 3D). After secondary immunisation with 2W1S/papain 3 weeks after the initial stimulation, we were able to identify 2W1S- responsive CD4+ T cells by tetramer staining [36] and FACS analysis. Tetramer+ CD4+ T cells in both the lung and medLN were largely CD44+ and Rorateal+ (Fig 3E), indicating that recently activated Th cells express Rora. As with the N. brasiliensis infection, in the S. mansoni model we found significant overlap of Foxp3 and Rora in activated T cells (S2F Fig in S1 File, Fisher test, p = 2*10−4). A higher fraction of FOXP3+ cells was also found in sorted Rora+/teal relative to RORA- cells, in the RWP model (Fig 3C). The gap is increased when the exposure to the pollen was longer. The same could also be seen after mice were exposed to 2W1S+papain. The Rora+/teal population contained 20% of the FOXP3+ cells, while the RORA-, only 4% of the cells.
Overall our data show that Rora is expressed in Th cells in worm immune models as well as in allergen type 2 immunity. We find Rora transcripts present in activated Th cells, not only in Tregs, but also in Th1 and Th2, as well as activated CD4+ T cells responding to a secondary immune challenge. Rora is also more correlated with cytokine- releasing cells, which suggests a function in effector cells.
Rora KO affects inflammation severity
To isolate the in vivo physiological function of Rora in Th cells, we created a conditional knockout mouse line, where exon 4 of the Rora gene was deleted only in CD4+ cells (CD4creRorafl/fl, verified by RNA-seq, S3A Fig in S1 File), resulting in a dysfunctional RORA protein. We noted no change in Th cell number when comparing CD4creRorafl/fl to control (S3B Fig in S1 File). We then infected CD4creRorafl/fl mice with N. brasiliensis and examined lung inflammation after 30 days. An unbiased image analysis algorithm was developed to detect emphysema and score stained lung sections (Fig 4A). Generally, uninfected samples score the lowest for emphysema, as expected (Fig 4B and 4C). During infection, single-allele Rora KO (CD4creRora+/fl) had a higher score of emphysema, and double-allele KO (CD4creRorafl/fl) even higher score (linear model fit, p = 0.018). We detected no differences in our negative controls, no-Cre KO mice (Rorafl/fl), and WT with Cre KO mice (CD4cre). As emphysema can be caused by eosinophil infiltration, in analogy with a previous model [13] we analyzed eosinophil infiltration rate and noted an increase but not statistically significant difference (S3C Fig in S1 File).
(a) Assessment of lung emphysema after N.brasiliensis infection, day 30 past infection, by microscopy and automated image analysis. Controls are included for no-Cre (against the Rora-Flox) and no infection (from left to right, n = 8, 17, 1, 14, 12, 9, 10, 29, 20, 39, 6, 6, 2, 2, 6, 6). (b-c) Lung lobes were fixed, embedded, sectioned and stained with haematoxylin and eosin. The pictures represent the least and most inflamed lung sections as detected by the algorithm, located in +/+ and Flox/Flox (control and KO).
We conclude that Rora expression in CD4+ T cells modestly promotes lung inflammation associated with N. brasiliensis infection. This may reflect a role for Rora in regulating the activity of T helper cells, and may also result from the modulation of Treg activity, as has previously been reported for Rora+ Tregs in skin [13].
Rora overexpression affects several key immune regulatory genes
To probe the role of Rora in T cells, we next proceeded to identify the molecular targets of Rora in Th2 cells by means of retroviral overexpression. Th2 cells were induced in vitro, infected the day after induction, and on day 5, gene expression was compared by RNA-seq (empty viral vector vs overexpression vector, Fig 5A and 5B). The in vitro differentially expressed (DE) genes were compared to the genes that are differentially expressed between skin Tregs from wt vs CD4creRorafl/fl mice (Fig 5C). Consistent with a previous report [13], we see a strong effect on tumor necrosis factor receptor superfamily member 25 (Tnfrsf25), which has been shown to be required for Th2 effector function, e.g., allergic lung inflammation [37]. The receptors C-C chemokine receptor type 2 (Ccr2), Ccr5, and Tnfrsf23 are in agreement, although with lower fold changes. We also find Ninjurin 2 (Ninj2) which has recently been implicated in endothelial inflammation [38].
To increase readability, only some gene names are written out (full DE gene data is in S2 File). (a) Design of the in vitro overexpression experiment. (b) Rora overexpression. (c) Previous data, DE genes in Treg from skin [13] (d) Design of the in vivo KO experiment. (e-h) DE genes in vivo, activated CD4+ T cells (n = 3+3), from lung, colon, siLP and spleen.
We also looked at a previous dataset of Rora small interfering RNA (siRNA)-treated human Th17 cells [39] and find that Tnfrsf25 is DE (p = 2.2*10−4, DESeq2). In another study of Tregs during steady-state in mouse [40], we find that Tnfrsf25 is higher in colon Tregs than in lymphoid tissue Tregs (p = 2*10−2, DESeq2). Thus four datasets, in different T helper cell types, confirm Tnfrsf25 as a downstream gene.
Overall the two in vitro systems strengthen the claim of Rora as regulator of Tnfrsf25. However the expression of Tnfrsf25 across Th2, Treg and Naive, but not Th1 (Fig 2I), shows that the function of Rora is not limited to Treg cells.
Rora affects activation in vivo
To see if Rora has additional functions in vivo, we generated bulk RNA-seq from activated (SELL-/CD44+) CD4+ T cells from control and CD4creRorafl/fl mice from lung, spleen, siLP and colon (30 days after N. brasiliensis infection), and compared the DE genes between the different tissues, the overexpression DE genes, and the previous Treg Rora KO data. Volcano plots are shown in Fig 5E (full list of fold changes in S2 File). From manual comparisons of these conditions, we have tried to find DE genes in common between the tissues, and if these in turn can explain the function of Rora. We generally do not see any DE cytokines of high statistical significance, but highlight some other types of genes here.
We find the transcription factor Aryl hydrocarbon receptor nuclear translocator-like (Arntl) to be differentially expressed in activated CD4+ T cells from siLP as well as in the overexpression analysis and in Tregs in skin (Fig 5B and 5C). Arntl is involved in the circadian responsiveness, and is expressed in CD4+ T cells [41]. The Arntl promoter contains ROR elements (RORE), and was found to be regulated by Rora in mice [42]. RORA promotes Arntl expression and thereby maintains the circadian rhythm. We find it to be expressed in all cell types (Fig 2I).
The S100 calcium-binding protein A4 (S100a4) is another gene which is differentially expressed in activated T cells from lung, colon, and to lesser extent in Tregs in skin. In the single cell dataset, it is expressed mostly in Tregs but generally in Th1/Th2 (Fig 2H). It exists in both intra- and extracellular forms and can activate NF-κB [43]. Using an antibody, S100a4 has been shown to mediate T cell accumulation, and Th1/Th2 polarization, at tumour sites [44]. Thus S100a4 might act downstream of Rora to control T cell activation.
C-X-C chemokine receptor type 6 (Cxcr6) is differentially expressed in siLP, spleen and in Tregs in skin, and is consistently regulated in other tissues. Cxcr6 is expressed in all Th cells but mostly in Th2 cells (Fig 2I). It might be involved in Th cell activation [45, 46].
The sialyltransferase St6galnac3 is downregulated in colon and overexpression (Fig 5B and 5F). This enzyme transfers sialic acids to glycolipids and glycoproteins. Along with the previously DE gene Alox8 (Fig 2I), we wondered if Rora may have any function on lipid metabolism. RORA has already been found to bind cholesterol [47] and is evolutionarily related to the nuclear receptor and lipid metabolism gene Peroxisome proliferator-activated receptor gamma (Pparg) [48]. We performed lipid LC-MS on in vivo CD4creRorafl/fl and WT cells from N. brasiliensis infected mice, as well as on in vitro WT naive and mature Th0/Th1/Th2 cells for reference (S4 Fig in S1 File, cells from negative magnetic bead selection as described in methods). The CD4creRorafl/fl presented small differences, and little overlap with the 85 lipids changed in vitro the first 6 days (cutoff p = 0.05). However, all the top DE lipids are phosphocholines or phospoethanolamines, and some phosphocholines can induce a type 2 immune response [49]. All the differentially expressed lipids are listed in S5 File.
In the spleen, Sell is downregulated in the KO (however it is still expressed as other cells are included as well), and opposite in the previous overexpression. This fits the idea of Rora driving activation and exit from the spleen, but why is it not DE in the other organs? We speculate that once a T cell has left the spleen or a lymph node, Rora might have already filled its purpose, and thus the effect on Sell will be smaller in any peripheral tissue.
We further noted the gene Emb (Embigin) in which is significantly DE in siLP, and in the same direction in overexpression, lung and spleen. It was also strongly DE in the previous skin Treg dataset. It has been shown to regulate cell motility in pancreatic cancer, and be controlled through the TGF‐β pathway [50]. It is possible that it could have a similar role in T cells.
Overall we see that Rora affects several genes in vivo, that are differentially expressed between different T helper cell types. Several of the downstreams genes may explain the effect of Rora, including for example Arntl, Cxcr6, St6galnac3, Emb and Sell. However, Rora abbrogration does not abolish T cells. It is possible that Rora acts as a helper factor, or tissue residence checkpoint, for T cell activation and migration.
Rora expression can be extrinsically regulated
A remaining question is which genes regulate Rora, and whether regulation is solely due to the T cell itself or from an external signal. More insights into when Rora is induced would provide understanding of its role. As previously seen, Rora is generally not expressed during in vitro activation, but is upregulated during in vivo activation. This discrepancy may be explained by an unknown signal from the environment that is present in vivo but not in vitro. To find candidate upregulators we screened for in vitro upregulated genes during 39 treatments, including cytokines and chemokines, with 15 unstimulated controls, in 3 biological replicates. In short, naive cells were taken from wild type mouse spleens, plated on CD3/CD28 coated plates, with addition of test cytokines/chemokines 24 hours later. 5 days after plating the cells were measured by bulk RNA-seq (Fig 6A).
The full list of stimuli, and the effect on any gene of interest, is provided in S2 File. All error bars are ±sample s.d. (a) Design of the in vitro experiment (b) Effect of all tested stimuli on Rora expression. (c) Global correlation of gene fold change in each stimuli, with fold change of gene expression during Rora overexpression. Lines have been drawn across panels b-c to highlight how similar a treatment is to directly perturbing Rora.
The stimulants and their effects are shown in Fig 6B (DE genes in S2 File). To validate our perturbation assay we looked at known marker genes. Of the genes having the strongest effect we find the expected genes: IL4 for Gata3 (Th2), Transforming growth factor beta 1 (TGFb1) for Foxp3 (Treg), IL6 and partially TGFb1 and IL21 for Rorc (Th17) [51] (S5a-S5c Fig in S1 File). We then looked for the effect on Rora expression. The most upregulating stimulants are CCL5, IL15 and IL33.
CCL5 is chemotactic for T cells, eosinophils, and basophils. CCL5 signaling through CCR3 has been reported to regulate Th2 (IL4+CD4+) cellular responses to promote metastasis of luminal breast cancers [52]. Interestingly also CCR5, another CCL5 receptor, was upregulated during overexpression, implicating a potential autocrine role for this factor. In vitro, CCR3 and CCR5 are expressed primarily in Th1, Th2 and Th17 [31].
IL15 is produced by non-lymphoid cells and important for the survival of several lymphoid subsets [53]. IL15 KO mice are markedly lymphopenic due to decreased proliferation and decreasing homing to peripheral lymph nodes [54]. Increase in IL15 thus fits our observation of Rora expression primarily in activated Th cells outside the lymph nodes. Il15ra is expressed in all murine CD4+ T cell subsets [31].
IL6 resulted in upregulation of Rora (statistical test significant if both batches of IL6 are pooled), consistent with previous IL6 in vitro culture data [31]. IL6 has several roles, beyond generating Th17 cells in vitro, both pro- and anti-inflammatory [55]. Il6ra is expressed in all CD4+ T cell subsets [31], though in our scRNA-seq data less in Th1 (Fig 2I).
IL33 also induces Rora, and is mainly associated with Th2 or ILC2 cells [56]. This agrees with our in vivo scRNA-seq where Rora is slightly elevated in Th2 cells compared to other Th subtypes (Fig 2I). Our scRNA-seq shows that the corresponding receptor Il1rl1 is expressed in all CD4+ T cell subsets but primarily in Th2 and Treg (Fig 2I).
Ccl7 upregulates Rora and can do so through several receptors; these include for example CCR1, CCR2, CCR3, CCR5, and CCR10 [57]. All of these receptors are present in CD4+ T cells to some extent, but Ccr5 is somewhat higher expressed than the rest [31].
The most downregulating stimulants are C-C motif chemokine 22 (CCL22) and stromal cell-derived factor 1 (SDF1a). CCL22 is secreted by dendritic cells and macrophages, and is suggested to induce Treg cell infiltration into the pleural space in patients with malignant effusion [58, 59]. In addition, its receptor CCR4 is required for CD4+ T cell migration to the skin [60]. In bulk in vitro data it is present in all CD4+ T cell subsets, especially after activation [31].
SDF1a (also known as CXCL12), a ligand for the chemokine receptor CXCR4, promotes bone marrow homing for T cells. It has been shown that Naive T cells downregulate CXCR4 to avoid this [61]. A further look at the Human Protein Atlas shows the tissues with expression of SDF1a to be particularly high in spleen, cervix/endometrium, and adipose tissue. In our study we compared with lung and gut, which have less Cxcl12 expressed. This is consistent with our screening data, suggesting SDF1a may be a negative regulator of Rora. Cxcr4 is expressed in all CD4+ T cell subsets at roughly equal level [31] according to bulk and our scRNA-seq data (Fig 2I).
To further validate if a treatment induces Rora, we also checked if a treatment induces the same downstream genes as in the Rora overexpression experiment. To do this we took genes from the Rora overexpression with a Log2 fold change > 2, and correlated the fold change with the corresponding fold change in a treatment (Fig 6C). Many of the most Rora influencing treatments also score similarly by this score. In particular, IL33 and CCL7 have a high level of agreement, followed by CCL22 (replicate a). Other treatments appear to have a less direct effect on Rora, for example IL6.
Based on our and previous data we conclude that Rora expression is mainly associated with an inflammatory state. Several cytokines may influence the expression of Rora. The effect on other genes can be studied using a similar screening approach, based on our supplementary data.
Conclusions
We have here investigated the impact of Rora on different subsets of CD4+ T cells, with an emphasis on Th2 cells, in several in vivo immune models. Rora is already known to be important for Th17 identity [12], and more recently it was discovered to regulate the function of Tregs in the skin [13]. Here we show that Rora is actually expressed in all activated T helper cells during different types of inflammation, including worm infections and allergy. We expand on the previous literature and suggest that Rora plays a wider role in the regulation of the immune response, in many T helper cell types—Th1, Th2 and secondary activated cells as well. Its function is likely to also involve further key downstream genes than just the previously shown Tnfrsf25. Our model of the Rora regulatory network is depicted in Fig 7.
Rora is expressed upon activation, and relies on additional stimuli from peripheral tissues for full induction of gene expression. Several cytokines have an impact on Rora expression level. Rora affects a diverse set of genes of importance for activated T cells, including Tnfrsf25 (negative regulation of Eosinophils during infection), Arntl (circadian rhythm), S100a4 (activation), Alox8 and Stgalnac3 (lipid metabolism).
Rora deletion did not alter the overall CD4+ T cell numbers, and therefore does not seem to be involved in T cell proliferation. This was also shown before in a study where Th2 cells were examined in RORα-deficient Staggerer mice [62]. In this study we found that Th2 cells from the RORα-deficient mice were still able to differentiate and proliferate. However, we noticed higher inflammation in the lungs of CD4creRorafl/fl mice, after N. brasiliensis infection, indicating a functional role in CD4+ T cells. A similar effect was shown in skin, when Rora was deleted in Tregs only [13]. The effect was demonstrated to be through Treg Tnfrsf25 signaling. Tnfrsf25 has been shown to be required also for Th2 cell function in lung inflammation [37]. Our analysis of the DE genes of activated Th cells, in the lung and colon of worm infected mice, as well as overexpression, and reanalysis of previously published data, confirms the importance of Tnfrsf25 in Tregs, but that its origin may include also other activated T helper cell subsets.
A higher fraction of activated T helper cells (both Th and Tregs) contained Rora, when located in non-lymphoid tissues (lung and siLP) compared to lymph node tissues. The same pattern was previously shown for skin-located versus LN-located Tregs [13]. The enrichment of Rora+/teal cells in peripheral tissues may reflect its correlation with T cell activation and effector function, or alternatively may relate to a role in influencing T cell migration. Our data is also unable to conclusively show if Rora is expressed before or after cytokines are expressed. Further research is required to see how Rora is involved in these regulatory programs.
One theory is that Rora has a complex interaction with circadian rhythm genes. T cells migrate in and out of the lymph node throughout the day [63], and Rora is already known to be involved in circadian rhythm. Our data strengthens this claim as the circadian gene Arntl is DE in both Tregs and other activated CD4+ T cells, in agreement with previous studies [42, 64]. It is possible that Rora is regulated through external signalling. We investigated this with a small screen and found, for example, SDF1a as a negative regulator. This could be interpreted in the context of leukocyte migration during circadian rhythms, as blocking SDF1a disables the circadian migration of leukocytes [63]. A caveat with the circadian rhythm influence is that, if the Rora RNA cycles over the day, so will the Rora protein—but with a delay induced by translation speed. This significantly complicates the interpretation of our data, and possibly helps explain any discrepancies in our RNA-seq vs FACS analysis. Further detailed work is required to account for any circadian effects.
Other theories can be brought forward based on DE genes downstream of Rora. If Sell is a direct downstream gene, it would explain most of the activation. However, other downstream genes such as S100a4, Emb and Cxcr6 have also been linked to activation and cell migration. As Rora is a transcription factor which may have thousands of direct targets, the effect may be through more than one downstream gene, making it hard to pinpoint a single molecular mechanism.
The previous Rora KO sg mice display a range of phenotypes, e.g., reduced body weight gain and neuronal disorders. Assuming the gene regulatory network is preserved in other cell types, our work provides hypothetical mechanisms to explore further in future experiments. Since Rora controls PPARa [65], it may be a possible intermediate, with our data suggesting effects on phosphocholines. In this work we did not investigate the effect of potential Rora ligands but these could further contribute to a feedback loop. As for the neuronal disorders, knowing that the SDF1a-CXCR4 interaction is also important for neuronal guidance [66], SDF1a-mediated regulation of Rora may also explain neuronal impacts.
In conclusion, Rora has been found to be important for several immune cell types, in different contexts. It is involved in multiple pathways and appears to link cell migration with cell cycle, but also affects the outcome of inflammation through, but not limited to, Tnfrsf25. In its natural context, Rora is mainly associated with non-lymphoid tissue and activated T cells. In vitro culture suggests that this is due to environmental stimuli, possibly SDF1a. By generating a large systematic transcriptomic characterization of cytokine effects on T helper cells, we have found several cytokines that may explain be responsible for Rora induction in other tissues. With the large number of conditions tested in this study we are one step closer to understanding the multiple functions of this gene.
Supporting information
S1 File. Collection of supplementary S1–S5 Figs.
https://doi.org/10.1371/journal.pone.0251233.s001
(PDF)
S2 File. RNA-seq count data and condition matrices.
https://doi.org/10.1371/journal.pone.0251233.s002
(ZIP)
S4 File. Biomark single cell RT-qPCR data, raw data and list of TaqMan probes.
https://doi.org/10.1371/journal.pone.0251233.s004
(ZIP)
Acknowledgments
We would like to thank Helen Jolin for help with mice, MRC Ares and Biomed staff, and Bee Ling Ng, Chris Hall, Jennie Graham, Maria Daly, Fan Zhang and Martyn Balmont for help with cell sorting. Ayesha Jinat provided technical assistance with RNA-seq.
References
- 1. Pulendran B, Artis D. New paradigms in type 2 immunity. Science. 2012;337: 431–435. pmid:22837519
- 2. Halim TYF. Group 2 innate lymphoid cells in disease. Int Immunol. 2016;28: 13–22. pmid:26306498
- 3. Gasteiger G, Fan X, Dikiy S, Lee SY, Rudensky AY. Tissue residency of innate lymphoid cells in lymphoid and nonlymphoid organs. Science. 2015;350: 981–985. pmid:26472762
- 4. Jetten AM. Retinoid-related orphan receptors (RORs): critical roles in development, immunity, circadian rhythm, and cellular metabolism. Nucl Recept Signal. 2009;7: e003. pmid:19381306
- 5. Solt LA, Burris TP. Action of RORs and their ligands in (patho)physiology. Trends Endocrinol Metab. 2012;23: 619–627. pmid:22789990
- 6. Jetten AM, Kang HS, Takeda Y. Retinoic acid-related orphan receptors α and γ: key regulators of lipid/glucose metabolism, inflammation, and insulin sensitivity. Front Endocrinol. 2013;4: 1.
- 7. Cook DN, Kang HS, Jetten AM. Retinoic Acid-Related Orphan Receptors (RORs): Regulatory Functions in Immunity, Development, Circadian Rhythm, and Metabolism. Nucl Receptor Res. 2015;2. pmid:26878025
- 8. Halim TYF, MacLaren A, Romanish MT, Gold MJ, McNagny KM, Takei F. Retinoic-acid-receptor-related orphan nuclear receptor alpha is required for natural helper cell development and allergic inflammation. Immunity. 2012;37: 463–474. pmid:22981535
- 9. Wong SH, Walker JA, Jolin HE, Drynan LF, Hams E, Camelo A, et al. Transcription factor RORα is critical for nuocyte development. Nat Immunol. 2012;13: 229–236. pmid:22267218
- 10. Lo BC, Gold MJ, Hughes MR, Antignano F, Valdez Y, Zaph C, et al. The orphan nuclear receptor RORα and group 3 innate lymphoid cells drive fibrosis in a mouse model of Crohn’s disease. Sci Immunol. 2016;1: eaaf8864.
- 11. Han Y-H, Kim H-J, Na H, Nam M-W, Kim J-Y, Kim J-S, et al. RORα Induces KLF4-Mediated M2 Polarization in the Liver Macrophages that Protect against Nonalcoholic Steatohepatitis. Cell Rep. 2017;20: 124–135. pmid:28683306
- 12. Yang XO, Pappu BP, Nurieva R, Akimzhanov A, Kang HS, Chung Y, et al. T helper 17 lineage differentiation is programmed by orphan nuclear receptors ROR alpha and ROR gamma. Immunity. 2008;28: 29–39. pmid:18164222
- 13. Malhotra N, Leyva-Castillo JM, Jadhav U, Barreiro O, Kam C, O’Neill NK, et al. RORα-expressing T regulatory cells restrain allergic skin inflammation. Sci Immunol. 2018;3. pmid:29500225
- 14. Walker JA, Clark PA, Crisp A, Barlow JL, Szeto A, Ferreira ACF, et al. Polychromic Reporter Mice Reveal Unappreciated Innate Lymphoid Cell Progenitor Heterogeneity and Elusive ILC3 Progenitors in Bone Marrow. Immunity. 2019;51: 104–118.e7. pmid:31128961
- 15. Henriksson J, Chen X, Gomes T, Ullah U, Meyer KB, Miragaia R, et al. Genome-wide CRISPR Screens in T Helper Cells Reveal Pervasive Crosstalk between Activation and Differentiation. Cell. 2019;176: 882–896.e18. pmid:30639098
- 16. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15: 550. pmid:25516281
- 17. Weigmann B, Tubbe I, Seidel D, Nicolaev A, Becker C, Neurath MF. Isolation and subsequent analysis of murine lamina propria mononuclear cells from colonic tissue. Nat Protoc. 2007;2: 2307–2311. pmid:17947970
- 18. Tucker MS, Karunaratne LB, Lewis FA, Freitas TC, Liang Y-S. Schistosomiasis. Curr Protoc Immunol. 2001; Chapter 19: Unit 19.1.
- 19. Mann VH, Morales ME, Rinaldi G, Brindley PJ. Culture for genetic manipulation of developmental stages of Schistosoma mansoni. Parasitology. 2010;137: 451–462. pmid:19765348
- 20. Pau G, Fuchs F, Sklyar O, Boutros M, Huber W. EBImage—an R package for image processing with applications to cellular phenotypes. Bioinformatics. 2010;26: 979–981. pmid:20338898
- 21. Bray NL, Pimentel H, Melsted P, Pachter L. Near-optimal probabilistic RNA-seq quantification. Nat Biotechnol. 2016;34: 525–527. pmid:27043002
- 22. Benaglia T, Chauveau D, Hunter D, Young D. mixtools: An R package for analyzing finite mixture models. J Stat Softw. 2009;32: 1–29.
- 23. Patro R, Duggal G, Love MI, Irizarry RA, Kingsford C. Salmon provides fast and bias-aware quantification of transcript expression. Nat Methods. 2017. pmid:28263959
- 24. McCarthy DJ, Campbell KR, Lun ATL, Wills QF. Scater: pre-processing, quality control, normalization and visualization of single-cell RNA-seq data in R. Bioinformatics. 2017;33: 1179–1186. pmid:28088763
- 25. Picelli S, Faridani OR, Björklund ÅK, Winberg G, Sagasser S, Sandberg R. Full-length RNA-seq from single cells using Smart-seq2. Nat Protoc. 2014;9: 171–181. pmid:24385147
- 26. Butler A, Satija R. Integrated analysis of single cell transcriptomic data across conditions, technologies, and species. bioRxiv. 2017. Available: https://www.biorxiv.org/content/early/2017/07/18/164889.abstract
- 27. Stubbington MJT, Lönnberg T, Proserpio V, Clare S, Speak AO, Dougan G, et al. T cell fate and clonality inference from single-cell transcriptomes. Nat Methods. 2016;13: 329–332. pmid:26950746
- 28. Bürglin TR, Henriksson J. FACSanadu: Graphical user interface for rapid visualization and quantification of flow cytometry data. bioRxiv. 2017. p. 201897.
- 29. Folch J, Lees M, Sloane Stanley GH. A simple method for the isolation and purification of total lipides from animal tissues. J Biol Chem. 1957;226: 497–509. pmid:13428781
- 30. Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43: e47. pmid:25605792
- 31. Stubbington MJT, Mahata B, Svensson V, Deonarine A, Nissen JK, Betz AG, et al. An atlas of mouse CD4+ T cell transcriptomes. Biol Direct. 2015;10: 14. pmid:25886751
- 32. Lönnberg T, Svensson V, James KR, Fernandez-Ruiz D, Sebina I, Montandon R, et al. Single-cell RNA-seq and computational analysis using temporal mixture modelling resolves Th1/Tfh fate bifurcation in malaria. Sci Immunol. 2017;2. pmid:28345074
- 33. Miragaia RJ, Gomes T, Chomka A, Jardine L, Riedel A, Hegazy AN, et al. Single-Cell Transcriptomics of Regulatory T Cells Reveals Trajectories of Tissue Adaptation. Immunity. 2019;50: 493–504.e7. pmid:30737144
- 34. Davidson S, Efremova M, Riedel A, Mahata B, Pramanik J, Huuhtanen J, et al. Single-cell RNA sequencing reveals a dynamic stromal niche within the evolving tumour microenvironment. bioRxiv. 2018. p. 467225.
- 35. Nemeth K, Keane-Myers A, Brown JM, Metcalfe DD, Gorham JD, Bundoc VG, et al. Bone marrow stromal cells use TGF-β to suppress allergic responses in a mouse model of ragweed-induced asthma. Proc Natl Acad Sci U S A. 2010;107: 5652–5657. pmid:20231466
- 36. Moon JJ, Chu HH, Pepper M, McSorley SJ, Jameson SC, Kedl RM, et al. Naive CD4(+) T cell frequency varies for different epitopes and predicts repertoire diversity and response magnitude. Immunity. 2007;27: 203–213. pmid:17707129
- 37. Fang L, Adkins B, Deyev V, Podack ER. Essential role of TNF receptor superfamily 25 (TNFRSF25) in the development of allergic lung inflammation. J Exp Med. 2008;205: 1037–1048. pmid:18411341
- 38. Wang J, Fa J, Wang P, Jia X, Peng H, Chen J, et al. NINJ2- A novel regulator of endothelial inflammation and activation. Cell Signal. 2017;35: 231–241. pmid:28431986
- 39. Castro G, Liu X, Ngo K, De Leon-Tabaldo A, Zhao S, Luna-Roman R, et al. RORγt and RORα signature genes in human Th17 cells. PLoS One. 2017;12: e0181868. pmid:28763457
- 40. Miragaia RJ, Gomes T, Chomka A, Jardine L, Riedel A. Single cell transcriptomics of regulatory T cells reveals trajectories of tissue adaptation. bioRxiv. 2017. Available: https://www.biorxiv.org/content/early/2017/11/22/217489.abstract
- 41. Bollinger T, Leutz A, Leliavski A, Skrum L, Kovac J, Bonacina L, et al. Circadian clocks in mouse and human CD4+ T cells. PLoS One. 2011;6: e29801. pmid:22216357
- 42. Akashi M, Takumi T. The orphan nuclear receptor RORα regulates circadian transcription of the mammalian core-clock Bmal1. Nat Struct Mol Biol. 2005;12: 441. pmid:15821743
- 43. Fei F, Qu J, Zhang M, Li Y, Zhang S. S100A4 in cancer progression and metastasis: A systematic review. Oncotarget. 2017;8: 73219–73239. pmid:29069865
- 44. Grum-Schwensen B, Klingelhöfer J, Beck M, Bonefeld CM, Hamerlik P, Guldberg P, et al. S100A4-neutralizing antibody suppresses spontaneous tumor progression, pre-metastatic niche formation and alters T-cell polarization balance. BMC Cancer. 2015;15: 44. pmid:25884510
- 45. Mandai Y, Takahashi D, Hase K, Ebisawa M, Nakagawa T, Sato T, et al. S1711 Distinct Roles for Cxcr6 and Cxcr6− CD4 T Cells in the Pathogenesis of Chronic Colitis. Gastroenterology. 2010;138: S–258.
- 46. Latta M, Mohan K, Issekutz TB. CXCR6 is expressed on T cells in both T helper type 1 (Th1) inflammation and allergen-induced Th2 lung inflammation but is only a weak mediator of chemotaxis. Immunology. 2007;121: 555–564. pmid:17437534
- 47. Bitsch F, Aichholz R, Kallen J, Geisse S, Fournier B, Schlaeppi JM. Identification of natural ligands of retinoic acid receptor-related orphan receptor alpha ligand-binding domain expressed in Sf9 cells—a mass spectrometry approach. Anal Biochem. 2003;323: 139–149. pmid:14622968
- 48. Zhang Z, Burch PE, Cooney AJ, Lanz RB, Pereira FA, Wu J, et al. Genomic analysis of the nuclear receptor family: new insights into structure, regulation, and evolution from the rat genome. Genome Res. 2004;14: 580–590. pmid:15059999
- 49. Harnett W, Harnett MM, Byron O. Structural/functional aspects of ES-62—a secreted immunomodulatory phosphorylcholine-containing filarial nematode glycoprotein. Curr Protein Pept Sci. 2003;4: 59–71. pmid:12570785
- 50. Jung DE, Kim JM, Kim C, Song SY. Embigin is overexpressed in pancreatic ductal adenocarcinoma and regulates cell motility through epithelial to mesenchymal transition via the TGF-β pathway. Mol Carcinog. 2016;55: 633–645. pmid:25773908
- 51. Weaver CT, Hatton RD, Mangan PR, Harrington LE. IL-17 family cytokines and the expanding diversity of effector T cell lineages. Annu Rev Immunol. 2007;25: 821–852. pmid:17201677
- 52. Yasuhara R, Irié T, Suzuki K, Sawada T, Miwa N, Sasaki A, et al. The β-catenin signaling pathway induces aggressive potential in breast cancer by up-regulating the chemokine CCL5. Exp Cell Res. 2015;338: 22–31. pmid:26363360
- 53. Khaled AR, Durum SK. Lymphocide: cytokines and the control of lymphoid homeostasis. Nat Rev Immunol. 2002;2: 817–830. pmid:12415306
- 54. Lodolce JP, Boone DL, Chai S, Swain RE, Dassopoulos T, Trettin S, et al. IL-15 Receptor Maintains Lymphoid Homeostasis by Supporting Lymphocyte Homing and Proliferation. Immunity. 1998;9: 669–676. pmid:9846488
- 55. Hunter CA, Jones SA. IL-6 as a keystone cytokine in health and disease. Nat Immunol. 2015;16: 448–457. pmid:25898198
- 56. Peine M, Marek RM, Löhning M. IL-33 in T Cell Differentiation, Function, and Immune Homeostasis. Trends Immunol. 2016;37: 321–333. pmid:27055914
- 57. Liu Y, Cai Y, Liu L, Wu Y, Xiong X. Crucial biological functions of CCL7 in cancer. PeerJ. 2018;6: e4928. pmid:29915688
- 58. Qin X-J, Shi H-Z, Deng J-M, Liang Q-L, Jiang J, Ye Z-J. CCL22 recruits CD4-positive CD25-positive regulatory T cells into malignant pleural effusion. Clin Cancer Res. 2009;15: 2231–2237. pmid:19318474
- 59. Vulcano M, Albanesi C, Stoppacciaro A, Bagnati R, D’Amico G, Struyf S, et al. Dendritic cells as a major source of macrophage-derived chemokine/CCL22 in vitro and in vivo. Eur J Immunol. 2001;31: 812–822. pmid:11241286
- 60. Campbell JJ, OConnell DJ, Wurbel M-A. Cutting Edge: Chemokine receptor CCR4 is necessary for antigen-driven cutaneous accumulation of CD4 T cells under physiological conditions. J Immunol. 2007;178: 3358–3362. pmid:17339428
- 61. Arojo OA, Ouyang X, Liu D, Meng T, Kaech SM, Pereira JP, et al. Active mTORC2 Signaling in Naive T Cells Suppresses Bone Marrow Homing by Inhibiting CXCR4 Expression. J Immunol. 2018. pmid:29934471
- 62. Halim TYF, Steer CA, Mathä L, Gold MJ, Martinez-Gonzalez I, McNagny KM, et al. Group 2 innate lymphoid cells are critical for the initiation of adaptive T helper 2 cell-mediated allergic lung inflammation. Immunity. 2014;40: 425–435. pmid:24613091
- 63. He W, Holtkamp S, Hergenhan SM, Kraus K, de Juan A, Weber J, et al. Circadian Expression of Migratory Factors Establishes Lineage-Specific Signatures that Guide the Homing of Leukocyte Subsets to Tissues. Immunity. 2018. pmid:30527911
- 64. Druzd D, Matveeva O, Ince L, Harrison U, He W, Schmal C, et al. Lymphocyte Circadian Clocks Control Lymph Node Trafficking and Adaptive Immune Responses. Immunity. 2017;46: 120–132. pmid:28087238
- 65. Kim K, Boo K, Yu YS, Oh SK, Kim H, Jeon Y, et al. RORα controls hepatic lipid homeostasis via negative regulation of PPARγ transcriptional network. Nat Commun. 2017;8: 162. pmid:28757615
- 66. Li M, Ransohoff RM. Multiple roles of chemokine CXCL12 in the central nervous system: a migration from immunology to neurobiology. Prog Neurobiol. 2008;84: 116–131. pmid:18177992