108–126 Nucleic Acids Research, 2022, Vol. 50, No. 1 Published online 10 December 2021 https://doi.org/10.1093/nar/gkab1182 Coordinated glucocorticoid receptor and MAFB action induces tolerogenesis and epigenome remodeling in dendritic cells Octavio Morante-Palacios1,2 , Laura Ciudad1 , Raphael Micheroli3 , Carlos de la Calle-Fabregat1 , Tianlu Li1 , Gisela Barbisan1 , Miranda Houtman3 , Sam G. Edalat3 , Mojca Frank-Bertoncelj3 , Caroline Ospelt3 and Esteban Ballestar 1,2,* 1 Epigenetics and Immune Disease Group, Josep Carreras Research Institute (IJC), 08916 Badalona, Barcelona, Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 Spain, 2 Germans Trias i Pujol Research Institute (IGTP), 08916 Badalona, Barcelona, Spain and 3 Center of Experimental Rheumatology, Department of Rheumatology, University Hospital Zurich, University of Zurich, Zurich, Switzerland Received September 20, 2021; Revised November 08, 2021; Editorial Decision November 10, 2021; Accepted November 16, 2021 ABSTRACT INTRODUCTION Glucocorticoids (GCs) exert potent anti- Dendritic cells (DCs) are a heterogeneous group of innate inflammatory effects in immune cells through immune cells with a central role not only in the response the glucocorticoid receptor (GR). Dendritic cells to threats, but also in the regulation of inflammatory re- (DCs), central actors for coordinating immune sponses and the induction of immune tolerance, which is responses, acquire tolerogenic properties in re- defined as the reduced responsiveness of the immune sys- tem to molecules with the potential to produce an immune sponse to GCs. Tolerogenic DCs (tolDCs) have response (1). Several types of human DCs exist in vivo, in- emerged as a potential treatment for various in- cluding conventional DCs (cDCs) and plasmacytoid DCs flammatory diseases. To date, the underlying cell (pDCs). Moreover, monocytes (MOs) can extravasate to type-specific regulatory mechanisms orchestrating tissues and differentiate to monocyte-derived macrophages GC-mediated acquisition of immunosuppressive (moMACs) or monocyte-derived dendritic cells (moDCs) properties remain poorly understood. In this study, (2,3). MOs cultured with GM-CSF and IL-4 are a widely we investigated the transcriptomic and epigenomic used in vitro model of DCs (4), with gene expression simi- remodeling associated with differentiation to DCs in larities with blood cDCs and in vivo moDCs (5). the presence of GCs. Our analysis demonstrates a Epigenetic determinants and transcription factors (TFs) major role of MAFB in this process, in synergy with critical to DC differentiation and function have been ex- GR. GR and MAFB both interact with methylcytosine tensively studied. In particular, STAT6 (6), early growth response 2 (EGR2) (7) and interferon regulatory factor 4 dioxygenase TET2 and bind to genomic loci that (IRF4) (5) have been linked to DC fate determination and undergo specific demethylation in tolDCs. We also specific DNA demethylation events through recruitment show that the role of MAFB is more extensive, of Ten-eleven translocation methylcytosine dioxygenase 2 binding to thousands of genomic loci in tolDCs. (TET2), which is the most expressed TET enzyme in MOs Finally, MAFB knockdown erases the tolerogenic (8). Moreover, aryl hydrocarbon receptor (AHR) can act properties of tolDCs and reverts the specific DNA as a molecular switch that enables monocyte differentiation demethylation and gene upregulation. The preemi- to moDC, via IRF4, whereas MAFB determines monocyte nent role of MAFB is also demonstrated in vivo for differentiation to moMAC (5,9,10). myeloid cells from synovium in rheumatoid arthritis DCs can acquire tolerogenic functions in vivo and in vitro following GC treatment. Our results imply that, once in response to several stimuli, including interleukin (IL)- 10, vitamin D3, rapamycin, and glucocorticoids (GCs) (1). directly activated by GR, MAFB plays a critical role In particular, DCs differentiated from MOs in vitro with in orchestrating the epigenomic and transcriptomic GM-CSF, IL-4 and GCs can suppress T cell proliferation remodeling that define the tolerogenic phenotype. in vitro, and display a high level of production of IL-10 and low levels of TNF␣ and IL-12p70 (11,12). Tolerogenic DCs * To whom correspondence should be addressed. Tel: +34 935572800; Fax: +34 934651472; Email: eballestar@carrerasresearch.org C The Author(s) 2021. Published by Oxford University Press on behalf of Nucleic Acids Research. This is an Open Access article distributed under the terms of the Creative Commons Attribution-NonCommercial License (http://creativecommons.org/licenses/by-nc/4.0/), which permits non-commercial re-use, distribution, and reproduction in any medium, provided the original work is properly cited. For commercial re-use, please contact journals.permissions@oup.com Nucleic Acids Research, 2022, Vol. 50, No. 1 109 (tolDCs), generated as indicated above, can be useful as a dexamethasone (Sigma-Aldrich) was also added to the treatment for autoimmune diseases. Several clinical trials medium. have yielded satisfactory results in diseases such as rheuma- For cell stimulation, LPS (10 ng/ml) and IFN␥ (20 toid arthritis (13,14), multiple sclerosis (15) and Crohn’s dis- ng/mL) were added to cell culture at day 5, 24 h before cell ease (16). collection. GCs are a family of steroid hormones that are ligands of the glucocorticoid receptor (GR), a nuclear receptor ex- CD8+ proliferation suppression assay pressed in most cell types that can trigger the expression Allogenic CD8+ were isolated by negative selection using of anti-inflammatory genes through direct DNA binding. Dynabeads Untouched Human CD8 T Cells Kit (Invitro- Moreover, GR also represses the action of inflammatory- gen) and labeled with carboxyfluorescein succinimidyl es- related TFs, such as the NF-␬B and AP-1 families, via ter (CFSE) CellTrace™ (Invitrogen), in accordance with the protein–protein interactions, in a process called transrepres- manufacturer’s instructions. sion (17). GR has been related to chromatin remodelers, in- Purified CD8+ cells were seeded in 96-well plates (200 000 cluding EP300 and BRG1 (18,19). The mechanisms under- Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 cells per well), and monocyte-derived cells (DCs or tolDCs) lying cell type-specific programs induced by GR upon lig- were added at different cell ratios (1:2, 1:3, 1:4 and 1:8). and binding, as well as the participation of TFs and epige- To stimulate CD8+ cells, 5 ␮l/ml of anti-CD3/CD28 Dyn- netic enzymes, remain to be fully determined, which is par- abeads (Invitrogen) were added to each well, except for the ticularly relevant in the case of the various innate immune negative control. cell types. For this experiment, cells were harvested in RPMI con- In this work, we have studied the transcriptional and epi- taining 10% fetal bovine serum (Gibco, ThermoFisher) genetic remodeling associated with tolDC differentiation, and 100 units/ml penicillin/streptomycin (Gibco, Ther- and identified a major role for MAFB in this process, in moFisher). Cells were cultured for 5 days, and the medium synergy with GR. We have shown that GR binds both the was changed on day 3. promoter and enhancer regions associated with MAFB, which is quickly upregulated and binds thousands of ge- Quantification of cytokine production nomic sites, correlating with widespread DNA demethyla- tion and gene upregulation. We demonstrate how MAFB is Cell culture supernatants were collected after 6 days and di- crucial for the acquisition of the transcriptomic and epige- luted appropriately. Enzyme-linked immunosorbent assays nomic remodeling that gives rise to the tolerogenic pheno- (ELISA) were performed, following the manufacturer’s in- type. This is achieved through a differentiation switch to structions: Human IL-10, Human IL-12p70 and Human macrophage-like cells with tolerogenic properties. The ma- TNF␣ from BioLegend, and Human IL-1␤ from Ther- jor role of MAFB in activating a tolerogenic expression moFisher. profile is also demonstrated in monocyte-derived cells from rheumatoid arthritis joints treated with GCs, in which an Flow cytometry expansion of cells with a transcriptomic signature similar To study cell-surface markers, cells were collected us- to MAFB-dependent tolDCs is shown. ing Versene, a non-enzymatic dissociation buffer (Ther- moFisher). Cells were resuspended in the staining buffer MATERIALS AND METHODS (PBS with 4% fetal bovine serum and 2 mM ethylenedi- CD14+ monocyte purification and culture aminetetraacetic acid (EDTA)). Cells were then incubated in ice with Fc block reagent (Miltenyi Biotec) for 10 min- Buffy coats were obtained from anonymous donors via utes, and stained with the viability dye LIVE/DEAD™ Fix- the Catalan Blood and Tissue Bank (CBTB). The CBTB able Violet (ThermoFisher), following the manufacturer’s follows the principles of the World Medical Association protocol. (WMA) Declaration of Helsinki. Before providing blood Cells were then stained to study the proteins of inter- samples, all donors received detailed oral and written in- est, using the following antibodies: CD16 (APC) (#130-113- formation and signed a consent form at the CBTB. 389, Miltenyi Biotec), CD14 (APC) (#130-110-520, Mil- PBMCs were isolated by density-gradient centrifuga- tenyi Biotec), CD163 (FITC) (#33618, BioLegend), CD1a tion using lymphocyte-isolation solution (Rafer). Pure MOs (PE) (#300106, BioLegend), CD80 (PE) (#H12208P, eBio- were then isolated from PBMCs by positive selection with Science), CD83 (APC) (#130-110-504, Miltenyi Biotec), magnetic CD14 MicroBeads (Miltenyi Biotec). Purity was CD86 (APC) (#130-113-569, Miltenyi Biotec), HLA-DR verified by flow cytometry, which yielded more than 90% of (PE) (#12-9956-42, eBioScience). CD14+ cells. After staining, cells were fixed with PBS + 4% MOs were resuspended in Roswell Park Memorial Insti- paraformaldehyde (Electron Microscopy Sciences) and tute (RPMI) Medium 1640 + GlutaMAX™ (Gibco, Ther- analyzed within 2 days using a BD FACSCanto™ II Cell moFisher) and immediately added to cell culture plates. Analyzer (BD Biosciences). Data were analyzed with the After 20 min, monocytes were attached to the cell culture FlowJo v10 software. plates, and the medium was changed with RPMI contain- ing 10% fetal bovine serum (Gibco, ThermoFisher), 100 Genomic DNA and total RNA extraction units/mL penicillin/streptomycin (Gibco, ThermoFisher), 10 ng/ml human GM-CSF (PeproTech) and 10 ng/ml hu- Genomic DNA and total RNA were extracted using the man IL-4 (PeproTech). In the case of tolDCs, 100 nM Maxwell RSC Cultured Cells DNA kit (Promega) and the 110 Nucleic Acids Research, 2022, Vol. 50, No. 1 Maxwell RSC simplyRNA Cells kit (Promega), respec- tion, covering 99% of the reference sequence (RefSeq) genes. tively, following the manufacturer’s instructions. Raw files (IDAT files) were provided for the Josep Carreras Research Institute Genomics Platform (Barcelona). Quality control and analysis of EPIC arrays were per- Gene expression microarrays formed using ShinyÉPICo (21), a graphical pipeline that RNA samples were processed in the Genomics Platform of uses minfi (22) for normalization, and limma (20) for dif- the Vall d’Hebron Research Institute (Barcelona). After en- ferentially methylated positions analysis. CpH and SNP suring the RNA quality using Bioanalyzer, samples were hy- loci were removed and the Noob + Quantile normaliza- bridized in Clariom™ S microarrays (ThermoFisher). tion method was used. Donor information was used as a Raw microarray data (CEL files) were analyzed using the covariate, and Trend and Robust options were enabled for oligo and limma packages of the Bioconductor project (20). the eBayes moderated t-test analysis. CpGs were considered First, raw data were normalized using the Robust Multichip differentially methylated when the absolute differential of Average algorithm (RMA), which is included as a function methylation was >20% and the FDR was <0.05. in the oligo package. We then performed an independent fil- Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 tering step, removing probes with fewer than three samples over the 75% percentile of the negative control probes. After Immunofluorescence collapsing the remaining probes by gene, using the aggre- gate and mean functions, we built a limma linear model us- Cells were fixed with PBS + 4% paraformaldehyde for 20 ing the condition and donor information, with the formula min and permeabilized with PBS + Triton X-100 0.5% ‘∼0 + condition + donor’. The eBayes function in limma for 10 min. Coverslips were washed twice with PBS and was then used in each pairwise comparison to obtain the blocked with PBS + 4% BSA for 1 h. Anti-MAFB antibody FDR and logFC of each gene. Genes were considered to be HPA005653 (Sigma-Aldrich) was diluted 1:100 in PBS and differentially expressed when the FDR was <0.05 and the incubated overnight with the samples in a humidity cham- absolute logFC was >0.5. ber. Then, cells were incubated with anti-rabbit Alexa Fluor 647 (ThermoFisher) in blocking solution (PBS + BSA 4% + 0.025% Tween 20) for 1 h. After four washes with PBS, Bisulfite pyrosequencing cells were stained with Alexa Fluor™ 488 Phalloidin (Ther- 500 ng of genomic DNA was converted using the EZ DNA moFisher) 1/200 and DAPI 2 ␮g/ml. Vectashield (Vector Methylation Gold kit (Zymo Research). PCR was per- Laboratories) was used for the final sample preparation in formed using the bisulfite-converted DNA as input and slides. Images were obtained with a Leica TCS-SL confocal primers designed for each amplicon (Supplementary Table microscope. S5). These primers were designed using the PyroMark As- say Design 2.0 software (Qiagen). PCR amplicons were py- rosequenced using the PyroMark Q48 system and analyzed Western blotting with PyroMark Q48 Autoprep software. Cytoplasmic and nuclear protein fractions were obtained using hypotonic lysis buffer (Buffer A; 10 mM Tris pH 7.9, 1.5 mM MgCl2 , 10 mM KCl supplemented with pro- Real-time quantitative reverse-transcribed polymerase chain tease inhibitor cocktail (Complete, Roche) and phosphatase reaction (qRT-PCR) inhibitor cocktail (PhosSTOP, Roche) to lyse the plasma 300 ng of total RNA were reverse-transcribed to cDNA membrane. Cells were visualized under the microscope to with Transcriptor First Strand cDNA Synthesis Kit ensure correct cell lysis. Nuclear pellets were resuspended in (Roche) following manufacturer’s instructions. qRT-PCR Laemmli 1× loading buffer. For whole-cell protein extract, was performed in technical triplicates for each biologi- cell pellets were directly resuspended in Laemmli 1× load- cal replicate, using LightCycler® 480 SYBR Green Mix ing buffer. (Roche), and 7.5 ng of cDNA per reaction. The average Proteins were separated by SDS-PAGE electrophoresis. value from each technical replicate was obtained. Then, the Immunoblotting was performed on polyvinylidene diflu- standard double-delta Ct method was used to determine the oride (PVDF) membranes following standard procedures. relative quantities of target genes, and values were normal- Membranes were blocked with 5% Difco™ Skim Milk (BD ized against the control genes RPL38 and HPRT1. Custom Biosciences) and blotted with primary antibodies. After primers were designed to analyze genes of interest (Supple- overnight incubation, membranes were washed three times mentary Table S5) for 10 min with TBS-T (50 mM Tris, 150 mM NaCl, 0.1% Tween-20) and incubated for 1 h with HRP-conjugated mouse or rabbit secondary antibody solutions (Thermo DNA methylation profiling Fisher) diluted in 5% milk (diluted 1/10000). Finally, pro- 500 ng of genomic DNA was converted using the EZ DNA teins were detected by chemiluminescence using Western- Methylation Gold kit (Zymo Research). Infinium Methyla- Bright™ ECL (Advansta). The following antibodies were tionEPIC BeadChip (Illumina) arrays were used to analyze used: Anti-MAFB (HPA005653, Sigma-Aldrich), Anti-GR DNA methylation, following the manufacturer’s instruc- (C15200010-50, Diagenode), Anti-GAPDH (2275-PC-100, tions. This platform allows around 850 000 methylation sites Trevigen), Anti-Lamin B1 (ab229025, Abcam), Anti-TET2 per sample to be interrogated at single-nucleotide resolu- (C15200179, Diagenode). Nucleic Acids Research, 2022, Vol. 50, No. 1 111 Co-immunoprecipitation (Co-IP) Chromatin immunoprecipitation Co-IP assays were performed using tolDCs differentiated After 4 or 24 h of cell culture, DCs and tolDCs, in bio- from CD14 + monocytes for 24 h. Cell extracts were pre- logical duplicates, were fixed with Pierce™ fresh methanol- pared in lysis buffer [50 mM Tris–HCl, pH 7.5, 1 mM free formaldehyde (ThermoFisher) for 15 min and prepared EDTA, 150 mM NaCl, 1% Triton-X-100, protease in- for sonication with the truChIP Chromatin Shearing Kit hibitor cocktail (cOmplete™, Merck)] with correspond- (Covaris), following the manufacturer’s instructions. Chro- ing units of Benzonase (Sigma) and incubated at 4◦ C for matin was sonicated for 18 min with the Covaris M220 in 4 h. 100 ␮l of supernatant was saved as input and di- 1 ml milliTubes (Covaris). Size distribution of the sonicated luted with 2× Laemmli sample buffer (5× SDS, 20% glyc- chromatin was checked by electrophoresis to ensure appro- erol, 1 M Tris–HCl (pH 8.1)). Supernatants were first pre- priate sonication, with a size of around 200 bp. cleared with PureProteome™ Protein A/G agarose suspen- Magna Beads Protein A + G (Millipore) was blocked sion (Merck Millipore) for 1 h. The lysate was then in- with PBS + BSA (5 mg/ml) for 1 h. Chromatin was pre- cubated overnight at 4◦ C with respective crosslinked pri- cleared with 25 ␮l of beads for 1.5 h and 10 ␮g of chromatin Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 mary antibodies. Cross-linking was performed in 20 mM were incubated overnight with each antibody: 10 ␮l Anti- dimethyl pimelimidate (DMP) (Pierce, ThermoFisher Sci- MAFB antibody HPA005653 (Sigma-Aldrich) and 5 ␮l entific, MA, USA) dissolved in 0.2 M sodium borate (pH Anti-GR antibody C15200010-50 (Diagenode), in a buffer 9.0). Subsequently, the beads were quenched with 0.2M with 1% Triton X-100 and 150 mM NaCl. of ethanolamine (pH 8.0) and resuspended at 4ºC in PBS Three washes were performed with the Low Salt Wash until use. Beads were then washed three times with lysis Buffer (0.1% SDS, 1% Triton X-100, 2 mM EDTA, pH 8.0, buffer at 4ºC. Samples were eluted by acidification using a 20 mM Tris–HCl pH 8.0, 150 mM NaCl), the High Salt buffer containing 0.2 M glycine (pH 2.3) and diluted with Wash Buffer (0.1%SDS, 1% Triton X-100, 2 mM EDTA, 2× Laemmli. Samples and inputs were denatured at 95◦ C pH 8.0, 20 mM Tris–HCl pH 8, 500 mM NaCl), and the in the presence of 1% ␤-mercaptoethanol. Anti-MAFB LiCl Wash Buffer (0.25 M LiCl, 1% Nonidet P-40, 1% De- HPA005653 (Sigma-Aldrich), Anti-GR C15200010-50, and oxycholate, 1 mM EDTA pH 8, 10 mM Tris–HCl), followed Anti-TET2 antibody ab124297 (Abcam) were used for Co- by a final wash with TE buffer (pH 8.0, 10 mM Tris–HCl, IP. 1 mM EDTA). Chromatin was eluted for 45 min at 65ºC with 100 ␮l of elution buffer (10 mM Tris–Cl, 1 mM EDTA, 1% SDS) and decrosslinked by adding 5 ␮l 5M NaCl and Transfection of primary human monocytes 5 ␮l 1 M NaHCO3 (2 h at 65ºC). Next, 1 ␮l of 10 mg/ml We used the #si19279 MAFB Silencer Select siRNA (Ag- proteinase K (Invitrogen) was added and the samples were ilent, ThermoFisher) to perform knockdown experiments incubated at 37ºC for 1 h. in MOs, using the Silencer Select Negative Control #1 (Ag- For DNA purification, iPure kit v2 (Diagenode) was ilent, ThermoFisher) as control. CD14+ MOs were cultured used, following the manufacturer’s instructions. in 12-well plates (1 million cells/well) and transfected with siRNAs (100 nM) using Lipofectamine 3000 Reagent (1.5 ␮l/well), following the manufacturer’s protocol. After 6 h, Chromatin immunoprecipitation sequencing (ChIP-seq) the medium was changed and DCs/tolDCs were cultured as analysis described in ‘CD14+ monocytes purification and culture’. ChIP-seq inputs and immunoprecipitated DNA from bio- To test transfection efficiency, siGLO Green Transfection logical duplicate samples were used to generate ChIP-seq li- Indicator (Horizon) was used, at the same concentration, braries in the Centre for Genomic Regulation (Barcelona), and checked by flow cytometry 24 h after transfection. using the TruSeq ChIP Library Preparation Kit (Illumina). Quality control of the libraries was performed using 2100 Bioanalyzer (Agilent). Libraries were sequenced with an Il- RNA-seq lumina HiSeq 2500, in 50-bp single-end, yielding between RNA-seq libraries of transfected DC/tolDCs were gener- 25 and 40 million reads per sample. ated and sequenced by BGI Genomics (Hong Kong), in Potential adapter contamination was trimmed from the 100-bp paired-end, with the DNBseq platform. More than raw reads using Cutadapt (https://doi.org/10.14806/ej.17.1. 30 million reads were obtained for each sample. Fastq files 200). Reads were aligned to the GRCh38 human genome were aligned to the hg38 transcriptome using HISAT2(23) assembly using the Burrows-Wheeler Alignment (BWA)- with standard options. Reads mapped in proper pair and MEM algorithm (27). Aligned reads were filtered by primary alignments were selected with SAMtools (24). MAPQ, removing alignments with MAPQ < 30, using Reads were assigned to genes with featureCounts (25). SAMtools (24). Aligned reads overlapping with the EN- Differentially expressed genes were detected with DE- CODE blacklist were also removed (28). Seq2 (26). The donor was used as a covariate in the Quality control of ChIP-seq data was performed using model. The Ashr shrinkage algorithm was applied and only the SPP package(29). The relative strand cross-correlation protein-coding genes with an absolute logFC >0.5 and an coefficient (RSC) was greater than or equal to 2 in all the FDR <0.05 were selected as differentially expressed. For immunoprecipitated samples. representation purposes, Variance Stabilizing Transforma- Bigwig files were generated for visualization, using the tion (VST) values and normalized counts provided by DE- bamCoverage function in the Deeptools package (30), Seq2 were used. with the bins per million mapped reads (BPM) method 112 Nucleic Acids Research, 2022, Vol. 50, No. 1 and ‘–binSize 20 –extendReads 150 –smoothLength 60 – (GSE40484, GSE117946 and GSE99056). Our expression centerReads’ options. Wiggletools was used to aggregate dataset and public expression data were combined by gene the bigwig duplicates (31). symbol and batch-corrected using ComBat (42). The 1000 MACS2 software with ‘–nomodel –extsize 200’ options most variable genes were used to plot the data in a T- was used for peak calling (32). Biological duplicates were distributed Stochastic Neighbor Embedding (t-SNE) rep- aggregated using the MSPC algorithm (33) with the options resentation, using the Rtsne package (https://www.jmlr.org/ ‘-r Biological -w 1E-6 -s 1E-12’. This approach allows us to papers/v9/vandermaaten08a.html). obtain more robust peaks, using the information from both biological duplicates. The resulting consensus peaks were used for the downstream analysis and figure panels. Synovial tissue extraction and processing Synovial biopsy and tissue analysis was approved by the Data analysis and representation Kantonal Ethic Commission Zurich, Switzerland (refs: 2019-00115 and 2016-02014). All patients signed informed Statistical analyses were performed in R 4.0. Gene expres- consent forms. Patient’s characteristics were as follows: Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 sion and DNA methylation heatmaps were created with the Rheumatoid Synovitis heatmap.2 function of the gplots package. The findMotif- Gender Site Diagnosis factor score Histology Therapy sGenome.pl function of HOMER (Hypergeometric Opti- 28 f MCP RA Negative 5 Lympho- Glucocorticoids mization of Motif EnRichment) was used to analyze known myeloid oral 10 mg/day and de novo motif enrichment. For ChIP-seq peaks, the 50 f Wrist RA Negative 8 Lympho- None parameter ‘-size 50’ was used, whereas the parameters ‘- myeloid size 200 -cpg) were used for methylation data. All EPIC array CpG coordinates were also used as background for One part of the synovial biopsies was embedded in the methylation data. GREAT software was used to cal- paraffin and histologically analyzed as previously described culate CpG-associated genes and gene ontology (GO) en- (43,44). For scRNA-seq synovial biopsies were mechani- richment (34). Gene set enrichment analysis (GSEA) and cally minced and enzymatically digested with 100 ug/ml GO enrichment of gene expression data were performed us- Liberase TL (Roche) and 100 ug/ml DNAseI (Roche) in ing the clusterProfiler package (35). ChIP-seq peaks files of pre-warmed RPMI 1640 containing 2 mM glutamine and histone marks from MO and DCs were downloaded from 25 mM HEPESGibco for 30 min at 37◦ C. The digestion was the BLUEprint webpage (http://dcc.blueprint-epigenome. stopped by adding 20% (v/v) fetal calf serum (FCS) and the eu). Consensus peaks of the different replicates were ob- dissociated tissue was filtered through a 40 um cell strainer. tained with the MSPC algorithm, using the options ‘-r Bio- Red blood cells were lysed with RBC lysis buffer (Milteny logical -w 1E-4 -s 1E-8 -c 3’. Biotec). Cell viability was analyzed with a LUNA-FL dual The chromatin state learning model for CD14 + mono- fluorescence cell counter (88% and 90%). cytes was downloaded from the Roadmap Epigenomics Project webpage, and chromatin state enrichments were cal- Single-cell RNA-seq processing and analysis culated using Fisher’s exact test. ChIP-seq overlaps and Venn diagrams were generated Dissociated cells were processed with the Chromium Sin- with the ChIPpeakAnno package (36). Genomic track plots gle Cell 3 GEM, Library & Gel Bead Kit v3 or v3.1 ac- were created with pyGenomeTracks (37). Methylated CpG cording to the manufacturer’s protocol (10x Genomics). Li- set enrichment analysis (mCSEA) (38) was used to calculate braries were sequenced on the NovaSeq 6000 (Illumina). CpG-set-specific DNA methylation modifications. Public Raw sequencing data was processed with cell ranger (v6.0.0, GR ChIP-seq datasets were extracted from the Remap2020 10× Genomics) using mkfastq and count with default set- database (39), excluding samples without glucocorticoid tings and the provided human reference (GRCh38-2020-A). treatment. Peak intersects were calculated using bedtools (40). Public peak-callings of histone marks were extracted Targeted cell number Estimated number of cells Reads/cell from the Blueprint database (41) and replicates were aggre- 28 10 000 4389 44 396 gated using the MSPC algorithm. Public DNAse-seq and 50 6000 5760 53 176 ATAC-seq bigwigs were aggregated using wiggletools. We used Seurat (4.0.3) to analyze scRNA-seq data. We first read and merge data from patients, and we filtered Student’s paired t test cells with counts lower than 200 or higher than 30 000 and Statistical analyses involved Student’s paired-samples t with RNA features lower than 100 or higher than 5000. tests, with which the means of matched pairs of groups were We also filtered cells with mitochondrial or ribosomal fea- compared, except where indicated otherwise. The levels of tures higher than 25%. We then integrated samples before significance were: ns, P > 0.05; *P < 0.05; **P < 0.01; clustering, using the standard Seurat pipeline. First, fea- ***P < 0.001. tures were normalized and variable features were found for each sample, independently. Next, functions FindInte- grationAnchors and IntegrateData were used to integrate Gene expression public data processing the data. The functions ScaleData (using RNA counts To compare gene expression with public data, expression and mitochondrial percentage as variables to regress), array matrices were downloaded from the GEO database RunPCA, FindNeigbors, FindClusters and RunUMAP Nucleic Acids Research, 2022, Vol. 50, No. 1 113 were used sequentially to perform the UMAP representa- ZBTB16 and FKBP5 in tolDCs, as expected (Supplemen- tion, with 30 dimensions. After selecting only monocyte- tary Figure S1G). derived clusters, scaling and representation was performed To predict potential additional TFs involved in tolDC again, using the same functions, with resolution 0.15 and transcriptome acquisition, we performed a Discrimi- 7 dimensions. nant Regulon Expression Analysis (DoRothEA) (47). DoRothEA uses a dataset of TF regulons, a collection of RESULTS TF - target genes interactions obtained from different types of evidence such as literature curated resources and ChIP- Dexamethasone modulates dendritic cell differentiation to a seq peaks and it calculates enrichment of these regulons in tolerogenic phenotype and drives transcriptome remodeling the condition of interest. This analysis revealed that eight associated with MAFB distinct TF regulons were enriched in tolDCs, although only To investigate the mechanisms underlying the glucocorti- two of these were associated with concomitant upregulation coid (GC)-mediated phenotypic remodeling of dendritic of their coding genes: MAFB and MYC (Figure 1G). MYC cells (DCs), as well as the involvement of myeloid-specific has been previously identified as a transcriptional regula- Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 TFs, monocytes (MOs) isolated from peripheral blood of tor in tolDCs (48). To assess the possible role of MAFB healthy donors were differentiated in vitro to DCs and in tolDC gene expression remodeling, we obtained a public tolerogenic DCs (tolDCs) for 6 days using GM-CSF and dataset of genes downregulated after the treatment with an IL-4 in the absence and presence of a GR ligand (dexam- siRNA that targets MAFB (9) in M2 macrophages. These ethasone), respectively (Figure 1A). CD8 + T cell prolifer- genes, positively regulated by MAFB, were more strongly ation assays in co-culture with both DCs and tolDCs re- expressed in tolDCs than in DCs, as can be observed in the vealed the immunosuppressive properties of the latter (Fig- GSEA (Supplementary Figure S1H). In addition, tolDCs ure 1B). ELISA assays also showed that LPS + IFN␥ - presented a transcriptome similar to those of in vitro M2 stimulated tolDCs produced higher levels of IL-10 and macrophages and in vivo mo-MACs (Figure 1H). MAFB smaller amounts of TNF␣, IL-12p70 and IL-1␤ than DCs has previously been identified as a key mediator in the dif- (Figure 1C). Moreover, in the steady-state, higher IL-10 and ferentiation of both cell types (5,9). lower IL-1␤ levels of production were observed in tolDCs, whereas TNF␣ and IL-12p70 were undetectable (Supple- Dexamethasone mediates the acquisition of a specific DNA mentary Figure S1A). tolDCs presented lower levels of ex- methylation pattern which inversely correlates with gene ex- pression of the costimulatory molecule CD86 and the mat- pression during tolDC differentiation uration marker CD83 (Figure 1D). Moreover, tolDCs pre- sented higher levels of CD14 and of the macrophage mark- In parallel, we obtained genome-wide DNA methylation ers CD163 and CD16 (Supplementary Figure S1B). profiles of MOs, DCs and tolDCs after 6 days of differ- We then profiled the transcriptome of DCs, tolDCs and entiation. Principal component analysis (PCA) of all dif- MOs. 981 genes were induced and 919 genes were re- ferentially methylated CpGs between groups showed non- pressed in tolDCs in comparison with DCs (FDR < 0.05, overlapping clustering of MOs, DCs and tolDCs along- logFC > 0.5) (Figure 1E and Supplementary Table S1). The side principal component 1 (Supplementary Figure S2A). transcriptomes of both tolDCs and DCs were notably dis- Demethylation occurring during MO-to-DC differentiation similar to that of MOs (Supplementary Figure S1C). was broadly inhibited in tolDCs (Supplementary Figure The Gene Ontology (GO) over-represented categories in S2B). tolDC-upregulated genes, including terms such as ‘negative The comparison of DCs and tolDCs (FDR < 0.05 and regulation of IL-12 production’, ‘regulation of complement absolute ß > 0.2) in relation to MOs revealed two clus- activation’ and ‘regulation of inflammatory response’ (Fig- ters of CpG sites: a group of CpGs that underwent spe- ure 1F). In tolDC-downregulated genes, terms such as ‘anti- cific demethylation in DCs and that was blocked in tolDCs gen processing and presentation of endogenous antigen’, (C1, 1353 CpGs), and a second group that was specifically ‘adaptive immune response’, ‘positive regulation of leuko- demethylated in tolDCs (C2, 411 CpGs) (Figure 2A and cyte activation’ and ‘response to interferon-gamma’ were Supplementary Table S2). over-represented (Supplementary Figure S1D). Both CpG clusters were enriched in introns and depleted tolDCs were enriched in several gene sets, including: in promoters (Figure 2B). CpGs were generally located far ‘mo-MAC signature’ (genes upregulated in mo-MACs in from CpG islands (Open Sea) (Figure 2C). The two clus- comparison with mo-DCs) (5), genes upregulated in M- ters were also enriched in enhancers and regions close to CSF macrophages relative to GM-CSF macrophages (45), active transcription start sites (TSSs) (Supplementary Fig- and genes downregulated in response to GM-CSF and IL- ure S2C). Looking at the enrichment in active enhancer hi- 4(46) (Supplementary Figure S1E). On the other hand, stone marks (H3K4me1 and H3K27ac) in DCs and MOs, tolDCs were depleted in inflammation-related gene sets an increase in the signal of these marks was noted in MO- such as Interferon-alpha response, TNF␣ Signaling via NF- to-DC differentiation, especially in the C1, which became ␬B genes upregulated in mo-DCs relative to mo-MACs) (5), specifically demethylated in DCs (Figure 2D). genes downregulated in M-CSF macrophages in compar- Employing the average signal of public MO DNAse- ison with GM-CSF macrophages, (45) and genes upregu- seq (Blueprint database) (41) and DC ATAC-seq triplicates lated in response to GM-CSF and IL-4 cytokines (46) (Sup- (49), we found that C1 and C2 CpGs had low accessibility plementary Figure S1F). We also confirmed the upregula- in MOs. We also observed that C1 had greater accessibil- tion of canonical GR targets, such as PDK4, TSC22D3, ity than C2 in DCs, demonstrating an inverse correlation 114 Nucleic Acids Research, 2022, Vol. 50, No. 1 A B Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 C D E F G H Figure 1. Phenotypic profiling of dexamethasone-mediated tolerogenic dendritic cells. (A) Schematic representation of the experimental approach, com- paring dendritic cell (DC) with tolerogenic dendritic cell (tolDC) differentiation. (B) DC and tolDC were cocultured with CD8+ cells for 5 days. The final Nucleic Acids Research, 2022, Vol. 50, No. 1 115 between chromatin accessibility and DNA demethylation Dynamics of GR and MAFB genomic binding in tolDC dif- (Figure 2E). ferentiation To detect potential TFs involved in the DNA methyla- Since MAFB was predicted to be a TF of importance in tion dynamics, we performed a TF motif enrichment analy- the dexamethasone-specific gene expression (Figure 1H) sis. As expected, C2 CpGs were enriched in Glucocorticoid and DNA methylation (Figure 2F) remodeling occurring in Response Elements (GREs) and very closely related mo- tolDCs, we investigated its role in the tolDC differentiation tifs (Androgen Receptor Elements, AREs and Progesterone process and its interplay with GR. Receptor Elements, PGR). The Maf recognition element MAFB expression was studied over time, where differ- (MARE) was also enriched, including several TFs from the ences between DCs and tolDCs could be observed from 1 MAF family: MAFA, MAFB, MAFF and MAFK (Figure h of differentiation (Supplementary Figure S3A). MAFB 2F). On the other hand, C1 CpGs were enriched in, among protein levels could be detected by western blot in tolDC nu- other binding motifs, AP-1 family, PU.1, IRF8, STAT6 and clei, at 12 and 24 h, concomitant with GR translocation to Egr1/Egr2 motifs (Supplementary Figure S2D), which are the nucleus (Figure 3A). Moreover, by immunofluorescence known to be associated with DC differentiation (6,7,50,51). Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 at 24 h of differentiation, MAFB was found to be localized Gene ontology (GO) analysis of C2 CpGs revealed en- in tolDC nuclei, whereas in DCs the level of expression was richment of functional categories associated with immune much lower (Figure 3B). system regulation, including ‘negative regulation of cy- We next generated ChIP-seq data of GR and MAFB in tokine production involved in inflammatory response’ and DCs and tolDCs, at both 4 and 24 h of differentiation, us- ‘regulation of inflammatory response’ (Figure 2G). In con- ing biological duplicates that we consolidated to obtain the trast, GO analysis of C1 CpGs included categories re- consensus peaks of each condition (see Materials and Meth- lated to inflammatory processes, such as ‘positive reg- ods). De novo motif discovery around GR and MAFB peaks ulation of MAPK cascade’, ‘myeloid leukocyte activa- in tolDC yielded very similar motifs to those of their respec- tion’ and ‘inflammatory response’ (Supplementary Fig- tive canonical binding sites. (Supplementary Figure S3B). ure S2E). Moreover, C2 CpGs are demethylated in M2 GR peak calling in tolDC resulted in hundreds of peaks macrophages (M-CSF), whereas C1 CpGs are demethy- with a strong signal (Figure 3C, Supplementary Figure lated in M1 macrophages (GM-CSF) (Supplementary Fig- S3C). By contrast, GR peak calling in DCs yielded a negli- ure S2F). gible number of significant peaks, at both 4 and 24 h, with We next associated each CpG with its nearest gene and weaker signals (Figure 3C, Supplementary Figure S3C). tested whether C1 and C2 associated genes were enriched On the other hand, MAFB presented binding in thou- in tolDC-induced or tolDC-repressed genes. We found a sands of sites in both DCs and tolDCs at 4h, with stronger strong enrichment of the induced genes over the C2 asso- signal and peak number in tolDCs (Figure 3C, Supplemen- ciated genes (FDR = 8.88e−41) and of the repressed genes tary Figure S3C). Moreover, MAFB DC peaks at 4h were over the C1 associated genes (FDR = 4.72e−26) (Figure highly enriched in the MAFB canonical binding site (Sup- 2H). In this regard, there was a significant inverse correla- plementary Figure S3D) This is consistent with the pre- tion between DNA methylation and gene expression (Fig- existing expression of MAFB in MOs, which is downregu- ure 2I), which is exemplified in some CpGs associated with lated during the MO to DC differentiation (5). In this re- genes linked to the biology of tolDCs and DCs (Figure 2J, gard, the MAFB peak number and signal in DCs at 24 h is Supplementary Figure S2G). minimal in comparison with tolDCs, concordantly to a very In relation to the potential association between MAFB low MAFB expression, and no motif compatible with the and our DNA methylation and expression data, we assessed canonical one was found (Figure 3C, Supplementary Fig- the potential upregulation of MAFB-positively regulated ure S3C, D). genes of M2 macrophages and C2-associated genes in other To distinguish potential specific features of early and late tolDC types. In addition to tolDCs, we found that both peaks of GR and MAFB, we classified the tolDC peaks as gene sets were upregulated in DC10 (DCs differentiated in ‘4 h-specific’ (present at 4h but not at 24 h), ‘24 h-specific’ vitro from MOs with GM-CSF/IL-4 and IL10), but not in (present at 24 h but not at 4 h), or ‘continuous’ (present rapamycin-treated or vitamin D3-treated tolDCs (Supple- at both times) (Figure 3D and Supplementary Table S3). mentary Figure S2H). In both TFs, the subset with the strongest binding was the ←−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− CFSE signal of CD8 + cells is shown (left panel). CD8+ with only CD3/CD28 T-activator beads (C+) or alone (C–) are also shown. In the right panel, the average proliferation is shown (mean ± standard error of the mean (SEM)) (n = 4). (C) IL-10, TNF␣, IL-12p70 and IL1-␤ production of DC and tolDC, after 5 days of differentiation and 24 h of LPS (10 ng/␮l) and IFNg (20 ng/␮l) stimuli. P-values of paired t-tests are shown (n = 4) (ns P > 0.05, * P < 0.05). (D) Box-plots of CD80, CD83, CD86 and HLA-DR surface expression (Median Fluorescence Intensity) in DCs and tolDCs in steady-state or stimulated with LPS (10 ng/␮l) and IFNg (20 ng/␮l) (n = 6) (ns P > 0.05, ** P < 0.01, *** P ≤ 0.001). (E) Gene expression heatmap of differentially expressed genes comparing tolDCs with DCs and also displaying the gene expression values of the precursor cell type (MO) (logFC > 0.5, FDR < 0.05). Scaled fluorescence values of expression arrays are shown, ranging from −2 (lower gene expression, green) to +2 (higher gene expression, orange) (n = 3). (F) Gene ontology (GO) over-representation of GO Biological Process categories. Fold change of tolDC induced genes over background and -log10(FDR) of Fisher’s exact tests are shown. (G) Discriminant regulon expression analysis (DoRothEA) of tolDC compared with DC. Only transcription factors with FDR <0.05 are shown. NES and logFC of transcription factor expression are depicted. (H) T-distributed stochastic neighbor embedding (t-SNE) plot of the aggregated and batch-corrected gene expression data from our study (MO, DC and tolDC) and two additional public datasets (GSE40484 (moMAC, moDC, cDC2, CM (classical MOs) and NCM (non-classical MOs) and GSE99056 (M-MAC (M2 macrophages) and GM-MAC (M1 macrophages)). The four different groups obtained using k-means clustering are represented with grey ellipses of multivariate t-distributions. 116 Nucleic Acids Research, 2022, Vol. 50, No. 1 A B Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 C D E F G H I J Figure 2. DNA methylation remodeling of dexamethasone-mediated tolerogenic dendritic cells. (A) DNA methylation heatmap of differentially methylated CpGs comparing DCs with tolDCs (␤ ≥ 0.2, FDR < 0.05) and showing the DNA methylation values of the precursor cell type (MO). Scaled ␤-values Nucleic Acids Research, 2022, Vol. 50, No. 1 117 ‘continuous’ peaks. The binding typically occurred within 24 h, GR was also bound to the MAFB gene promoter (Fig- open chromatin in MOs, even though some GR and MAFB ure 3G). This rapid binding, together with the rapid MAFB peaks occurred in low-accessibility regions (Supplementary RNA and protein upregulation (Figure 3A, Supplementary Figure S3E). Figure S3A), suggests a direct, GR-mediated, mechanism Whereas GR ‘4h-specific’ and ‘continuous’ peaks were lo- of MAFB regulation. cated mainly in intergenic and intronic regions and were enriched in MO enhancers and active transcription start GR and MAFB binding are correlated with gene upregulation sites (TSS), ‘24 h-specific’ peaks were most frequently found and DNA demethylation in tolDCs in promoters, and were more enriched in active TSS but not enriched in enhancers (Figure 3E, Supplementary Fig- We measured the association between each ChIP-seq peak ure S3D). On the other hand, no notable differences were and its nearest gene in order to examine the link between described in the location of MAFB ‘4h-specific’, ‘24 h- TF binding and gene expression remodeling. Overall, both specific’ and ‘continuous’ peaks, all of them being found MAFB and GR ‘4h-specific’, ‘24 h-specific’, and ‘contin- mostly in intergenic and intronic regions, where they were uous’ associated genes were more strongly expressed in Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 enriched in monocyte active TSS and enhancers (Figure 3E, tolDC than in DCs, indicating that both TFs may be in- Supplementary Figure S3D). volved in upregulating tolDC -specific genes (Figure 4A, B). Overall, there were many more MAFB than GR peaks, It is notable that, in both TFs, the most frequently associ- at 4 and 24 h. Remarkably, GR binding, which was ubiqui- ated subset was the ‘continuous’ associated genes. tously expressed in most cell types, significantly overlapped To determine the relationship between TF binding and with the peaks observed in a panel of 53 GR ChIP-seqs DNA methylation, we first profiled the initial methyla- (Supplementary Figure S3G), suggesting a common role of tion state of CpGs surrounding GR and MAFB peaks, GR across these cell types. using public whole-genome bisulfite sequencing data of We also inspected overlapping GR and MAFB peaks, CD14+ MOs (41). Both TFs bound methylated regions, which involved around half of all the GR peaks (Figure 3F). even though they were located mainly in non-methylated These GR∩MAFB peaks can be classified into four groups, regions (Supplementary Figure S4A). However, GR 24 h- depending on their temporal association between TFs: ‘GR specific peaks were found almost exclusively in non- first’, where GR peaks at 4 h, and MAFB peaks at 24 h, but methylated regions, concordantly with its enrichment in not at 4 h; ‘MAFB first’, where MAFB peaks at 4h, and MO active TSS (Supplementary Figure S3F). GR peaks at 24 h, but not at 4 h; ‘constant’, where GR and We then performed co-immunoprecipitation of TET2, a MAFB peaks at both 4 and 24 h; and ‘transient’, overlap- key mediator of active demethylation in myeloid cells, re- ping peaks not included in any of the former groups (Figure vealing its interaction with both GR and MAFB (Figure 3G). In general, ‘GR first’, ‘constant’ and ‘transient’ peaks 4C) in tolDCs. In this regard, GR ChIP-seq signal in tolDCs were enriched in GREs. In contrast, ‘MAFB first’ peaks, at both 4 and 24 h was found around some CpGs specifically in which MAFB is bound before GR, are enriched in ETS demethylated in tolDCs (C2 CpGs), with a strong signal family motifs, including PU.1, and AP-1 motifs, but not around a small subset of total C2 CpGs (Figure 4D). More- MAREs. (Supplementary Figure S3H). All groups were en- over, MAFB signal in tolDCs at both 4 and 24 h was found riched in MO enhancers and TSSs, except for ‘MAFB first’ more generally around C2 CpGs. Interestingly, MAFB sig- peaks, which were not enriched in enhancers (Supplemen- nal was also present in DC at 4 h, at a lower level than that tary Figure S3I). of tolDCs, but was not present at 24 h (Figure 4E). To address the mechanism of the GC-induced upregula- C1 CpGs were associated with EGR2, among other TFs tion of MAFB, we checked the binding of GR around the (Supplementary Figure S2D), which prompted us to use MAFB promoter. At both 4h and 24h, GR was bound in a public dataset of EGR2-FLAG ChIP-seq of DCs (7), an enhancer in MOs whose closest gene was MAFB. At which revealed a strong signal in C1 CpGs but not in C2 ←−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− are shown (lower and higher DNA methylation levels in blue and red, respectively) (n = 4). On the right side, violin plots of Cluster 1 (C1) and Cluster 2 (C2) depict scaled DNA methylation data. (B) Bar-plot of genomic features, percentages of C1 and C2 CpGs in comparison with background CpGs (Bg). (C) Bar-plot of CpG island contexts, percentages of C1 and C2 CpGs in comparison with background CpGs. (D) Accessibility (ATAC-seq) data of C1 CpGs (red) and C2 CpGs (blue) in MOs and DCs. The average DNAse-seq from an MO duplicate (BLUEprint) and ATAC-seq triplicate from DCs (GSE100374) were used in the representation. (E) ChIP-seq data of H3K27ac and H3K4me1 of CD14+ monocytes and DCs were downloaded from the BLUEPRINT database. Odds ratios of histone marks enrichment were calculated for bins of 10bp, 2000bp upstream and downstream in relation to C1 and C2 CpGs. CpGs annotated in the EPIC array were used as background. (F) Bubble scatter-plot of TF-binding motif enrichment for C2 CpGs. The x-axis shows the percentage of windows containing the motif and the y-axis shows the factor of enrichment of the motif over the EPIC background. Bubbles are colored according to the TF family. FDR is indicated by bubble size. (G) GO (Gene Ontology) over-represented categories in C2 CpGs. Fold change in comparison with background (EPIC array CpGs) and −log10(FDR) are shown. (H) C1 and C2 CpGs were associated with the nearest gene and the enrichment of both gene sets (C1 CpGs- and C2 CpGs- associated genes) over the tolDC-induced and tolDC-repressed genes were calculated using Fisher’s exact tests. Odds ratios ±95% confidence intervals are shown. (I) DNA methylation of differentially methylated CpGs were correlated with gene expression of differentially expressed genes in the tolDC vs. DC comparison. LogFC of expression is represented in the y-axis, where a higher number represents a higher level of expression in tolDC and a lower number a higher level of expression in DC. DNA methylation is depicted on the x-axis as ␤, where a lower number represents a lower level of methylation in tolDC, and a higher number a lower level of methylation in DC. Points are colored according to their genomic context. A significant negative correlation between methylation and expression is observed (R = −0.55, P < 2.2e−16). (J) The fold-changes of DC and tolDC expression (with respect to MO) of some examples of C2 CpGs associated genes are shown, calculated from fluorescence values of expression arrays. Below each gene, the methylation (␤-values) of an example of associated CpG is indicated. Three and four biological replicates are provided, respectively (n = 3–4). 118 Nucleic Acids Research, 2022, Vol. 50, No. 1 Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 G Figure 3. Delineation of GR and MAFB binding during DC and tolDC differentiation. (A) Western blot of glucocorticoid receptor (GR) and MAFB proteins in cytoplasm and nuclei of DC and tolDC. GAPDH and LaminB proteins were used as loading controls for cytoplasm and nuclei, respectively. (B) Immunofluorescence of MAFB in DC and tolDC after 24 h of differentiation. Fluorescent signal of MAFB (green), actin filaments (Phalloidin, red), nuclei (DAPI, blue) and a composite image with the sum of the three fluorescences are shown. (C) Average signal in bins per million mapped reads (BPMs) of GR and MAFB binding in DCs and tolDCs at 4 and 24 h of differentiation, after peak consolidation between two biological replicates of each condition. (D) Venn diagrams showing GR and MAFB peaks at 4, 24 h and at both times (continuous). (E) Bar-plot of genomic features percentages of ‘4 h-specific’, ‘24 h-specific’, and ‘continuous’ peaks of GR and MAFB. (F) Venn diagrams showing GR and MAFB peaks at 4 and 24 h, and the concatenation of the two (total peaks). (G) Heatmap of total peaks with GR/MAFB overlap (GR MAFB). Peaks were classified in ‘GR first’ (GR bound at 4 h and MAFB at 24 h), ‘MAFB first’ (MAFB bound at 4 h and GR at 24 h), ‘Constant’ (GR and MAFB peaks at 4 and 24 h), and ‘Transient’ (overlapping peaks not included in any of the former groups). (H) Representation of the GR ChIP-seq signal (BPMs) close to the MAFB gene. 15-states ChromHMM of MOs in peaks are depicted in orange. Significant peaks are shown against a grey background. Nucleic Acids Research, 2022, Vol. 50, No. 1 119 A B C Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 D E F G H Figure 4. Integration of GR and MAFB binding with DNA methylation and gene expression changes in tolDCs. (A) Gene set enrichment analysis (GSEA) of tolDC vs. DC, with the genes closest to GR ‘continuous’, ‘4 h-specific’ and ‘24 h-specific’ peaks. The running enrichment score is presented with the normalized enrichment score (NES) is shown above (FDR < 0.01). (B) Gene set enrichment analysis (GSEA) of tolDC versus DC, with the genes closest to the top500 MAFB ‘continuous’, ‘4 h-specific’ and ‘24 h-specific’ peaks (500 peaks with the highest joint P-value). The running enrichment score is represented in the y-axis and the normalized enrichment score (NES) is shown above (FDR < 0.01). (C) Western blot of the co-immunoprecipitation of TET2, showing the signal of MAFB, GR and TET2 proteins. On the right, the image signal quantifications of three independent western blots are shown for each protein. P-values of paired t-tests are shown (n = 3) (*P < 0.05, **P < 0.01). (D, E) Heatmaps of the GR and MAFB ChIP-seq signal (BPMs) around C1- and C2-CpGs in DCs and tolDCs at 4 h and 24 h of differentiation. (F, G) Representation of the GR and MAFB ChIP-seq signal (BPMs) close to the SIAH2 and C1QB genes. 15-state ChromHMM of MOs in peaks are depicted in orange. Significant peaks are shown against a grey background. (H) SIAH2 and C1QB time-course gene expression obtained from qRT-PCR (relative arbitrary units) (n = 3) and DNA methylation (␤-values) of associated differentially methylated CpGs obtained from pyrosequencing (n = 3). P-values of paired t-tests are shown (n = 3) (* P < 0.05 ** P ≤ 0.01, *** P ≤ 0.001). 120 Nucleic Acids Research, 2022, Vol. 50, No. 1 CpGs (Supplementary Figure S4B). Strikingly, the EGR2 ing the role of MAFB in the tolDC demethylation pro- gene was downregulated in tolDCs (Supplementary Figure cess (Figure 5C). In contrast, no differences were observed S2G), which may partially explain the blockage of C1 CpG in C1 CpGs, corresponding to the absence of MAREs in demethylation. In addition, the genes closest to the EGR2 the cluster (Supplementary Figure S2B) and the weaker peaks were more strongly expressed in DCs, and downreg- signal of MAFB in the ChIP-seq (Figure 4D). Based ulated in tolDCs, suggesting an association between a po- on their overlap with MAFB and GR ChIP-seq peaks, tential loss-of-function of EGR2 in tolDCs and the specific C2 CpGs were then divided into three groups: MAFB- downregulation of genes (Supplementary Figure S4C). specific, GR/MAFB, and GR-specific CpGs. In siMAFB- We also studied the temporal relation between DNA treated tolDCs, MAFB-specific and GR/MAFB CpGs methylation and gene expression changes. We selected were more methylated, whereas GR-specific CpGs were CpGs with GR and/or MAFB binding associated with not affected by MAFB inhibition (Figure 5D). This sug- tolDC-induced or tolDC-repressed genes (Figure 4F, G, gests that both MAFB and GR can direct demethyla- Supplementary Figure S4D−F). In some loci, DNA tion to C2 CpGs, consistent with their interaction with demethylation was concomitant with or preceded gene up- TET2. Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 regulation (Figure 4H). However, there were also examples Surface markers CD14, CD16 and CD163 were signifi- where gene upregulation clearly anticipated DNA demethy- cantly reduced in tolDCs with the MAFB inhibition (Fig- lation (Supplementary Figure S4G). ure 5E), providing evidence of the functional role of MAFB in the tolerogenic phenotype. Moreover, the inhibition of MAFB also reduced IL-10 MAFB downregulation reverts dexamethasone-induced ex- production at steady-state and after stimulation of the cells pression and DNA methylation remodeling, damping the with LPS. TNF␣ production in stimulated cells was sig- tolerogenic phenotype of tolDCs nificantly higher in siMAFB-treated tolDCs, showing that Given the association between MAFB binding, gene upreg- MAFB inhibition not only reduced the tolerogenic features ulation and DNA demethylation in tolDCs, we performed of tolDCs but also boosted an increase in some proinflam- MAFB knockdown, using small-interfering RNAs (siR- matory traits, concordantly with RNA-seq data (Figure NAs) targeted against MAFB (siMAFB) or non-targeting 5F). (siCTL). Under our conditions, we achieved around 75% Consequently, the suppression of CD8+ T cell prolifera- transfection efficiency (Supplementary Figure S5A), >50% tion, a main feature of tolDC, was reduced with the MAFB RNA reduction (Supplementary Figure S5B) and a drastic inhibition (Figure 5G). Overall, these data indicate that decrease in MAFB protein (Supplementary Figure S5C). MAFB is a key player in the acquisition of tolDC tolero- Gene expression profiling of tolDCs transfected with genesis, and in the transcriptomic and epigenomic events siMAFB or siCTL was performed by RNA-seq, obtain- driving that phenotype. ing 222 downregulated genes and 259 upregulated genes (FDR < 0.05, logFC > 0.5) (Figure 5A and Supple- Glucocorticoids skew an MAFB-associated, monocyte- mentary Table S4). Among the downregulated genes, rele- derived cell differentiation program in rheumatoid arthritis vant tolDC-induced genes were found, including RNASE1, joints CCL18, LGMN, MERTK, IL10, IFIT3 and SLCO2B1. Moreover, the upregulated genes included tolDC-repressed In vitro tolDCs and DCs derived from MOs resembled in genes such as GPT, CD1C, TNF, EGR2, CSF2RB, CD1B vivo cell types, moMACs and moDCs, respectively (Figure and FLT3. 1H). We therefore performed single-cell RNA-seq from un- GO Biological Process-over-represented categories in sorted cells of the synovium of rheumatoid arthritis (RA) siMAFB-downregulated genes included terms such as joints from one patient treated with GCs and one treatment- ‘monocyte chemotaxis’, ‘regulation of tolerance induc- naive patient with a similar histology type and synovi- tion’ and ‘negative regulation of interferon-gamma produc- tis score to further explore the in vivo effects of GCs in tion’ (Supplementary Figure S5D). Furthermore, among monocyte-derived cell populations. siMAFB-upregulated genes, terms such as antigen process- From the total single-cell transcriptomes, we exclusively ing and presentation via MHC class Ib, positive regulation selected putative MOs and MO-derived cells, based on of leukocyte-mediated immunity and leukocyte differentia- the expression of CD14, S100A9, and MRC1 (2), and the tion were over-represented (Supplementary Figure S5E). absence of THY1, CD248, CD27, IGHN, CD3G, CD3E, In this regard, tolDC-induced genes were, in gen- CD34, KLRD1 and NKG7 expression (Supplementary Fig- eral, significantly downregulated with the MAFB siRNA, ure S6A−C). Clustering of these cells yielded five distinct whereas tolDC-repressed genes were upregulated (Figure subpopulations. We then excluded clusters containing fewer 5B). Moreover, genes associated with MAFB ‘continuous’ than 50 cells and annotated the three remaining clusters as and ‘4h-specific’ ChIP-seq peaks were linked to siMAFB M1, M2 and M3 (Figure 6A). downregulation, whereas those GR peaks were not related Looking at the gene markers with a high level of expres- to downregulation or upregulation (Supplementary Figure sion in each group, we found some monocytic and proin- S5F). flammatory genes, such as VCAN, LYZ, S100A8, S100A9 We then tested the effect of MAFB downregulation and IL1B, in the M1 cluster, macrophagic markers (C1QB, on the differentially methylated CpGs. CpGs specifically FCGR3A, CD163, RNASE1 and FOLR2) in the M2 clus- demethylated in tolDCs (C2 CpGs) were more methy- ter, and dendritic cell-related markers (FCER1A, IRF4, lated in tolDC when MAFB was downregulated, confirm- CD1C, CCR7 and CCL17) in the M3 cluster (Figure 6B). Nucleic Acids Research, 2022, Vol. 50, No. 1 121 A B Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 C D E F G Figure 5. Effects of MAFB knockdown during tolDC differentiation. (A) Volcano plot comparing tolDCs treated with control siRNA (siCTL) and MAFB siRNA (siMAFB). Dashed lines indicate significance thresholds (FDR < 0.05, absolute logFC > 0.5) (n = 4). tolDC-induced and tolDC-repressed genes 122 Nucleic Acids Research, 2022, Vol. 50, No. 1 In addition, module scores of blood CD14+ MO, mo- sociated with a small fraction of the epigenomic and tran- MAC and mo-DC signatures displayed an increased signal scriptomic remodeling produced. Instead, after the glu- in M1, M2 and M3 clusters, respectively (Figure 6C). cocorticoid stimulus, GR binds to an enhancer close to We then analyzed the relative proportions of clusters in MAFB, and the upregulation of the gene is observed af- the treatment-naive and the GC-treated patient joints. In ter as little as 1h of differentiation, pointing to a direct, the GC-treated patient joints, M1 and M3 clusters were de- GR-mediated regulation mechanism. MAFB, in turn, binds pleted, whereas the M2 cluster was increased (Figure 6D). thousands of genomic loci in tolDCs, and is involved not In this regard, when M1, M2 and M3 clusters expression only in DNA demethylation and gene upregulation, but were analyzed in bulk, GC-upregulated genes were more also in the acquisition of the tolerogenic phenotype, as expressed in in vitro tolDCs and downregulated with the demonstrated by its knockdown. In this regard, MAFB acts MAFB siRNA (Figure 6E). as a surrogate to induce tolerogenesis in DCs on behalf Moreover, in vitro tolDC-induced genes presented a of GR. higher level of expression in the M2 cluster, whereas tolDC- MAFB is a known downstream target of the IL- repressed genes were more strongly expressed in M3 (Fig- 10/STAT3 signaling pathway (54), and additional evidence Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 ure 6F). Furthermore, the module scores of the genes as- of the major role of MAFB in our model comes from the sociated with MAFB ‘continuous’ ChIP-seq peaks and the finding that in vitro DC-10 (tolDCs generated with IL-10) siMAFB-downregulated genes for each cluster showed an present upregulation of MAFB target genes and C2 CpG- increased signal in the M2 cluster for both modules (Figure associated genes, despite the absence of active GR. 6G), and MAFB is more expressed in the M2 than the M1 A central question about GR biology has been how a TF and M3 clusters (Supplementary Figure S6D), supporting that is ubiquitously expressed across most tissues can have the involvement of MAFB in the transcriptomic signature different cell-type-specific functions(55–57). Chromatin ac- of the M2 cluster. cessibility has been described as pre-determining GR bind- ing and shaping its differential binding across various cell types (58). Our data also indicate a preference for GR to DISCUSSION bind preexisting open chromatin. However, we found a sig- In this study, we demonstrate a fundamental role for nificant overlap of tolDC GR peaks across very different MAFB, in combination with GR, in the acquisition of cell lines. Our findings support that, upon GC exposure, glucocorticoid-mediated tolerogenesis by tolDCs. We show MAFB confers the cell type specificity required to acquire a coordinated action of MAFB and GR in their binding the tolerogenic phenotype of DCs. to genomic sites, DNA methylation and gene expression DNA methylation (5mC) is generally considered a repres- changes, which lead to the establishment of a tolerogenic sive epigenetic mark that is associated with gene downregu- phenotype, in which MAFB plays a predominant role. This lation(59). TET2, the most strongly expressed TET in MOs, major role for MAFB is confirmed in vivo, by examining has been linked to active demethylation events in these cells single-cell data from the synovium of rheumatoid arthritis in several terminal myeloid differentiation models (6–8,60). patients treated with GC. We prove that GR and MAFB both interact with TET2 in GR, as a ubiquitously expressed nuclear factor, is dis- tolDCs, indicating that the demethylation process is proba- tributed fundamentally in the cytoplasm of cells, until GCs bly triggered by MAFB- or GR-driven TET2 recruitment. induce a conformational change and promote its translo- However, we cannot rule out the possibility that other TET cation to the nucleus and binding to the genome in a mat- proteins are involved. ter of minutes (52). This, together with the previously de- DNA methylation changes occurring during the MO-to- scribed pioneer capacity of GR(53), suggested a direct and tolDC differentiation is potentially both a cause and a con- major role in the acquisition of the tolDC phenotype and sequence of the reshaping of gene expression. Some exam- methylome. On the other hand, MAFB has been previously ples show that DNA demethylation can precede gene up- linked to MO differentiation (10) and the M2 macrophage regulation, or can occur at the same time (SIAH2, C1QB), phenotype (9). Strikingly, we show that GR has an impor- whereas other genes, such as CD163 and ETS2, present a tant but limited direct role in this context, given that its noticeable gene upregulation before DNA demethylation, binding occurs in a few hundred genomic loci and is as- as shown in other immune contexts(61). The functional ←−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− are shown in blue and orange, respectively. (B) Gene set enrichment analysis (GSEA) of tolDCs (siCTL) versus tolDCs (siMAFB), using tolDC-induced and tolDC-repressed gene sets. The running enrichment score is represented and the normalized enrichment score (NES) is shown above (FDR < 0.01). (C) DNA methylation heatmap of previously obtained differentially methylated CpGs (C1-CpGs and C2-CpGs) in tolDCs (siCTL) and tolDCs (siMAFB). Scaled ␤-values are shown (lower DNA methylation levels in blue and higher methylation levels in red). On the right side, violin plots of Cluster 1 (C1) and Cluster 2 (C2) depict ␤-values (n = 4) (ns P > 0.05, *** P ≤ 0.001). (D) Methylated CpG set enrichment analysis (mCSEA) of tolDCs (siCTL) versus tolDCs (siMAFB), using MAFB-only CpGs, GR/MAFB CpGs and GR-only CpGs as CpG-sets (depending on the overlap of CpGs with GR or MAFB peaks). The running enrichment score is represented and the normalized enrichment score (NES) and FDR are shown above. (E) Box-plots of median fluorescence intensity (MFI) of CD14, CD16, CD163 and CD1a flow cytometry data from DCs (siCTL), tolDCs (siCTL) and tolDCs (siMAFB) (n = 7) (ns P > 0.05, * P < 0.05, ** P ≤ 0.01). (F) Box-plots of supernatant concentration from DCs (siCTL), tolDCs (siCTL) and tolDCs (siMAFB) (n = 7) of IL-10 in steady-state and stimulated conditions (LPS 10 ng/␮l and IFN␥ 20 ng/␮l) and IL-12p70 and TNF␣ under stimulated conditions (pg/mL). TNF␣ and IL-12p70 in steady state were not detected. (ns P > 0.05, * P < 0.05, ** P ≤ 0.01) (G) DC (siCTL), tolDC (siCTL) and tolDC (siMAFB) were cocultured with CD8 + cells for 5 days (n = 4). The final CFSE signal of CD8+ cells is shown (left panel). CD8+ with only CD3/CD28 T-activator beads (C+) or alone (C-) are also shown. On the right panel, the average proliferation of the quadruplicate is shown (mean ± standard error of the mean (SEM)) (** P ≤ 0.01, *** P ≤ 0.001). Nucleic Acids Research, 2022, Vol. 50, No. 1 123 Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 Figure 6. In vivo effects of glucocorticoid treatment in myeloid cells from RA synovium. (A) UMAP of putative monocyte-derived clusters from synovial tissues identified in the scRNA-seq analysis. A total of one treatment-naive and one GC-treated RA patient is shown in the UMAP. (B) Heatmap of the top 20 genes most differentially expressed from the M1, M2 and M3 clusters. Relevant cluster markers and the total number of genes identified in each cluster are shown. (C) UMAP heatmap of module scores of gene sets associated with blood CD14+ cells, mo-MACs and mo-DCs (67). (D) UMAP is divided depending on the treatment of RA patients (treatment-naive or GC-treated). In the right panel, proportions of each cluster in each group are shown. (E) Gene set enrichment analysis (GSEA) of tolDCs vs. DCs and tolDCs (siCTL) vs. tolDCs (siMAFB), using the genes upregulated with glucocorticoids in M1, M2 and M3 cluster cells from patients as the gene set. The running enrichment score is represented and the normalized enrichment score (NES) is shown above (FDR < 0.01). (F) UMAP heatmap of module scores of tolDC-induced and tolDC-repressed genes. (G) Violin plots of module scores of genes associated with MAFB ‘continuous’ peaks and siMAFB-downregulated genes in the M1, M2 and M3 clusters (**** P ≤ 0.0001). role of the latter is enigmatic. Increasing evidence suggests cytes, and that an increase of phagocytosis is a typical fea- that active demethylation intermediates, such as 5hmC, may ture of tolDCs (63,64). Here, we show that glucocorticoids be epigenetic marks with regulatory functions (62). DNA skew MO-to-DC differentiation through MAFB, resem- demethylation could also have a role later on, stabilizing the bling an in vivo cell type. In this regard, we have shown de- phenotype or fine-tuning the immune response after suc- pletion of cells with an expression pattern similar to moDCs cessive inflammatory stimuli. This is compatible with IL10 and an increase of cells similar to moMACs in RA patients gene behavior; whereby dexamethasone-mediated demethy- treated with GCs. Since in vivo moDCs are involved in the lation precedes a higher level of production of IL-10 after pathogenesis of several inflammatory diseases (3,65,66), the an inflammatory stimulus. GC-mediated remodeling of monocyte-derived populations The transcriptome of GM-CSF/IL-4 DCs is very sim- in the synovium could be a significant process that modi- ilar to the in vivo mo-DCs described in ascitic and syn- fies their proportions and modulates the inflammation pro- ovial fluids (5). Intriguingly, tolDCs are more similar to in duced by monocyte-derived cells in tissues. However, given vivo moMACs and to in vitro M2 macrophages. MAFB up- the limited number of RA samples included in this study, regulation is a hallmark of these three cell types. Several the in vivo effects of GC treatment in MO differentiation studies have produced evidence to suggest that there is a should be further explored in follow up studies with a more macrophagic phenotype of dexamethasone-treated mono- clinical perspective. 124 Nucleic Acids Research, 2022, Vol. 50, No. 1 Our results shed light on the regulatory mechanisms of 4. Sallusto,F. and Lanzavecchi,A. (1994) Efficient presentation of GC-induced tolDC differentiation, identifying the critical soluble antigen by cultured human dendritic cells is maintained by granulocyte/macrophage colony-stimulating factor plus iuterleukin 4 role of MAFB, which takes over GR to fulfil the main and downregulated by tumor necrosis factor ␣. J. Exp. Med., 179, roles of tolerogenesis induction. We have described a new 1109–1118. mechanism of action of GCs in MOs, consisting of MAFB 5. Goudot,C., Coillard,A., Villani,A.C., Gueguen,P., Cros,A., induction, overriding the MO-to-DC differentiation pro- Sarkizova,S., Tang-Huau,T.L., Bohec,M., Baulande,S., Hacohen,N. gram, rendering macrophage-like tolerogenic cells. More- et al. (2017) Aryl hydrocarbon receptor controls monocyte differentiation into dendritic cells versus macrophages. Immunity, 47, over, in MO-derived cells from synovial tissues, we have 582–596. also shown a concordant GC-mediated depletion of mo- 6. Vento-Tormo,R., Company,C., Rodrı́guez-Ubreva,J., de la Rica,L., DCs and an increase of mo-MACs. By improving the un- Urquiza,J.M., Javierre,B.M., Sabarinathan,R., Luque,A., Esteller,M., derstanding of the molecular mechanisms underlying the Aran,J.M. et al. (2016) IL-4 orchestrates STAT6-mediated DNA demethylation leading to dendritic cell differentiation. Genome Biol., tolDC generation mediated by GCs and the effects of GCs 17, 4. in MOs in vivo, our results can provide insights to enhance 7. Mendes,K., Schmidhofer,S., Minderjahn,J., Glatz,D., Kiesewetter,C., the in vitro generation of tolDCs and to create more specific Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 Raithel,J., Wimmer,J., Gebhard,C. and Rehli,M. (2021) The anti-inflammatory therapies. epigenetic pioneer EGR2 initiates DNA demethylation in differentiating monocytes at both stable and transient binding sites. Nat. Commun., 12, 1556. DATA AVAILABILITY 8. Klug,M., Schmidhofer,S., Gebhard,C., Andreesen,R. and Rehli,M. (2013) 5-Hydroxymethylcytosine is an essential intermediate of active DNA methylation, expression, ChIP-seq data and single- DNA demethylation processes in primary human monocytes. cell RNAseq data for this publication have been deposited Genome Biol., 14, R46. in the NCBI Gene Expression Omnibus and are accessible 9. Cuevas,V.D., Anta,L., Samaniego,R., Orta-Zavalza,E., Vladimir de through GEO SuperSeries accession number GSE180542 la Rosa,J., Baujat,G., Domı́nguez-Soto,Á., Sánchez-Mateos,P., Escribese,M.M., Castrillo,A. et al. (2017) MAFB determines human and GSE181082. Flow Cytometry data has been deposited macrophage anti-inflammatory polarization: relevance for the to the FlowRepository under repository ID FR-FCM- pathogenic mechanisms operating in multicentric carpotarsal Z4NS (https://flowrepository.org/). osteolysis. J. Immunol., 198, 2070–2081. 10. Kelly,L.M., Englmeier,U., Lafon,I., Sieweke,M.H. and Graf,T. (2000) MafB is an inducer of monocytic differentiation. EMBO J., 19, SUPPLEMENTARY DATA 1987–1997. 11. Garcı́a-González,P., Morales,R., Hoyos,L., Maggi,J., Campos,J., Supplementary Data are available at NAR Online. Pesce,B., Gárate,D., Larrondo,M., González,R., Soto,L. et al. (2013) A short protocol using dexamethasone and monophosphoryl lipid A ACKNOWLEDGEMENTS generates tolerogenic dendritic cells that display a potent migratory capacity to lymphoid chemokines. J. Transl. Med., 11, 128. We are very grateful to Prof. Angel Corbı́ and his team 12. Unger,W.W., Laban,S., Kleijwegt,F.S., van der Slik,A.R. and for useful discussions and feedback. We thank the CERCA Roep,B.O. (2009) Induction of Treg by monocyte-derived DC modulated by vitamin D3 or dexamethasone: differential role for Programme/Generalitat de Catalunya and the Josep Car- PD-L1. Eur. J. Immunol., 39, 3147–3159. reras Foundation for institutional support. 13. Bell,G.M., Anderson,A.E., Diboll,J., Reece,R., Eltherington,O., Harry,R.A., Fouweather,T., MacDonald,C., Chadwick,T., McColl,E. et al. (2017) Autologous tolerogenic dendritic cells for rheumatoid FUNDING and inflammatory arthritis. Ann. Rheum. Dis., 76, 227–234. E.B. was funded by the Spanish Ministry of Science and 14. Kurochkina,Y., Tikhonova,M., Tyrinova,T., Leplina,O., Sizikov,A., Sulutian,A., Chumasova,O., Ostanin,A. and Chernykh,E. (2018) Innovation (MICINN) [SAF2017-88086-R, PID2020- SAT0212 The safety and tolerability of intra-articular injection of 117212RB-I00; AEI/10.13039/501100011033]; O.M.-P. tolerogenic dendritic cells in patients with rheumatoid arthritis: the holds an i-PFIS PhD Fellowship [IFI17/00034] from preliminary results. In Annals of the Rheumatic Diseases. BMJ, 77, Acción Estratégica en Salud 2013–2016 ISCIII, co- 966.2–967. financed by Fondo Social Europeo; C.O. and M.F.-B. 15. Zubizarreta,I., Flórez-Grau,G., Vila,G., Cabezón,R., España,C., Andorra,M., Saiz,A., Llufriu,S., Sepulveda,M., Sola-Valls,N. et al. were supported by the Swiss National Science Foun- (2019) Immune tolerance in multiple sclerosis and neuromyelitis dation [320030 176061]; M.H. was supported by the optica with peptide-loaded tolerogenic dendritic cells in a phase 1b Iten-Kohaut Foundation and the MLR foundation. Fund- trial. Proc. Natl. Acad. Sci. U.S.A., 116, 8463–8470. ing for open access charge: Spanish Ministry of Science 16. Jauregui-Amezaga,A., Cabezón,R., Ramı́rez-Morros,A., España,C., Rimola,J., Bru,C., Pinó-Donnay,S., Gallego,M., Masamunt,M.C., and Innovation (MICINN) [PID2020-117212RB-I00]. Ordás,I. et al. (2015) Intraperitoneal administration of autologous Conflict of interest statement. None declared. tolerogenic dendritic cells for refractory crohn’s disease: a phase I study. J. Crohns. Colitis, 9, 1071–1078. 17. Adcock,I.M., Nasuhara,Y., Stevens,D.A. and Barnes,P.J. (1999) REFERENCES Ligand-induced differentiation of glucocorticoid receptor (GR) 1. Morante-Palacios,O., Fondelli,F., Ballestar,E. and trans-repression and transactivation: preferential targetting of Martı́nez-Cáceres,E.M. (2021) Tolerogenic dendritic cells in NF-kappaB and lack of I-kappaB involvement. Br. J. Pharmacol., autoimmunity and inflammatory diseases. Trends Immunol., 42, 127, 1003–1011. 59–75. 18. McDowell,I.C., Barrera,A., D’Ippolito,A.M., Vockley,C.M., 2. Coillard,A. and Segura,E. (2019) In vivo differentiation of human Hong,L.K., Leichter,S.M., Bartelt,L.C., Majoros,W.H., Song,L., monocytes. Front. Immunol., 10, 1907. Safi,A. et al. (2018) Glucocorticoid receptor recruits to enhancers and 3. Segura,E., Touzot,M., Bohineust,A., Cappuccio,A., Chiocchia,G., drives activation by motif-directed binding. Genome Res., 28, Hosmalin,A., Dalod,M., Soumelis,V. and Amigorena,S. (2013) 1272–1284. Human inflammatory dendritic cells induce Th17 cell differentiation. Immunity, 38, 336–348. Nucleic Acids Research, 2022, Vol. 50, No. 1 125 19. Hoffman,J.A., Trotter,K.W., Ward,J.M. and Archer,T.K. (2018) 40. Quinlan,A.R. and Hall,I.M. (2010) BEDTools: a flexible suite of BRG1 governs glucocorticoid receptor interactions with chromatin utilities for comparing genomic features. Bioinformatics, 26, 841–842. and pioneer factors across the genome. Elife, 7, e35073. 41. Fernández,J.M., De La Torre,V., Richardson,D., Royo,R., 20. Carvalho,B.S. and Irizarry,R.A. (2010) A framework for Puiggròs,M., Moncunill,V., Fragkogianni,S., Clarke,L., oligonucleotide microarray preprocessing. Bioinformatics, 26, Consortium,BLUEPRINT, Flicek,P. et al. (2016) The BLUEPRINT 2363–2367. Data Analysis Portal. Cell Syst., 3, 491–495. 21. Morante-Palacios,O. and Ballestar,E. (2021) shinyÉPICo: a graphical 42. Leek,J.T., Johnson,W.E., Parker,H.S., Jaffe,A.E. and Storey,J.D. pipeline to analyze Illumina DNA methylation arrays. Bioinformatics, (2012) The sva package for removing batch effects and other 37, 257–259. unwanted variation in high-throughput experiments. Bioinformatics, 22. Aryee,M.J., Jaffe,A.E., Corrada-Bravo,H., Ladd-Acosta,C., 28, 882–883. Feinberg,A.P., Hansen,K.D. and Irizarry,R.A. (2014) Minfi: a flexible 43. Krenn,V., Morawietz,L., Burmester,G.R., Kinne,R.W., and comprehensive Bioconductor package for the analysis of Mueller-Ladner,U., Muller,B. and Haupl,T. (2006) Synovitis score: Infinium DNA methylation microarrays. Bioinformatics, 30, discrimination between chronic low-grade and high-grade synovitis. 1363–1369. Histopathology, 49, 358–364. 23. Kim,D., Langmead,B. and Salzberg,S.L. (2015) HISAT: A fast 44. Humby,F., Lewis,M., Ramamoorthi,N., Hackney,J.A., Barnes,M.R., spliced aligner with low memory requirements. Nat. Methods, 12, Bombardieri,M., Setiadi,A.F., Kelly,S., Bene,F., DiCicco,M. et al. Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 357–360. (2019) Synovial cellular and molecular signatures stratify clinical 24. Li,H., Handsaker,B., Wysoker,A., Fennell,T., Ruan,J., Homer,N., response to csDMARD therapy and predict radiographic progression Marth,G., Abecasis,G. and Durbin,R. (2009) The Sequence in early rheumatoid arthritis patients. Ann. Rheum. Dis., 78, 761–772. Alignment/Map format and SAMtools. Bioinformatics, 25, 45. Sierra-Filardi,E., Puig-Kröger,A., Blanco,F.J., Nieto,C., Bragado,R., 2078–2079. Palomero,M.I., Bernabéu,C., Vega,M.A. and Corbı́,A.L. (2011) 25. Liao,Y., Smyth,G.K. and Shi,W. (2014) FeatureCounts: an efficient Activin A skews macrophage polarization by promoting a general purpose program for assigning sequence reads to genomic proinflammatory phenotype and inhibiting the acquisition of features. Bioinformatics, 30, 923–930. anti-inflammatory macrophage markers. Blood, 117, 5092–5101. 26. Love,M.I., Huber,W. and Anders,S. (2014) Moderated estimation of 46. Rutella,S., Bonanno,G., Procoli,A., Mariotti,A., de Ritis,D.G., fold change and dispersion for RNA-seq data with DESeq2. Genome Curti,A., Danese,S., Pessina,G., Pandolfi,S., Natoni,F. et al. (2006) Biol., 15, 550. Hepatocyte growth factor favors monocyte differentiation into 27. Li,H. (2013) Aligning sequence reads, clone sequences and assembly regulatory interleukin (IL)-10++IL-12low/neg accessory cells with contigs with BWA-MEM. arXiv doi: dendritic-cell features. Blood, 108, 218–227. https://arxiv.org/abs/1303.3997v2, 26 May 2013, preprint: not peer 47. Garcia-Alonso,L., Holland,C.H., Ibrahim,M.M., Turei,D. and reviewed. Saez-Rodriguez,J. (2019) Benchmark and integration of resources for 28. Amemiya,H.M., Kundaje,A. and Boyle,A.P. (2019) The ENCODE the estimation of human transcription factor activities. Genome Res., Blacklist: Identification of Problematic Regions of the Genome. Sci. 29, 1363–1375. Rep., 9, 9354. 48. Garcı́a-González,P.A., Maggi,J., Schinnerling,K., 29. Kharchenko,P.V., Tolstorukov,M.Y. and Park,P.J. (2008) Design and Sepúlveda-Gutiérrez,A., Soto,L., Neira,O., Mehdi,A.M., Nel,H.J., analysis of ChIP-seq experiments for DNA-binding proteins. Nat. Pesce,B., Aravena,O. et al. (2019) Regulation of tolerogenic features Biotechnol., 26, 1351–1359. on dexamethasone-modulated MPLA-activated dendritic cells by 30. Ramı́rez,F., Ryan,D.P., Grüning,B., Bhardwaj,V., Kilpert,F., MYC. Front. Immunol., 10, 1171. Richter,A.S., Heyne,S., Dündar,F. and Manke,T. (2016) deepTools2: 49. Johnson,J.S., Lucas,S.Y., Amon,L.M., Skelton,S., Nazitto,R., a next generation web server for deep-sequencing data analysis. Carbonetti,S., Sather,D.N., Littman,D.R. and Aderem,A. (2018) Nucleic Acids Res., 44, W160–W165. Reshaping of the dendritic cell chromatin landscape and interferon 31. Zerbino,D.R., Johnson,N., Juettemann,T., Wilder,S.P. and Flicek,P. pathways during HIV infection. Cell Host Microbe, 23, 366–381. (2014) WiggleTools: parallel processing of large collections of 50. Bakri,Y., Sarrazin,S., Mayer,U.P., Tillmanns,S., Nerlov,C., Boned,A. genome-wide datasets for visualization and statistical analysis. and Sieweke,M.H. (2005) Balance of MafB and PU.1 specifies Bioinformatics, 30, 1008–1009. alternative macrophage or dendritic cell fate. Blood, 105, 2707–2716. 32 Zhang,Y., Liu,T., Meyer,C.A., Eeckhoute,J., Johnson,D.S., 51. Sichien,D., Scott,C.L., Martens,L., Vanderkerken,M., Van Gassen,S., Bernstein,B.E., Nusbaum,C., Myers,R.M., Brown,M., Li,W. et al. Plantinga,M., Joeris,T., De Prijck,S., Vanhoutte,L., (2008) Model-based analysis of ChIP-Seq (MACS). Genome Biol., 9, Vanheerswynghels,M. et al. (2016) IRF8 transcription factor controls R137. survival and function of terminally differentiated conventional and 33. Jalili,V., Matteucci,M., Masseroli,M. and Morelli,M.J. (2015) Using plasmacytoid dendritic cells, respectively. Immunity, 45, 626–640. combined evidence from replicates to evaluate ChIP-seq peaks. 52. Htun,H., Barsony,J., Renyi,I., Gould,D.L. and Hager,G.L. (1996) Bioinformatics, 31, 2761–2769. Visualization of glucocorticoid receptor translocation and 34. McLean,C.Y., Bristor,D., Hiller,M., Clarke,S.L., Schaar,B.T., intranuclear organization in living cells with a green fluorescent Lowe,C.B., Wenger,A.M. and Bejerano,G. (2010) GREAT improves protein chimera. Proc. Natl. Acad. Sci. U.S.A., 93, 4845–4850. functional interpretation of cis-regulatory regions. Nat. Biotechnol., 53. Johnson,T.A., Chereji,R.V., Stavreva,D.A., Morris,S.A., Hager,G.L. 28, 495–501. and Clark,D.J. (2018) Conventional and pioneer modes of 35. Yu,G., Wang,L.G., Han,Y. and He,Q.Y. (2012) ClusterProfiler: an R glucocorticoid receptor interaction with enhancer chromatin in vivo. package for comparing biological themes among gene clusters. Omi. Nucleic Acids Res., 46, 203–214. A J. Integr. Biol., 16, 284–287. 54. Gemelli,C., Zanocco Marani,T., Bicciato,S., Mazza,E.M., 36. Zhu,L.J., Gazin,C., Lawson,N.D., Pagès,H., Lin,S.M., Lapointe,D.S. Boraschi,D., Salsi,V., Zappavigna,V., Parenti,S., Selmi,T., and Green,M.R. (2010) ChIPpeakAnno: a Bioconductor package to Tagliafico,E. et al. (2014) MafB is a downstream target of the annotate ChIP-seq and ChIP-chip data. BMC Bioinformatics, 11, 237. IL-10/STAT3 signaling pathway, involved in the regulation of 37. Lopez-Delisle,L., Rabbani,L., Wolff,J., Bhardwaj,V., Backofen,R., macrophage de-activation. Biochim. Biophys. Acta, 1843, 955–964. Grüning,B., Ramı́rez,F. and Manke,T. (2021) pyGenomeTracks: 55. Evans,R.M. (1988) The steroid and thyroid hormone receptor reproducible plots for multivariate genomic datasets. Bioinformatics, superfamily. Science, 240, 889–895. 37, 422–423. 56. Franco,L.M., Gadkari,M., Howe,K.N., Sun,J., Kardava,L., 38. Martorell-Marugán,J., González-Rumayor,V. and Carmona-Sáez,P. Kumar,P., Kumari,S., Hu,Z., Fraser,I.D.C., Moir,S. et al. (2019) (2019) MCSEA: detecting subtle differentially methylated regions. Immune regulation by glucocorticoids can be linked to cell Bioinformatics, 35, 3257–3262. type-dependent transcriptional responses. J. Exp. Med., 216, 384–406. 39. Chèneby,J., Ménétrier,Z., Mestdagh,M., Rosnet,T., Douida,A., 57. Hardy,R.S., Raza,K. and Cooper,M.S. (2020) Therapeutic Rhalloussi,W., Bergon,A., Lopez,F. and Ballester,B. (2020) ReMap glucocorticoids: mechanisms of actions in rheumatic diseases. Nat. 2020: a database of regulatory regions from an integrative analysis of Rev. Rheumatol., 16, 133–144. human and Arabidopsis DNA-binding sequencing experiments. 58. John,S., Sabo,P.J., Thurman,R.E., Sung,M.H., Biddie,S.C., Nucleic. Acids. Res., 48, D180–D188. Johnson,T.A., Hager,G.L. and Stamatoyannopoulos,J.A. (2011) 126 Nucleic Acids Research, 2022, Vol. 50, No. 1 Chromatin accessibility pre-determines glucocorticoid receptor macrophages that share functional and phenotypical aspects with binding patterns. Nat. Genet., 43, 264–268. erythroblastic island macrophages. Haematologica, 103, 395–405. 59. Luo,C., Hajkova,P. and Ecker,J.R. (2018) Dynamic DNA 64. van der Goes,A., Hoekstra,K., van den Berg,T.K. and Dijkstra,C.D. methylation: In the right place at the right time. Science (80-)., 361, (2000) Dexamethasone promotes phagocytosis and bacterial killing 1336–1340. by human monocytes/macrophages in vitro. J. Leukoc. Biol., 67, 60. de la Rica,L., Rodrı́guez-Ubreva,J., Garcı́a,M., Islam,A.B.M.M.K., 801–807. Urquiza,J.M., Hernando,H., Christensen,J., Helin,K., 65. Kamada,N., Hisamatsu,T., Okamoto,S., Chinen,H., Kobayashi,T., Gómez-Vaquero,C. and Ballestar,E. (2013) PU.1 target genes Sato,T., Sakuraba,A., Kitazume,M.T., Sugita,A., Koganei,K. et al. undergo Tet2-coupled demethylation and DNMT3b-mediated (2008) Unique CD14 intestinal macrophages contribute to the methylation in monocyte-to-osteoclast differentiation. Genome Biol., pathogenesis of Crohn disease via IL-23/IFN-gamma axis. J. Clin. 14, R99. Invest., 118, 2269–2280. 61. Pacis,A., Mailhot-Léonard,F., Tailleux,L., Randolph,H.E., 66. Zaba,L.C., Fuentes-Duculan,J., Eungdamrong,N.J., Abello,M.V., Yotova,V., Dumaine,A., Grenier,J.-C. and Barreiro,L.B. (2019) Gene Novitskaya,I., Pierson,K.C., Gonzalez,J., Krueger,J.G. and activation precedes DNA demethylation in response to infection in Lowes,M.A. (2009) Psoriasis is characterized by accumulation of human dendritic cells. Proc. Natl. Acad. Sci., 116, 6938–6943. immunostimulatory and Th1/Th17 cell-polarizing myeloid dendritic 62. He,B., Zhang,C., Zhang,X., Fan,Y., Zeng,H., Liu,J., Meng,H., cells. J. Invest. Dermatol., 129, 79–88. Downloaded from https://academic.oup.com/nar/article/50/1/108/6459126 by guest on 19 January 2023 Bai,D., Peng,J., Zhang,Q. et al. (2021) Tissue-specific 67. Tang-Huau,T.L., Gueguen,P., Goudot,C., Durand,M., Bohec,M., 5-hydroxymethylcytosine landscape of the human genome. Nat. Baulande,S., Pasquier,B., Amigorena,S. and Segura,E. (2018) Human Commun., 12, 4249. in vivo-generated monocyte-derived dendritic cells and macrophages 63. Heideveld,E., Hampton-O’Neil,L.A., Cross,S.J., van Alphen,F.P.J., cross-present antigens through a vacuolar pathway. Nat. Commun., 9, van den Biggelaar,M., Toye,A.M. and van den Akker,E. (2018) 2570. Glucocorticoids induce differentiation of monocytes towards