K B -Ova Littermate

Littermate Suv39h1

Leticia Laura Niborski

Paul Guegen

Mengliang Ye

Allan Thiolat

Rodrigo Nalio Ramos

Pamela Caudana

Jordan Denizeau

Christel Goudot

Jean-Michel Luccarini

Anne Soudé

Bruno Bournique

Pierre Broqua

Luigia Pace
email: luigia.pace@iigm.it

Sylvain Baulande

Christine Sedlik

Jean-Pierre Quivy

Geneviève Almouzni
email: genevieve.almouzni@curie.fr

José L Cohen

Elina Zueva

Joshua J Waterfall
email: joshua.waterfall@curie.fr

Sebastian Amigorena
email: sebastian.amigorena@curie.fr

Eliane Piaggio
email: eliane.piaggio@curie.fr

Paul Gueguen

Christina Metoikidou

Thomas Dupic

Myriam Lawand

Sonia Lameiras

Olivier Lantz

Nicolas Girard

Agathe Seguin-Givelet

Marine Lefevre

Thierry Mora

Aleksandra M Walczak

Epigenetic control of CD8 + T cell responsiveness to anti-PD-1 by Suv39h1

RESEARCH ARTICLE

◥

IMMUNOLOGY

The epigenetic control of stemness in CD8 + T cell fate commitment Luigia Pace, 1,2,3 * † Christel Goudot, 1,2 Elina Zueva, 1,2 Paul Gueguen, 1,2 Nina Burgdorf, 1,2 Joshua J. Waterfall, 1,4,5 Jean-Pierre Quivy, 1,6,7 Geneviève Almouzni, 1,6,7 * Sebastian Amigorena 1,2 * After priming, naïve CD8 + T lymphocytes establish specific heritable transcription programs that define progression to long-lasting memory cells or to short-lived effector cells. Although lineage specification is critical for protection, it remains unclear how chromatin dynamics contributes to the control of gene expression programs. We explored the role of gene silencing by the histone methyltransferase Suv39h1. In murine CD8 + T cells activated after Listeria monocytogenes infection, Suv39h1-dependent trimethylation of histone H3 lysine 9 controls the expression of a set of stem cell-related memory genes. Single-cell RNA sequencing revealed a defect in silencing of stem/memory genes selectively in Suv39h1-defective T cell effectors. As a result, Suv39h1-defective CD8 + T cells show sustained survival and increased long-term memory reprogramming capacity. Thus, Suv39h1 plays a critical role in marking chromatin to silence stem/memory genes during CD8 + T effector terminal differentiation.

M

emory T lymphocytes provide lifelong protection against pathogens and cancer (1). In contrast to naïve and effector T cells, memory cells possess unique properties of "stemness," enabling longterm survival and plasticity to replenish effector pools after renewed antigen challenges (2). Understanding the lineage relationships among naïve, effector, and memory T cells, as well as the molecular pathways that regulate gene expression during the transitions from one to another of these distinct states, is essential for the rational design of vaccines and the development of new immune-therapeutic protocols (3).

Although many studies have characterized the transcription factors that control the differentiation of T cells, the corresponding epigenetic states and associated chromatin dynamics involved in the establishment and maintenance of CD8 + T cell memory and effector identities is still incompletely understood (4-7). Several epigenetic pathways, including trimethylated histone H3 Lys 9 (H3K9me3)/HP-1/Suv39h1 and Polycomb repressive complexes, can contribute to establishing or maintaining transcriptional silencing (7)(8)(9)(10)(11). The H3K9me3 modification is considered to be a repressive mark, a hallmark of both constitutive and facultative heterochromatin (12), most often associated with silent gene loci. Mouse Suv39h1 and Suv39h2, two H3K9 site-specific histone methyltransferases (HMTs), are critical heterochromatin regulators (10,13). Suv39h1 is involved in heterochromatin organization, gene silencing, and lineage stability (10,11,14). It also limits somatic reprogramming of differentiated cells into induced pluripotent stem cells (15 ). In B cells and CD4 + T cells, Suv39h1 is involved in gene silencing and lineage plasticity (16,17). Although the mechanisms underlying the induction of genes critically involved in effector and memory T cell generation have been extensively analyzed, the impact of heterochromatin-dependent gene expression silencing on the fates of T lymphocytes during differentiation has not been addressed. Here, we explore the role of Suv39h1-dependent gene silencing in the establishment and maintenance of memory CD8 + T cell stemness, plasticity, and transition to terminally differentiated effectors.

Long-term protection against Listeria monocytogenes infection requires Suv39h1

To investigate the role of Suv39h1 in T cell responses to infectious agents, we infected Suv39h1knockout (KO) and wild-type littermate mice with OVA-expressing L. monocytogenes (LM-OVA) (Fig. 1A). The bacterial burden in the spleen and liver w a sm e a s u r e do nd a y3( p e a ko fi n f e c t i o ni nt h e spleen) and day 7 (resolution). After primary infection, we observed similar bacteremia in the spleen and liver of Suv39h1-KO and littermate control mice, both at the peak and after the resolution of infection (Fig. 1B). In contrast, upon a secondary challenge with LM-OVA 48 days after the primary infection, most littermate control mice showed complete protection, whereas high levels of LM-OVA were detected in the liver and/ or spleen in more than 85% of Suv39h1-KO mice (Fig. 1C). Because protection against secondary L. monocytogenes infection is primarily mediated by CD8 + Tcells (18), these results suggested a defect in the CD8 + T cell response against L. monocytogenes in Suv39h1-KO mice.

To further investigate the role of Suv39h1 in antigen-specific CD8 + T cell responses to LM-OVA, we measured the antigen-specific CD8 + T cell response using SIINFEKL (Ser-Ile-Ile-Asn-Phe-Glu-Lys-Leu)-H-2K b multimers (K b -OVA;Fig. 1D).We observed a factor of 2 reduction in the percentage and absolute numbers of K b -OVA + CD8 + T cells in the blood (Fig. 1, D andE) and a factor of 9 decrease in the numbers of K b -OVA + CD8 + T cells in the spleen (Fig. 1F) on day 7 in Suv39h1-KO mice relative to littermates. An analysis of K b -OVA + CD8 + T cells on days 5, 6, and 7 after LM-OVA infection suggests that the reduced antigen-specific response was not due to differences in the survival of Suv39h1-KO cells (fig. S1). Upon ex vivo restimulation, the percentages of interferon-g + (IFN-g + ) and granzyme B + CD8 + T cells in Suv39h1-KO mice were lower than those of littermates (Fig. 1,D,G, and H, and fig. S2). To address whether the reduced K b -OVA + CD8 + TcellnumbersinLM-OVA-infected Suv39h1-KO mice was T cell-intrinsic, we adoptively transferred Suv39h1-proficient or -deficient T cell receptor transgenic OT-I CD8 + T cells (which recognize OVA) to wild-type mice. Five days after infection with LM-OVA, the numbers of Suv39h1-KO OT-I cells were reduced relative to wild-type OT-I (Fig. 1I). Although the percentage of IFN-g + OT-I cells was similar, we observed a lower frequency of granzyme B + T cells (fig. S3, A to C), confirming that the defect observed in Suv39h1-KO mice was intrinsic to CD8 + T cells. In line with these in vivo results, lower numbers of Suv39h1-KO effector OT-I CD8 + T cells were recovered in vitro under effector-polarizing conditions, relative to numbers of Suv39h1-proficient T cells. Notably, the defective expansion in Suv39h1-KO T cells was overcome by the overexpression of SUV39H1 by a retroviral expression vector (fig. S4, A andB). Thus, CD8 + T cell responses to LM-OVA are impaired in Suv39h1-KO mice because of a CD8 + T cell-intrinsic defect.

Transcriptional silencing by Suv39h1

To investigate whether the reduced K b -OVA + CD8 + T cell response results from defective gene expression programming, we analyzed RNA profiles by Affymetrix microarrays of purified naïve and K b -OVA + CD8 + T cells, which were wild-type or Suv39h1-KO and isolated by fluorescence-activated cell sorting (FACS) 7 days after LM-OVA infection (Fig. 2A and fig. S5A). The overall number of transcripts up-regulated in K b -OVA + CD8 + T cells from infected mice, as compared to naïve CD8 + T cells, was highly similar between littermates (1571) and Suv39h1-KO (1433) T cells (Fig. 2, B andC). In contrast, the number of transcripts significantly down-regulated in Suv39h1-KO mice (1108) was lower relative to littermates (1738) ( F i g .2 ,Ba n dC ) .T h ea b s e n c eo fS u v 3 9 h 1e x p r e ssion impaired the transcriptional silencing of app r o x i m a t e l y9 9 7g e n e s( F i g. 2C), including genes involved in critical T cell functions such as Il7r (CD127), Sell (CD62L), Ccr7,a ndCxcr4 (Fig. 2D) (19), consistent with the known repressive role of Suv39h1 in gene expression.

To examine the phenotype of the antigenspecific CD8 + T cells that accumulate in Suv39h1-KO mice, we used gene set enrichment analysis (GSEA) by pairwise comparison between differentiated Suv39h1-KO and wild-type K b -OVA + CD8 + T cells isolated 7 days after LM-OVA infection (Fig. 2E). Suv39h1-KO cells showed significant enrichment in both memory and naïve gene signatures, whereas the effector signature was not significantly enriched in either wild-type or Suv39h1-KO K b -OVA + CD8 + T cells (Fig. 2E). The interrogation of a series of public gene signatures revealed that Suv39h1-KO K b -OVA + CD8 + T cells were also significantly enriched for a "lymphoid-stem cell" signature (20)(Fig. 2E, right). For further investigation of the nature of these differences, the core lymphoid-stem cell gene signature (20) H3K9me3 enriched islands are represented. Transcribed control genes are also shown in the right panels. (D) H3K9me3 enrichment is plotted against mRNA expression of stem cell/memory signature genes. H3K9me3 island deposition in in vitro differentiated effectors correlates with stem/memory gene silencing in wild-type dump -K b -OVA + CD8 + T cells and is impaired in Suv39h1-KO effectors. Pearson linear regression is displayed. (E)Littermate and Suv39h1-KO polyclonal naïve (CD44 lo CD62L + CD127 + )CD8 + T cells, central memory (CD44 hi CD62L + CD127 + )CD8 + T cells, and effector (dump - CD44 hi CD127 lo/-KLRG1 + K b -OVA + )CD8 + Tcells isolated 7 days after LM-OVA infection were analyzed by mChIP-qPCR using antibodies to H3K9me3. qPCR was performed with primers specific for the promoters of Il7r and Sell.

(F) Naïve and effector (dump -CD44 hi CD127 lo/-KLRG1 + K b -OVA + )CD8 + Tcells from day 7 LM-OVA infected mice were analyzed by mChIP-qPCR using antibodies to H3K9ac. Data are means ± SEM of at least three independent experiments. *P <0 . 05,**P < 0.01; for (E) and (F), analysis of variance, Tu k e y t e s t .

RESEARCH | RESEARCH ARTICLE

on January 21, 2021 http://science.sciencemag.org/ Downloaded from to be commonly enriched both in memory CD8 + T cells and long-term hematopoietic stem cells (21), overexpressed in mouse and human CD8 + T stem cell-like memory subset (22)andalsoretrieved from signaling pathways that regulate stem cell pluripotency from the KEGG database. We designated this as a "stem cell-like memory" signature (composed of 86 genes; table S1). K b -OVA + CD8 + TcellsfromLM-OVA-infected Suv39h1-KO mice were strongly enriched for this signature ( F i g .2 Fa n df i g.S 5 ,Ba n dC ) .T h ee x p r e s s i o n levels of stem cell-like and memory markers Il7r, Ly6a (Ly6a/e, Sca-1), Fas (CD95), and Pou6f1 were all increased in Suv39h1-KO K b -OVA + CD8 + cells relative to cells from control littermates (table S1 and fig. S5C). The CD8 + T-stem cell-like memory signature was not enriched in naïve Suv39h1-KO CD8 + T cells relative to wild-type littermates (Fig. 2G, BubbleGUM representation, top) (23). Thus, by comparison to wild-type cells, Suv39h1-KO K b -OVA + CD8 + T cells overexpress a series of genes associated with stem cells and memory functions.

To search for these stem cell/memory CD8 + T cells in wild-type animals, we isolated the following wild-type subpopulations by FACS (from dump -CD44 hi K b -OVA + CD8 + Tce ll s,da y7aft er LM-OVA infection): (i) CD127 hi KLRG1 -memory precursors, (ii) CD127 lo/-KLRG1 + effectors, (iii) CD127 + KLRG1 + double-positive cells, and (iv) CD127 -KLRG1 -double-negative cells (Fig. 2G, left), as described previously (5,19). The RNA expression profile of these four FACS-sorted populations was analyzed using Affymetrix microarrays. We tested each of the four subsets (test class) against a combination of the three other FACS-sorted subpopulations, using high-throughput GSEA for the CD8 + T naïve, memory, effector, and stem cell-like/memory signatures (Fig. 2G, lower right). As expected, sorted memory precursors were enriched for both memory and naïve signatures, whereas sorted effectors were enriched for the effector signature (Fig. 2G, right). Notably, only sorted memory precursors from wild-typemicewereenrichedfortheCD8 + Tstem cell-like signature when compared to the three other sorted subpopulations (Fig. 2G, lower right). Thus, the stem cell-like memory signature, which is specific to memory precursors in wild-type T cells, is broadly enriched in bulk K b -OVA + CD8 + T cells in Suv39h1-KO mice.

H3K9me3 deposition at stem/memory-associated loci

To elucidate the epigenetic states (11)associated with the stem/memory and effector (19,21)gene loci (tables S1 and S2), we first analyzed a set of public chromatin immunoprecipitation sequencing (ChIP-seq) data sets, generated with naïve, memory precursor, and effector OT-I CD8 + T cells isolated 8 days after LM-OVA infection in vivo (24). We identified the regions significantly enriched for monomethylated histone H3 Lys 4 (H3K4me1), trimethylated histone H3 Lys 4 (H3K4me3), acetylated histone H3 Lys 27 (H3K27ac), and trimethylated histone H3 Lys 27 (H3K27me3) histone marks. For each histone mark, the significantly enriched islands were associated with the nearest gene (25). Differences in histone mark relative enrichment were calculated by pairwise comparison of terminal effectors to naïve cells (Fig. 3A, top) or to memory precursors (Fig. 3A, bottom). The regions associated with the stem cell-like memory genes showed an enrichment in H3K4me1 and H3K4me3 in both naïve cells and memory precursors, as compared to terminal effectors. The enrichment in these active marks, H3K4me1 and H3K4me3, corresponds to the transcription level of the stem cell-like memory genes in both wild-type naïve and memory precursors, as compared to effectors (fig. S5B). In contrast, the profile of H3K27ac deposition, typically distal, was more difficult to associate with target genes, including certain stem/memory genes and most effector signature genes in effectors (Fig. 3A); of note, cell cycle genes were excluded from this signature (table S2). In line with these results, H3K4me1 and H3K27ac marked the effector genes mostly in the effector subset (Fig. 3A). The enrichment of the repressive mark H3K27me3 on stem cell-like memory genes corresponded with the silencing of the gene signature in the effector CD8 + T cells (Fig. 3A and fig. S5B). Thus, posttranslational chromatin modifications corresponding to transcriptionally active and repressive marks follow the stem cell-like memory and effector gene expression patterns defined in memory precursors and effector T cell populations, respectively.

To investigate the contribution of Suv39h1dependent chromatin changes, we performed ChIP-seq for the repressive mark H3K9me3 in naïve and in vitro differentiated effector CD8 + T cells from wild-type and Suv39h1-KO mice (fig. S6A). The total numbers of H3K9me3 islands and targeted genes were higher in wild-type effectors than in wild-type naïve CD8 + T cells (Fig. 3B and fig. S6B). Most de novo H3K9me3 islands (present in effector but not in naïve T cells) were distal to the transcription start site (TSS), with equal proportions of intra-and intergenic locations (Fig. 3B, top), and with highest enrichment b e t w e e n1 0a n d1 0 0k bf r o mt h eT S S( f i g .S 6 C ) . Thus, H3K9me3 is broadly distributed in both wild-type naïve and effector CD8 + T cells, and the number of H3K9me3 islands increases in effector CD8 + T cells relative to naïve CD8 + T cells.

Among the 997 genes that are less efficiently silenced in vivo in Suv39h1-KO cells relative to littermates (Fig. 2C), 145 genes were decorated by H3K9me3 (fig. S6C). Several of these genes encode immune and stem/memory-related proteins, including CD127 or CD62L (Fig. 3C). In effectors, the Il7r gene acquires H3K9me3 at sites both proximal and distal to the TSS, including the promoter, previously described enhancers, introns, and intergenic regions (Fig. 3C, upper left) (26). H3K9me3 deposition is also increased at the Sell locus, upstream the promoter region, in effectors as compared to naïve CD8 + T cells (Fig. 3C, lower left). No significantly enriched H3K9me3 islands were found for actively transcribed housekeeping and control genes (Fig. 3C,right). These results suggest that H3K9me3 silences genes linked to naïve versus effector differentiation, including Il7r and Sell.

To evaluate whether these stem cell-like memory gene loci are direct targets of Suv39h1, we compared the relative (fold change) enrichment of H3K9me3 islands detected in effectors to stem cell/memory gene expression in wild-type and Suv39h1-KO K b -OVA + CD8 + T cells (Fig. 3D). H3K9me3 enrichment negatively correlated with mRNA expression of stem cell/memory signature genes in wild-type but not in Suv39h1-KO K b -OVA + CD8 + T cells (Fig. 3D). In contrast, none of the H3K9me3 domains were associated with effector signature genes. Thus, the silencing of stem cell-like memory genes correlated with significant H3K9me3 enrichment in wild-type but not in Suv39h1-KO T cells, in which increased gene expression correlated with reduced or absent H3K9me3 deposition.

To validate and refine these results, we used mChIP-quantitative polymerase chain reaction (qPCR) analysis of H3K9me3 at critical gene loci in naïve (CD44 lo CD62L + CD127 + ), central memory (CD44 hi CD62L + CD127 + ), and effector (dump - CD44 hi CD127 lo/-KLRG1 + K b -OVA + )C D 8 + Tc e l l s purified 7 days after LM-OVA infection. Naïve cells showed barely detectable H3K9me3 at both the Il7r and Sell promoters (Fig. 3E and fig. S6D), whereas in effector T cells, the levels of H3K9me3 were increased at both loci, again in correlation with gene silencing (Fig. 3E). Consistent with the expression results, central memory cells, like naïve cells, have low levels of H3K9me3 at both promoters, whereas the reduced expression of Il7r and Sell in effector cells correlated with H3K9me3 enrichment.

Suv39h1-KO effector CD8 + T cells did not show a significant increase of H3K9me3 at the Il7r or Sell loci, consistent with impaired silencing of Il7r and Sell expression (Figs. 2D and3E). Likewise, CD8 + central memory T cells and naïve cells presented low levels of H3K9me3 at both the Sell and Il7r promoters (Fig. 3E). Also consistent with these results, the level of the alternative active mark H3K9ac is increased in effectors at both Il7r and Sell loci in Suv39h1-KO CD8 + T cells alone, which was again consistent with gene e x p r e s s i o np r o f i l e s( F i g .3 F ) .T h e s er e s u l t si n d icate that Suv39h1 dynamically decorates genes encoding important regulators of CD8 + Tc e l l stem/memory fate with H3K9me3, and that these genes are silenced in wild-type CD44 hi CD127 lo/- KLRG1 + K b -OVA + CD8 + Tef fectors.InSuv39h1-KO antigen-specific CD8 + T cells, these stem/ memory genes fail to acquire the repressive mark, resulting in defective silencing.

Long-lasting memory and short-lived effector CD8 + T cell differentiation These gene expression results suggest that the Suv39h1 defect may affect memory versus effector differentiation in vivo. No major differences between wild-type and Suv39h1-KO were observed in terms of CD44, CD122, and PD1 expression patterns (fig. S7A) in blood K b -OVA + CD8 + T cells 7 days after LM-OVA infection. In contrast, CD127 expression was increased in both memory precursors (KLRG1 -) and effector cells (KLRG1 + )i nSuv39h1-KO mice relative to littermates (Fig. 4A). The increased proportion of CD127 + cells was due to decreased numbers of CD127 lo/-KLRG1 + effector T cells (Fig. 4, A and B, right). The absolute numbers of CD127 lo cells (both KLRG1 + and KLRG1 -) were reduced, whereas the number of CD127 high memory precursor cells was unchanged (Fig. 4, A and B, left). In Suv39h1-KO mice, K b -OVA + CD8 + T cells also included a lower proportion of CD27 lo cells and impaired down-regulation of CD62L in a subpopulation of cells (fig. S7, A andB). Enhanced memory differentiation of wild-type/ Suv39h1-KO 1:1 mixed bone marrow chimeras upon LM-OVA infection resulted in a significantly higher proportion of effectors among wild-type cells, whereas the proportion of memory precursors was increased among Suv39h1-KO K b -OVA + T cells (Fig. 4C and fig. S8). Similar resultswereobtainedaftertheadoptivetransfer of wild-type and Suv39h1-KO OT-I CD8 + T cells (fig. S3, D andE). Thus, Suv39h1-KO mice develop increased proportions of memory T cells in response to L. monocytogenes infection as the result of a T cell-intrinsic defect in effector differentiation.

We next sought to determine the impact of the Suv39h1 defect on long-term memory versus effector differentiation and persistence. The numbers of K b -OVA + CD8 + T memory cells and effectors were analyzed 7, 18, and 47 days after LM-OVA infection, and 3 days after rechallenge on day 51. As expected from previous results (Fig. 4B), the number of effectors was reduced at all times, including at the peak of the rechallenge memory response in the blood (Fig. 4D) and spleen (fig. S9, right). The numbers of memory cells did not change between wild-type and Suv39h1-KO at the peak of the response, during contraction, or after rechallenge (Fig. 4E An analysis of master regulators involved in memory and effector differentiation showed that the proportion of LM-OVA-specific T cells expressing T-bet is reduced in both IFN-g + and IFN-g -Suv39h1-KO cells as compared to littermates (Fig. 4I and fig. S11). Similarly, the expression of Eomes, Blimp1, and Bcl6 was also reduced (fig. S11). Thus, Suv39h1-KO CD8 + Tcells display a central memory-like phenotype but express reduced levels of both effector and memory transcription master regulators (4,5).

Previous adoptive transfer experiments showed that only CD127 + memory precursors give rise to long-term memory cells and confer protective immunity (29). To evaluate their in vivo memory self-renewal and differentiation properties, we isolated wild-type and Suv39h1-KO CD45.2 K b -OVA + CD8 + T cells 7 days after LM-OVA infection and adoptively transferred them at low numbers into naïve congenic CD45.1 recipients (Fig. 5A and fig. S12A). Forty days after challenge with LM-OVA, as expected, few wild-type K b -OVA + T cells persisted or responded to the LM-OVA challenge. In contrast, donor Suv39h1-KO CD45.2 K b -OVA + T cells were clearly present and responded to the infection (Fig. 5, B and C), with an increased proportion of memory precursor and effector cells (Fig. 5, B and D, left). To further evaluate the self-renewal properties of wildtype and Suv39h1-KO central memory T cells, we isolated total CD45.2 + CD44 hi CD62L + CD127 + KLRG1 -K b -OVA -central memory T cells, harvested 7 days after LM-OVA infection, and transferred them into naïve CD45.1 congenic mice (Fig. 5A and fig. S12A). Thirty-nine days after adoptive transfer, similar numbers of wild-type and Suv39h1-KO donor CD8 + T cells, which had maintained a central memory phenotype, were present in the blood (Fig. 5G, left, and fig. S12B). Four days after LM-OVA rechallenge, however, both the percent-ageandthetotalnumbersofSuv39h1-KO donor CD8 + T cells were increased relative to control Suv39h1-sufficient cells (Fig. 5, E and G, right). Donor Suv39h1-KO CD8 + T cells displayed higher expression of both CD127 and CD62L, lower levels of KLRG1 and CD44, a slight decrease in the expression of SCA-1, and similar levels of CD122 (Fig. 5, E andF). The donor memory subset was increased in the animals adoptively transferred with Suv39h1-KOTcells(Fig. 5H).Thus,Suv39h1-KO K b -OVA + and central memory CD8 + T cells have superior self-renewal and repopulation potential relative to their wild-type counterparts.

Stemness gene silencing in terminal effectors requires Suv39h1

Increased proportions of T cells with a central memory phenotype and stem cell-like properties were found in Suv39h1-KO mice. This phenotype could be due to the accumulation of a defined population of stem cell-like memory T cells, or to the expression of stem cell-related genes across different T cell subpopulations. We used singlecell RNA sequencing (scRNA-seq) to explore and dissect the heterogeneity of wild-type and Suv39h1-KO K b -OVA + CD8 + T lymphocytes. Purified naïve and K b -OVA + CD8 + T cells from wild-type and Suv39h1-KO mice were isolated by FACS 7 days after LM-OVA infection and processed for scRNAseq (fig. S13). For naïve cells, we sequenced 1102 and 991 cells from wild-type and Suv39h1-KO mice, respectively. For antigen-specific K b -OVA + CD8 + T cells, we processed two technical replicates for both wild-type and Suv39h1-KO infected mice (approximately 1200 and 1000 cells, respectively) and an additional biological replicate (from different mice, 404 wild-type and 283 Suv39h1-KO K b -OVA + CD8 + T cells). The cells from the two technical replicates were pooled and used for the rest of the analysis. A comparison between technical and biological replicates within wild-type and Suv39h1-KO mice showed a strong correlation (figs. S13 to S15).

A principal components analysis (PCA) of differentially expressed genes by wild-type and Suv39h1-KO naïve and K b -OVA + CD8 + T cells, is visualized as a set of t-distributed stochastic neighbor embedding (t-SNE) plots in Fig. 6A. In wild-type cells, unsupervised clustering of naïve and K b -OVA + CD8 + T cells revealed eight clusters in wild-type T cells and six clusters in Suv39h1-KO antigen-specific K b -OVA + CD8 + Tc e l l s( f i g . S14, A and B, and tables 3 and 4). On the basis of distinct transcription profiles, we grouped the unsupervised clusters into four major subset categories: naïve, memory precursor, effector, and cycling cells (Fig. 6, B to D, and fig. S14C). As expected, the naïve cells grouped in a homogeneous category characterized by the highest expression of Sell, Ccr7, and Tcf7 (Fig. 6B). The memory precursors were enriched in Il7r, Cxcr3, Cd27, Cd28,a n dLy6a expression (Fig. 6, B and E); the effectors were characterized by Zeb2, Klrg1, and granzyme B (Fig. 6,BandE,andfig.S16);and the cycling subsets were defined by cell cycle genes including Pcna and Mcm5 (Fig. 6, B and E, right). Although effectors represented the majority of multimer-positive CD8 + T cells when analyzed by FACS (CD127 -KLRG1 + ), the terminally differentiated effectors only represented approximately 30% of K b -OVA + CD8 + T cells in the categories defined by scRNA-seq (Fig. 6D). This is likely due to the resolution of the cycling cells, which also express high levels of different effector markers, including KLRG1 and granzyme A and B (Fig. 6, B and E, and fig. S16). Consistent with flow cytometric analysis, the proportion of Suv39h1-KO effector K b -OVA + CD8 + T cells was decreased (Fig. 6D and fig. S14). Thus, unsupervised scRNA-seq allows resolution of the expected memory precursors and effector populations among K b -OVA + CD8 + T cells. The scRNA-seq analysis also reveals a population expressing high levels of genes involved in the cell cycle, as well as markers of memory precursors and effectors.

Having defined the different populations of antigen-specific CD8 + T cells, we next sought to analyze the expression of stem cell/memory genes in wild-type and Suv39h1-KO cells. As expected, wild-type cells showed enrichment of stem cell markers (i.e., Cxcr6, Rnf138, Il18r1,andTraf1)in t hemem o ryprec u rso rs(F ig.6F ).I nSuv39h1-KO T cells, the expression of these markers was also high in memory precursors, but, in contrast to wild-type cells, expression was also detected in effector and cycling populations. Similar increased expression of the stem/memory signature in Suv39h1-KO effectors could be visualized on t-SNE plots where the number of genes from the stem/memory signature per cell is represented by a color scale (Fig. 6G). These results suggest that Suv39h1-defective effectors express higher levels of certain stem/memory-related genes, consistent with impaired silencing.

To further analyze the expression of the stem cell-like memory and effector signature genes at the single-cell level among the different populations, we used density scatterplots in which the numbers of genes from each signature are represented against each other (Fig. 6H). In most naïve cells from both wild-type and Suv39h1-KO mice, we detected an average of seven stem/ memory signature genes per cell, and fewer than four genes from the effector signature. As expected, in the wild type, memory precursors, relative to naïve cells, expressed higher numbers of genes from the stem/memory signature (mean of 9 genes per cell; 43% with more than 10 genes per cell) as well as high numbers of genes from the effector signature (mean of 18 genes per cell). The proportion of effector cells expressing high numbers (>10 per cell) of stem/memory genes was reduced (16%) as compared with memory cells (43%), whereas there was only a modest increase in the number of genes from the effector signature (mean 20 genes). Of note, the gene expression levels from the effector signature were strongly increased in effectors as compared to memory cells (Fig. 6B and fig. S16). Notably, cycling cells coexpressed either low or high numbers of genes from both signatures. Thus, scRNA-seq analysis reveals the alternative expression of the stem/memory and effector signatures in the two cell types, respectively, and the concomitant low or high expression of these two signatures in the cycling cells. These results suggest that cycling cells may represent bipotent differentiation intermediates expressing both effector and stem/memory potential. Furthermore, the commitment to effector differentiation paths appears to be acquired by the silencing of stem/memory genes. The total number of unique molecular identifiers (UMI) measured in each subset category did not differ between wild-type and Suv39h1-KO cells (fig. S17, A andB). Naïve, cycling, and memory Suv39h1-KO cells bear similar patterns of gene expression signatures as compared to wild-type cells. In contrast, a significant difference was observed in effector cells, in which the numbers of stem/ memory genes per cell were increased relative to effector cells from wild-type mice. The proportion of cells expressing more than 10 genes from the stem/memory signature was augmented from 16%inthewildtypeto34%ineffectorSuv39h1-KO cells. Thus, rather than a specific subpopulation of stem/memory cells accumulating in the Suv39h1-KO mice, the expression of stem/memory-related genes was derepressed mainly in Suv39h1-KO effector T cells.

Conclusions

We argue that after priming, cycling CD8 + T lymphocytes reprogram both self-renewing and effector gene expression profiles (Fig. 6I). These cycling cells may represent bipotent intermediates, which would then repress either the effector or stem cell/memory programs while they differentiate to memory precursors or effectors, respectively (Fig. 6I). The silencing of the stem cell/memory gene expression program is under the control of Suv39h1 by imposing the H3K9me3 modification on chromatin at the corresponding l o c i .I nd o i n gs o ,S u v 3 9 h 1 / H 3 K 9 m e 3w o u l de stablish an epigenetic barrier on the stem/memory gene expression program, preventing effector re-programing into memory cells (Fig. 6I These gene expression programs are intimately connected to precise control of chromatin dynamics, which determines the accessibility of genes and regulatory elements to transcription regulators [START_REF] Yadav | Chromatin plasticity: A versatile landscape that underlies cell fate and identity[END_REF]. Schematically, chromatin permissiveness to transcription is associated with chromatin modifiers including so-called writers and erasers, which biochemically modify histone tails, allowing the recruitment of effectors (readers) that modify the structure of the chromatin.

Suv39h1, a histone methyl transferase that promotes trimethylation of lysine 9 of H3 (H3K9me3) [START_REF] Rea | Regulation of chromatin structure by site-speci®c histone H3 methyltransferases[END_REF], plays critical roles in lymphocyte differentiation. In CD4+ T cells, Suv39h1 is associated with silencing of the IFNγ locus during Th2 differentiation [START_REF] Allan | An epigenetic silencing pathway controlling T helper 2 cell lineage commitment[END_REF]. In the absence of Suv39h1 In this study, we analyze the role of Suv39h1 in CD8+ T cell exhaustion and re-invigoration by PD-1 blockade in mouse tumors. Using both genetic and pharmacological approaches, we show that Suv39h1 inhibition overcomes tumor resistance to anti-PD-1 treatment. Immune phenotype and single cell RNAseq analysis show that Suv39h1-deficient TILs display phenotypic characteristics of exhausted cells, including high expression of multiple inhibitory checkpoints. However, the Suv39h1-defective TILs also express a strong IFN-I signature, and respond vigorously to PD-1 blockade, inducing potent tumor rejection and showing broader chromatin accessibility, in particular around genes linked to effector functions, as compared to wild type cells. We conclude that Suv39h1 expression in CD8+ T cells is instrumental for establishing an irreversible, anti-PD-1 resistant exhausted phenotype in tumors. As inhibition of Suv39h1 reverses T cell dysfunction and causes tumor rejection, Suv39h1 represents an actionable "epigenetic checkpoint" for T cell functions in cancer.

Results

Suv39h1 deficiency enhances solid tumor rejection in combination with anti-PD-1 Ab

To interrogate the role of Suv39h1 in T cell responses to cancer, we used a partially anti-PD-1 blockade resistant melanoma model, B16F10-OVA (in which the tumor cells express the surrogate tumor-antigen ovalbumin) grafted in WT or Suv39h 1-KO mice. Tumor bearing littermate WT and Suv39h 1-KO mice were randomly assigned to two groups, treated or not with anti-PD-1 monoclonal Ab after establishment of the tumors (Fig. 1A). As expected, anti-PD-1 treatment induces a delay in tumor growth, but not complete rejection in WT mice (Fig. 1B, 1C and 1D). In Suv39h 1-KO littermates, the growth of the tumor is similarly delayed, with a small number of tumor rejections (2/27), while treatment with anti-PD-1 induces complete tumor rejection in almost one-third of the mice (5/16, Fig. 1B, 1C and 1D). We conclude that Suv39h1 deficiency in the host not only delays tumor growth to levels similar to anti-PD-1 administration, but also synergizes with anti-PD-1 treatment.

The Suv39h 1-KO mice used here are constitutive knockouts affecting all cell types. To investigate the specific role of Suv39h1-deficient T cells in the observed tumor rejection phenotype, we adoptively transferred in vitro activated control WT and Suv39h 1-KO OVAspecific T cell receptor (TCR) transgenic OT-1 cells to B6 mice bearing established EL4-OVA tumors (Fig. S1A). Mice that received Suv39h 1-KO OT-1 T cells controlled tumor growth better than mice that received littermate WT OT-1 T cells (Fig. S1B), demonstrating that Suv39h1 deficiency in CD8+ T cells augments their anti-tumor potential in the absence of any other Suv39h1-defective cell type.

To further understand the implication of Suv39h1-deficient CD8+ T cells in tumor-rejection, we investigated how Suv39h1-defective T cells behave in a Graft vs Host Disease/ Graft vs Leukemia (GvHD/GVL) model. Irradiated haplotype-matched B6xDBA/2 F1 mice (Hk2bxd; termed B6D2) received a transplant of C57BL/6 (H2kb) bone marrow, with or without Suv39h1-knockout (KO) or wild-type littermate (WT) T cells, before being challenged with P815-GFP murine mastocytoma cells (H2kd) (Fig. S1C). Both Suv39h 1-KO and littermate WT T cells induced a modest increase in survival (10-30%), as compared to the control group that did not receive T cells (Fig. S1D). To determine if the mice died due to tumor or to GvHD development, tumor growth (Fig. S1E and1F) and GvHD severity (Fig. S1G and 1H) were measured. While mice injected with WT T cells showed a high tumor incidence (more than 60% had detectable tumor cells) and developed mild GvHD, the mice injected with Suv39h 1-KO T cells showed instead a lower tumor incidence (20%) and developed more severe GvHD. Indeed, while mice injected with WT T cells died from tumor development, the ones injected with Suv39h1-deficient T cells died from GvHD (78%) (Fig. S1I). Therefore, in this adoptive T cell transfer model Suv39h1-deficient T cells develop into much more potent effectors (that reject the tumor and attack the host) than their Suv39h1-proficient counterparts.

These two latter experiments (Fig. S1) suggest that the increased tumor control observed in Suv39h1-defective mice (Fig. 1A-D) could be attributed, at least partially, to an increased anti-tumor activity of Suv39h 1-KO CD8+ T cells. Along these lines, the percentages of tumor-infiltrating myeloid cells (CD11b+ cells, cDC1, and cDC2 cells), CD4+ T cells, and NK and were similar in untreated WT and Suv39h 1-KO animals, and were not significantly modified by the anti-PD-1 treatment (Fig. S2A); while the percentage and absolute number of CD8+ tumor-infiltrating lymphocytes (TILs) were increased in the Suv39h 1-KO mice as compared to WT littermates, and were expanded further after treatment with anti-PD-1 Ab (Fig. 1E). The increase in CD8+ T cells was restricted to the tumor, as no such changes were observed in the tumor-draining lymph nodes (DLNs) (Fig. S2B). These results highlight the link between increased CD8+ T-cell infiltration and tumor rejection in anti-PD-1 Ab treated Suv39h 1-KO mice.

It has been previously described that Suv39h1 is a critical regulator of the transition of CD8+ T cells between the memory and effector differentiation states (Pace et al., 2018). Thus, we analyzed the proportions of naïve (CD44loCD62Lhi), central memory (CD44hiCD62Lhi) and effector (CD44hiCD62Llo) CD8+ T cell subsets in blood (during tumor growth progression, at day 13), and in spleen, DLNs, and tumor (at sacrifice, day 20). We observed that in general, compared to WT mice, Suv39h 1-KO littermates treated or not with anti-PD-1, showed lower proportions of naïve CD8+ T cells and higher proportions of effector and central memory subsets, in the three locations (Fig. S2C-E). In contrast, in the tumor, virtually no naïve CD8+ T cells were detected, and the main effect was observed in PD-1 treated Suv39h 1-KO mice, in which the proportion of memory cells (CD44hiCD62Lhi) was significantly lower, while both the proportion and the absolute number of effector CD8+ T cells were significantly higher (Fig. 1F). Therefore, the absence of Suv39h1 expression not only leads to higher CD8+ T-cell tumor-infiltration upon anti-PD-1 treatment, but also a change in the phenotype of the infiltrating cells, away from memory and towards a more pronounced effector phenotype.

Suv39h1-KO CD8+ TILs display enhanced functional capacity

A vast body of literature shows that, upon tumor invasion, CD8+ T cells differentiate into an exhausted and unresponsive state in late-stage tumors (Phillip et al., 2017). Because tumors are rejected more efficiently in Suv39h1-defective mice treated or not with anti-PD-1, we hypothesized that the absence of Suv39h1 could delay or prevent exhaustion programs and together with checkpoint blockade could overcome dysfunctional states. To further investigate if Suv39h1-defective TILs are exhausted, we first analyzed the expression of negative-immune checkpoints, which represent a hallmark of exhaustion (Wherry and Kurachi, 2015). In Suv39h 1-KO mice, treatment with anti-PD-1 caused a marked increase in the proportion of TILs expressing PD-1, TIM-3, LAG-3 and 2B4 (Fig. 2A) and of all combinations of these receptors (Fig. 2B), as well as other surface molecules associated with exhaustion, such as CD39, CD38 and CD101 (Fig. S3A), suggesting more advanced differentiation (Philip et al., 2017;[START_REF] Hudson | Proliferating Transitory T Cells with an Effector-like Transcriptional Signature Emerge from PD-1+ Stem-like CD8+ T Cells during Chronic Infection[END_REF][START_REF] Zander | CD4+ T Cell Help Is Required for the Formation of a Cytolytic CD8+ T Cell Subset that Protects against Chronic Infection and Cancer[END_REF]. This increase was not observed in WT littermates treated in the same way.

To further study the effect of anti-PD-1 treatment on the activation state of Suv39h1-defective CD8+ TILs, we analyzed their effector functions, proliferation and survival, which have been previously described as impaired in exhausted T cells [START_REF] Paley | Progenitor and Terminal Subsets of CD8+ T Cells Cooperate to Contain Chronic Viral Infection[END_REF][START_REF] Pauken | Overcoming T cell exhaustion in infection and cancer[END_REF]. We observed that upon ex vivo stimulation CD8+ TILs from Suv39h 1-KO mice expressed the highest levels of GZMb, a marker of cytotoxic function (Fig. 2C), readily produced IFNγ, as detected by ELISPOT (Fig. 2D), and responded to IFN, as indicated by the high expression of ISG15 (a hallmark gene of IFN-response) (Fig. 2E).

Furthermore, around 13% of CD8+ TILs from PD-1-treated Suv39h 1-KO mice expressed higher levels of Ki67 (a marker of cycling cells, Fig. 2F). Nevertheless, increased cycling was not associated with augmented apoptosis, as suggested by the increased expression of BCL-2 (an anti-apoptotic molecule) (Fig. 2G) and small proportions of annexinV/7ADD + cells (Fig. 2H), suggesting higher function and conserved viability. To further understand this intriguing observation, we analyzed the expression of T-box transcription factors EOMES and TBET, both crucial markers in effector and memory functions in T cells [START_REF] Paley | Progenitor and Terminal Subsets of CD8+ T Cells Cooperate to Contain Chronic Viral Infection[END_REF][START_REF] Stelekati | Long-Term Persistence of Exhausted CD8 T Cells in Chronic Infection Is Regulated by MicroRNA-155[END_REF][START_REF] Khan | TOX transcriptionally and epigenetically programs CD8+ T cell exhaustion[END_REF], as well as the HMG-box transcription factor TOX, primary regulator of exhaustion epigenetic program [START_REF] Alfei | TOX reinforces the phenotype and longevity of exhausted t cells in chronic viral infection[END_REF]. We observed similar levels of TOX in CD8+ TILs from all mice conditions, and increased levels of EOMES in anti-PD-1 treated Suv39h 1-KO mice (Fig. S3B). Of interest, as cells get activated and enter into an exhaustion program, they upregulate the expression of Eomes; however, co-expression of EOMES with TBET has been associated with a rescue from a terminally exhausted state, as PD-1+EOMES+TBET+ cells retain effector functions. We observed that CD8+ PD-1+ TILs co-expressing EOMES and TBET where enriched in anti-PD-1 treated Suv39h 1-KO mice (Fig. 2I). These results indicate the ability to rewire the transcriptional control of TBET and EOMES in CD8+ TILs in Suv39h1 deficient mice after immunotherapy, independently of TOX expression. Thus, although Suv39h1-defective CD8+ TILs display features of exhausted cells (eg express high levels of negative checkpoints and EOMES), they are found in increased numbers (Fig. 1E), produce IFNγ and GZMb and importantly are not apoptotic, indicating enhanced poly-functionality and survival.

It has been recently described that anti-PD-1 Ab treatment may work not only by reinvigorating exhausted TILs, but also by inducing the proliferation of progenitor exhausted cells with memory/stem cell-like characteristics, as well as their differentiation into cytotoxic and short-lived, terminally exhausted cells resistant to anti-PD-1 therapy [START_REF] Miller | Subsets of exhausted CD8+ T cells differentially mediate tumor control and respond to checkpoint blockade[END_REF]Chen et al., 2019). We quantified the amount of progenitor and terminally exhausted T cells, identified as TCF1+TIM3-or TCF1-TIM3+ cells respectively, among the memory CD8+PD-1+ TILs compartment [START_REF] Miller | Subsets of exhausted CD8+ T cells differentially mediate tumor control and respond to checkpoint blockade[END_REF]. We observed that in the tumors of Suv39h 1-KO mice, the proportion of stem-like progenitor exhausted CD8+ T cells was lower than in their WT littermates, and the percentage of terminally exhausted T cells was significantly higher in anti-PD-1 treated Suv39h 1-KO CD8+ T cells compared with the WT mice (Fig. 2J), suggesting that in the absence of Suv39h1 activity CD8+ TILs differentiate from a stem celllike progenitor exhausted T cell to a more effector exhausted one. Similar results were obtained when using the cell surface marker SLAMF6 instead of the transcription factor TCF-1 [START_REF] Miller | Subsets of exhausted CD8+ T cells differentially mediate tumor control and respond to checkpoint blockade[END_REF] to quantify the precursor and the later effector populations (Fig. S3C).

Recent studies in chronic viral infection demonstrated that the expression of the type I transmembrane glycoprotein CD101 subdivide the terminally exhausted population into two subgroups [START_REF] Hudson | Proliferating Transitory T Cells with an Effector-like Transcriptional Signature Emerge from PD-1+ Stem-like CD8+ T Cells during Chronic Infection[END_REF]. Stem cell-like progenitor exhausted TCF1+CD8+ T cells first differentiate into a transitory CD101-TIM3+ state marked by an effector-like profile displaying chemokine receptor CX3CR1, cytokines and GZMb expression, contribute to viral control and eventually progress to a terminally differentiated and dysfunctional CD101+ TIM3+ state. Notably, Hudson et al. observed that PD-1 blockade increases the transitory cytolytic effector-like CD101-TIM3+ cells compared to stem-like CD8+ T cells. In the cancer context, we observed that anti-PD-1 treatment induced a subtle increase in the percentages of terminally exhausted CD101+TIM3+ cells in both WT and Suv39h1KO CD8+ TILs. The increase of transitory CD101-TIM3+ cells, in contrast, was stronger in Suv39h1-defective mice relative to the control littermates (Fig. 2K). These results indicate that in the absence of SUV39H1, anti-PD-1 treatment increases the proliferative and effector T cell differentiation by enhancing the conversion to transitory exhausted CD101-TIM3+ T cells.

scRNAseq analysis of tumor-infiltrating CD8+ T cells

To further characterize the impact of Suv39h1 on TIL programing, we performed single cell RNA sequencing (scRNA-seq) of CD8+ TILs isolated from B16F10-OVA tumors of Suv39h 1-KO and littermate WT mice treated or not with anti-PD-1 (Fig. 3A). In total, 21646 cells (1524 to 5953 cells per sample) were merged, and 8 clusters were defined, which are visualized using UMAP in Fig. 3B (see methods and Fig. S4 for bioinformatics pipeline details). Figure 3 shows the signatures used to define cluster identities (Fig. 3C), as well as key exemplifying genes (violin plots, Fig. 3D; and heatmap, Fig. 3E). Fig S4D -E shows additional information used to define each cluster. In more detail, cluster 1 and 2 were composed of cells bearing high memory signatures: "C1-memory 1" characterized by high expression of L ef 1 and Il7r ; and "C2-memory 2" by Id3, Ikzf 2/H eli os and killer cell lectinreceptors (Klr a1 and Klr a6); cluster 3, "C3-early activated" included cells with a signature of recent TCR-engagement, as featured by the expression of N r 4a1 and D usp2; cluster 4, "C4-IFN-a response" expressed an IFN-a-induced gene signature; cluster 5, "C5-effector cytolytic", shared similarities with the previously described short-lived effector T cell cytolytic cluster required to control chronic viral infection containing the chemokine receptor Cx3cr1, genes associated with cytotoxicity (G zm a, Klr g 1 and G zm b

)

Suv39h1-deficient TILs show an IFN-a response signature and a higher cytolytic potential upon PD-1 blockade

To further investigate how treatment with anti-PD-1 Ab affects reprograming of CD8+ TILs in WT and Suv39h 1-KO mice, we then explored the transcriptomic differences of the clusters between conditions. All clusters were present in each of the four analyzed conditions, but in different proportions. Thus, we first performed downsampling, so as to analyze the same amount of cells in each condition (randomly sampling 4097 cells from each condition to match the smallest one, n=16388 total cells after downsampling), and then quantified the changes in the proportions of cells belonging to each cluster from treated and untreated, littermate WT and Suv39h 1-KO mice. As shown in Figure 4A and 4B, CD8+ TILs from littermate WT untreated mice contained the highest proportions of memory and progenitor exhausted cells compared to all the other conditions. Upon anti-PD-1 treatment, these populations seemed to evolve into an exhausted state. Distinctively, in untreated conditions Suv39h 1-defective TILs were characterized by higher proportions of "IFN-responding" and "exhausted" cells than WT cells. Upon anti-PD-1 treatment, Suv39h 1-KO CD8+T cells comprised the highest frequencies of early activated, effector cytolytic and IFN-responsive cells. These results suggest that along the activation/differentiation program of anti-tumor T cells, Suv39h1 dependent silencing of the IFN-I response and of pathways of the exhaustion program, imposes an epigenetic barrier that prevents re-programing of TILs towards a terminally exhausted program by anti-PD-1. Moreover, our data indicate that in the absence of Suv39h1, CD8+ TILs are poised to more effectively respond to TCR-activation and to become highly cytolytic upon PD-1-blockade, underlying the increased tumor rejection compared to WT mice.

The anti-PD-1 response signature of Suv39h1-deficient CD8+ T cells discriminates response to anti-PD-1 in melanoma patients

We then asked whether the signature of Suv39h 1-KO CD8+ T cells that increase after anti-PD-1 in the mice is relevant to human cancer. To address this question, we first interrogated whether our data correlated with available published signatures from scRNAseq of CD8+ T cells from melanoma patients responding or not to immune checkpoint blockade (Sade-Feldman et al., 2018). To that end, we first confirmed that the "responder" and "nonresponder" signatures from the original paper correspond to distinct subpopulations of cells in the human dataset (Fig. S5A, S5B and Fig. 4C). To compare mice and human anti-PD-1 responses, we computed ortholog-matched signatures of each experimental group (littermate WT, littermate WT + anti-PD-1, KO, and KO + anti-PD-1) and quantified them on the reprocessed human melanoma cohort UMAP (Fig. 4D). We observed that the KO + anti-PD-1 signature best matched with the "responder" signature from the original study, and that they were strongly positively correlated (R = 0,54, p < 2.2e-16) (Fig. 4E). Interestingly, the KO + anti-PD-1 signature was also anti-correlated with the "non-responder" signature (R = -0,39, p < 2.2e-16). Finally, we also noticed an enrichment of the KO + anti-PD-1 signature coming from our mice data in responding patients (p=4.8 x 10 -12 ) (Fig. S5C), suggesting that our KO + anti-PD-1 signature could be a potential biomarker of response to immune checkpoint blockage (ICB) therapy.

We next used our mouse scRNAseq data to identify the specific subtype of CD8+ T cell response induced by the treatment in the responder patients. To that end, we first plotted signatures of our different conditions (littermate WT, littermate WT + anti-PD-1, KO, KO + anti-PD-1) (Fig. 4F) 4G). While each of these four previous studies differ with regards to detailed patient clinical characteristics, RNAseq profiling techniques and analytical approaches, we observed that across the four studies, the human signatures of response to immunotherapy were consistently associated to the cytolytic effector, TCR activation and IFN type I clusters, which are also those upregulated in the mouse KO + anti-PD-1 CD8+ T cells, and to a lesser degree to an increase of a fraction of "terminally exhausted" cells.

and
Overall, the transcriptional program induced by anti-PD-1 in the absence of Suv39H1 activity is associated with human T cells responding to PD-1 blockade and clearly links the emergence of effector CD8+ T cell subpopulations to the effectiveness of immune checkpoint blockade.

Suv39h1-deficiency enhances chromatin accessibility induced by anti-PD-1

To investigate the effect of anti-PD-1 reprograming on chromatin accessibility in WT and Suv39h 1-KO T cells, we next used Assay f or Tr an sposase-Ac c essi b le Ch r om ati n wi th h i g h th r oug h put seq uen c i n g (ATAC-seq). CD8+ TILs were FACS-sorted from tumors growing in WT and Suv39h 1-KO mice, treated or not with anti-PD-1. Accessible regions were first identified by calling ATAC-peaks for each condition. Tens of thousands of peaks were detected in each sample, with numbers systematically higher for anti-PD-1-treated cells, both in WT and Suv39h 1-KO CD8+ T cells, relative to untreated samples (Fig. 5A). To get quantitative insight into changes caused by anti-PD-1 treatment, we first compared treated samples directly to their untreated genotype-matched counterparts. An ATAC peak (identified in either condition) was considered modulated if it harbored twice or more sequencing reads (normalized for sequencing depth) as compared to the corresponding counterpart (Fig. 5B, left panel). Strikingly, with this approach no regions emerge as enriched in anti-PD-1 untreated relative to treated samples. On the contrary, anti-PD-1 treatment induces a widespread increase in chromatin accessibility in both WT and Suv39h 1-KO populations. In WT cells, almost 6000 regions at least double their ATAC-signal, as compared to their untreated counterparts. In Suv39h 1-KO cells, the effect of treatment is even more dramatic, with 35000 genomic regions gaining accessibility (almost half of all identified ATAC peaks in KO cells). In WT CD8+ TILs, the majority (68%, 4009/5866) of treatment-induced regions (T-peaks) are shared with Suv39h 1-KO T-peaks (Fig. 5B, left panel). ATAC signal density plots generated for 20% of the most variable T-peaks (Fig. 5B, right panel) reveal welldefined clusters: (i) Suv39h 1-KO -specific, (ii) WT-specific, and (iii) common to both genotypes. To cross-compare accessibilities of identified treatment-specific peaks in all four conditions, we plotted the distribution of ATAC signal in peaks for each cluster separately. As expected, average signals in WT or Suv39h 1-KO treatment-specific peaks are higher, as compared to the corresponding untreated counterparts (Fig. 5C). Importantly, a stronger signal in Suv39h 1-KO cells on shared T-peaks suggests that they are more accessible in a higher proportion of cells, as compared to WT. These results indicate that anti-PD-1 treatment promotes chromatin accessibility, and that this effect is more pronounced in the absence of Suv39h1.

To reveal biological signatures associated with anti-PD-1 treatment, we next assigned treatment-specific peaks to their nearest genes. We observed that the majority of genes assigned to T-peaks in WT cells are shared with Suv39h 1-KO (Fig. 5D, left panel), but conversely, the majority of T-peaks in Suv39h1-KO cells are not differentially accessible in WT cells. Signature enrichment analysis showed that anti-PD-1 treatment leads to an increase in chromatin accessibility in genomic loci linked to T cell memory, pluripotency and activation (TCR pathway, IL-2 pathway, IFNγ-response, and exhaustion), in both WT and Suv39h 1-KO CD8+ T cells (Fig. 5D, right panel). The progenitor signature, however, is more enriched in WT, while the IFNγ-and IFN-I signatures are more enriched in KO cells. Moreover, multiple effector signatures (cytolytic, terminally differentiated, as well as cell cycle) are uniquely associated with T-peaks in Suv39h 1-defective cells. Figure 5D, lower panel shows loci of representative TCR mediated activation gene EG r 1, effector genes (Cx3c r 1 and G zm b ) and IFN-stimulated gene M x1, that become more accessible in Suv39h 1-KO cells after anti-PD-1 treatment as compared to other conditions. These results suggest that altogether, anti-PD-1 treatment favors chromatin opening in genomic loci linked to T cell activation, memory and pluripotency, but in a Suv39h 1-KO context cells acquire preferential accessibility in functional IFN-related and cytolytic effector loci as compared to WT cells.

To investigate the differences between WT and Suv39h 1-KO populations after anti-PD-1 treatment, we next compared chromatin opening in the two populations. We only retained ATAC regions harboring twice or more reads and, importantly, called as peaks uniquely open in one, but not the other population after the treatment (Fig. 5E, upper panel). Of these peaks, around 1500 in WT and 1700 in KO are identified as condition-specific (see density plots in Fig. 5E with selected adjacent genes). Peaks uniquely accessible in treated WT cells (when compared to treated KO cells), are equally accessible in untreated WT population, (see ATAC signal distribution in Fig. 5E, lower panel). Only 3% of them are modulated by the treatment, suggesting that these are mostly stably opened regions in WT TILs. Peaks specific to treated Suv39h 1-KO cells (compared to treated WT cells), are also accessible in untreated KO samples, although over half of them (58%) respond to the treatment with a strong increase in accessibility (Fig. 5E, lower panel). Therefore, in the absence of Suv39h1, a series of chromatin regions, mostly treatment-sensitive, fail to remain closed. As ATAC-regions are typically associated with active or poised regulatory elements [START_REF] Buenrostro | Transposition of native chromatin for fast and sensitive epigenomic profiling of open chromatin, DNA-binding proteins and nucleosome position[END_REF], these results collectively suggest that Suv39h1 may contribute to the regulation of chromatin opening at cis-regulatory elements that support differentiation of T cells into functional cytolytic effectors upon re-activation.

Pharmacological inhibition of Suv39h1 phenocopies the genetic defect

Since genetic deficiency for Suv39h1 increases the response to anti-PD-1 treatment, we next evaluated the effect of ETP-69, a pharmacological inhibitor of Suv39h1 activity. ETP-69, which is not strictly specific for Suv39h1 (it also inhibits G9A, another mono-and dimethylase of H3K9) has been previously reported to function in mice [START_REF] Overman | Preparation of Epipolythiodioxopiperazine ETP derivatives for treatment of cancer[END_REF][START_REF] Snigdha | H3K9me3 Inhibition Improves Memory, Promotes Spine Formation, and Increases BDNF Levels in the Aged Hippocampus[END_REF]. Here, we evaluated its effect in combination therapy with PD-1blockade. For this, C57BL/6 mice bearing established B16F10-OVA tumors were treated daily with ETP-69 combined or not with anti-PD-1 (Fig. 6A). ETP-69 and anti-PD-1 administered as monotherapies induced a clear delay in tumor growth, and the combination of both was more effective than either of the single agents (Fig. 6B and6C). The increased efficiency of the combination treatment was ETP-69 dose-dependant (Fig. S6A andS6B). We conclude that pharmacological inhibition of Suv39h1 augments the efficacy of anti-PD-1 blockade in this tumor model.

We observed that frequencies and absolute numbers of CD8+ TILs were increased in ETP-69 treated mice, as compared to control groups, especially after treatment with anti-PD-1 (Fig. 6D). Similar to the results obtained in Suv39h1-defecient mice, ETP-69 alone or in combination with anti-PD-1 induced an increase in the proportion of effector CD8+ TILs, paralleled by a reduction in the proportion of memory TIL populations (Fig. 6E). In blood and secondary lymphoid organs, ETP-69 treatment lead to a reduction of naïve CD8+ T cells and increased proportions of central and effector memory subsets (Fig. S6A). Furthermore, CD8+ TILs from ETP-69 treated mice expressed multiple inhibitory checkpoint receptors (Fig. 6F and6G), and higher proportions of GZMb and IFNγ upon restimulation (Fig. 6H), and of Ki67+ cells (Fig. 6I). Finally, as observed for Suv39h1-KO mice treated with anti-PD-1, a higher proportion of CD8+ PD-1+ TILs from ETP-69 + anti-PD-1 treated mice coexpressed EOMES and TBET (Fig. 6J). We conclude that, similar to the Suv39h1-genetic defect, treatment of WT tumor-bearing mice with the Suv39h1 inhibitor ETP-69 reduces tumor growth, especially in combination with anti-PD-1. The inhibitor also enhances tumor infiltration by effector CD8+ T cells that express high levels of multiple inhibitory checkpoints but with an increased effector cytotoxic phenotype. Therefore, pharmacological inhibition phenotypically and functionally phenocopies the genetic defect of Suv39h1.

Discussion

In the tumor microenvironment, conditioned by hypoxia, nutrient shortage and tumor factors, "chronic" activation of CD8+ T lymphocytes progressively leads to loss of their effector functions and unresponsiveness to further stimulation, a process often referred to as "final exhaustion". Although the signaling pathways involved in the induction of exhaustion are incompletely understood, the transcriptional programs expressed in exhausted cells have been analyzed in some details. Induction and silencing of hundreds of genes are under the control of critical transcription factors such as IRF4, cMAF, NFAT, TOX, EOMES and T-bet and others [START_REF] Martinez | The Transcription Factor NFATPromotes Exhaustion of Activated CD8+ T Cells[END_REF][START_REF] Chihara | Induction and transcriptional regulation of the co-inhibitory gene module in T cells[END_REF]. Similar to any gene expression program, transcriptional regulation of exhaustion and re-invigoration upon immune checkpoint blockade can only be promoted under a precise control of chromatin accessibility at the corresponding loci. Here, we identify Suv39h1, the main H3K9 tri-methylase and a hallmark of facultative and constitutive heterochromatin [START_REF] Nielsen | Heterochromatin formation in mammalian cells: interaction between histones and HP1 proteins[END_REF][START_REF] Wiencke | Differentially expressed genes are marked by histone 3 lysine 9 trimethylation in human cancer cells[END_REF], as a key epigenetic positive regulator of functional exhaustion.

Based on our results, we propose a working model in which, Suv39h1 contributes to progression of TILs along the exhaustion pathway through heterochromatin-mediated silencing of IFN-I gene expression programs and pathways downstream of TCR-signaling, inhibiting the deployment of an effector program. Thus, as cells become exhausted, silencing of these gene expression programs establishes an epigenetic barrier that prevents effective reprograming by PD-1 blockade. In the absence of Suv39h1, this epigenetic barrier is incomplete and TCR-triggering in memory/effector cells leads to an alternative activation pathway that allows survival and differentiation into highly cytolytic effectors without reaching final exhaustion. This in turn allows more effective re-invigoration of exhausted TILs by anti-PD-1 Ab, and more efficient tumor rejection How then does Suv39h1 repress TCR activation, terminal differentiation and ISG expression programs? Our previous work showed that several stem cell/memory related genes display reduced levels of H3K9me3 deposition in Suv39h 1-KO T cells, suggesting that they are direct targets of Suv39h1 (Pace et al., 2018). Chromatin analyses in TILs, however, are complicated because of reduced cell numbers and heterogeneity among individual mice. Interestingly, enhanced IFN-I signature genes are also overexpressed in hematopoietic stem cells [START_REF] Wu | Intrinsic Immunity Shapes Viral Resistance of Stem Cells[END_REF], suggesting a shared link between the mechanisms that control the two gene expression programs. The link could be direct, such as through deposition of the same histone marks (including H3K9me3) at memory and ISGs loci. H3K9me3 could also control ISGs indirectly, for example through regulation of expression of transposable elements (TEs), which in turn drive IFN-signaling. Indeed, Suv39h1 was shown to control expression of several families of TEs in fibroblasts and RNA from TEs can be sensed in the cytosol and induce activation of STING and IFN-I [START_REF] Sidler | A role for SUV39H1-mediated H3K9 trimethylation in the control of genome stability and senescence in WI38 human diploid lung fibroblasts[END_REF].

To investigate the mechanisms of Suv39h1 action, we sequenced single cell transcriptomes from over 21000 infiltrating CD8+ T cells from B16F10-OVA tumors in WT, or Suv39h1defective mice, treated or not-treated by anti-PD-1 Ab. Our single cell analysis captured a high heterogeneity of CD8+ TILs, encompassing lymphocytes undergoing different states of activation: cycling, memory, memory effectors, progenitor exhausted, recently-TCR activated, IFN-responding, and cytolytic effectors subpopulations; but as expected no naïve cells. The frequency of cycling cells detected by this technique was not increased in anti-PD-1-treated mice (WT or Suv39h 1-KO, Fig. 4B), while, at least in the KO, anti-PD-1 induced increased proportions of Ki67-positive cells. This could be due to limitations of the available transcriptomic signatures of cycling, as in our single cell analysis, the cluster of cycling cells includes cells in different phases of the G2M and S phases, which seem not to change in a systematic way among the different conditions (Fig. S4B andS4C). The changes observed in the relative proportion of cells belonging to the other clusters indicate that during the antitumor response, Suv39h1 represses IFN-I signaling, and part of the exhaustion program, restraining the cells in a more memory/progenitor exhausted state, and T cell reinvigoration induced by PD-1 blockade drives these cells into a terminally exhausted program. However, in the absence of this epigenetic barrier imposed by Suv39h1, TCR-mediated reactivation is restored, and leads to highly cytolytic effectors with anti-tumor potential. Consistent with this working model, signatures from scRNAseq CD8+ T cells from melanoma patients that respond to immune checkpoint blockers correlated strongly with the transcriptomic signature of PD-1-treated Suv39h 1-KO CD8+ T cells. Reciprocally, the PD-1-treated Suv39h 1-KO signature was enriched in responding patients, indicating that this signature could potentially be used as biomarker of response in melanoma patients. Furthermore, the high granularity of CD8+ T cell populations characterized in scRNAseq analysis of the mouse data, allowed us to identify which were the specific CD8+ T subpopulations responding to ICB in patients, using as input published melanoma ICB response signatures obtained from single cell or bulk RNA sequencing. Of note, this analysis identified that, as in PD-1-treated Suv39h 1-KO mice, response to immune checkpoint blockade was consistently associated to the emergence of cytolytic effector and TCR activation-; a feature not previously highlighted in patients-, and also IFN type I signatures, as previously described [START_REF] Liu | Integrative molecular and clinical modeling of clinical outcomes to PD1 blockade in patients with metastatic melanoma[END_REF].

Our results identify a critical epigenetic regulator of exhaustion. They also establish a clear link between IFN-I, early activation, cytolytic effectors and re-programing of progenitor exhausted cells. This link may be particularly relevant to the multiple immunotherapy approaches based on inducing IFN-I responses in the tumor environment, including cytolytic viruses, STING agonists, and DNA-demethylating agents. Since Suv39h1 genetic deficiency or pharmacological inhibition both re-activate anti-tumor immune responses, and promote reprograming by anti-PD-1, the role of Suv39h1 in this process must be at least partially nonredundant. By analogy to other non-redundant immunosuppressive proteins whose inhibition unleashes anti-tumor immune responses (often referred to as "immune checkpoints"), Suv39h1 can be considered as an "epigenetic immune checkpoint". Should the role of Suv39h1 be conserved in human T cells, its blockade would open new perspectives for the epigenetic manipulation of anti-tumor T cell responses in the clinic.

supported by a fellowship from the International PhD program of Institut Curie. The authors declare no competing financial interests.

INTRODUCTION

Lung cancer is the leading cause of cancer-related death worldwide. Non-small cell lung cancer (NSCLC) represents 85% of lung cancer diagnoses (1), and tumor infiltration by lymphocytes is associated with a favorable survival prognosis (2) and a better clinical response to immune checkpoint blockade (ICB) (3). ICBs are thought to "reprogram" CD8 + tumor-infiltrating lymphocytes (TILs) to produce antitumor responses (4,5) by targeting inhibitory receptors such as programmed cell death protein 1 (PD-1), which are highly expressed by most TIL populations. How this reprogramming is achieved and which TIL subpopulations are targeted by ICBs are still open questions.

Most studies have distinguished early precursors and terminally differentiated CD8 + T cell populations (6,7) in chronic viral infections and tumors. Early precursors present characteristics of memory cells, including memory/stem markers and regenerative capacity. They are characterized by high expression of the chemokine recep-tor CXCR5, which is also expressed by B cells and follicular helper T cells (Tfh) (8,9). These precursors are not fully committed and can therefore be "reprogrammed" by ICB (9). Terminally differentiated populations are clonally expanded and express higher levels of effector markers and immune checkpoints (ICPs)-including PD-1, T cell immunoglobulin and mucin domain-3 (TIM3), and CD39 (6, 9)-and share many characteristics with "exhausted" or "dysfunctional" CD8 + T cells. Dysfunctional TILs were initially described in chronic viral infection mouse models by poor effector function, expression of inhibitory receptors, low proliferation, and distinct transcriptional and epigenetic states compared with effector or memory T cells (10,11). Several groups have recently analyzed T cell transcriptional programming in cancer, revealing a strong heterogeneity among TILs (12)(13)(14)(15)(16)(17), similar to chronic viral infections, with distinct subpopulations of exhausted progenitors and terminally differentiated T cells. The idea that terminally differentiated T cells in cancer are dysfunctional, however, is still debated, especially in human tumors.

Here, we used a combination of single-cell RNA sequencing (scRNA-seq) and single-cell T cell receptor (TCR) sequencing (scTCR-seq) in tumors, normal tissues adjacent to the tumor (juxtatumor), and peripheral blood to delineate the processes of CD8 + T cell differentiation in patients with untreated, primary NSCLC. We show that precursor, memory-like CD8 + TILs are composed of two main populations: one is also present in the blood (circulating precursors) and the other is also present in juxtatumor tissue and bears markers of memory-resident T cells (resident precursors). Both precursor subtypes differentiate into a main population of terminal effectors through a similar "transitional" stage. Terminal effectors are not observed in blood or juxtatumor tissue, are more clonally expanded, and express signatures of exhaustion. A substantial proportion of transitional and terminal effectors also express cell cycle signatures and Ki67, suggesting that clonal expansion that occurs in situ is part of the terminal differentiation process.

RESULTS

Transcriptional states of CD3 + TILs in NSCLC

We performed scRNA-seq on CD3 + T cells from 11 primary, untreated, and early-stage resected patients with NSCLC (Fig. 1A and table S1). For six of these patients, matched scTCR-seq information was also collected. Combined analysis of all patient samples identified a large source of variation associated with 3′ versus 5′-oriented single-cell chemistries (fig. S1,A to C) that was resolved using computational methods (18). Rigorous assessment of intrasample heterogeneity, including TCR, showed that integration did not overcorrect the sample expression profiles. After removing contaminating cell types, we collected data on 28,936 cells over the 11 patients (Fig. 1B). Uniform manifold approximation and projection (UMAP) and unsupervised graph-based clustering partitioned cells into 21 clusters based on their transcriptome (Fig. 1, B andC). Differential gene expression analyses (fig. S1D and table S2) and analyses of canonical marker genes (Fig. 1D) revealed the identities of the major cell clusters. We used gene signatures (table S3) (17,(19)(20)(21)(22)(23)(24) to reduce the impact of sparsity in the scRNA-seq values and assigned identities to different clusters, including naïve, memory, regulatory, helper, and effector T cells (Fig. 1E). Among these clusters, CD4 + regulatory T cells (T regs ) (CD4-IL--32-T regs ) and CD4 + memory (CD4-CD69-activated memory) are the most abundant, whereas the most abundant CD8 + subsets are CD8 + -circulating precursors [CD8-Kruppel-like factor 2 (KLF2)] and transitional CD8 + [CD8-granzyme H (GZMH)] (Fig. 1C and fig. S1E). We identify 10 clusters of CD4 + cells [including C-C motif chemokine receptor 8 (CCR8)-T regs , interleukin-32 (IL-32)-T regs , GZMA effectors, SELL-naïve, CD69-activated memory, IL-7R-memory, myelin and lymphocyte protein (MAL)-T regs , tumor necrosis factor receptor superfamily member 18 (TNFRSF18)-Tfh, sestrin 1 (SESN 1)-Tfh, and heat shock protein family H member 1 (HSPH1)-memory], 7 clusters of CD8 + cells [including FCGR3A effectors, KLF2-circulating precursors, GZMK-circulating precursors, XCL1-resident precursors, LAYN-terminally differentiated, GZMH-transitional, and SLC4A10-mucosal-associated invariant T cells (MAIT)], a cluster of gd Τ cells (TRDC), and 3 clusters of both CD4 + and CD8 + cells: one characterized by high expression of interferon (IFN)-related genes [CD4/8-IFN-stimulated gene 15 (ISG15)] and the other two characterized by strong expression of cycling genes (CD4/8-MCM5 and CD4/8-TOP2A). Analysis of ICP molecules expression showed that CD4-CCR8 cells express both inhibitory and stimulatory receptors, with particularly high levels of stimulatory molecules TNFRSF18 and TNFRSF4, whereas CD8-LAYN cells show higher expression of inhibitory molecules such as hepatitis A virus cellular receptor 2 (HAVCR2)/ TIM3 and lymphocyte activating 3 (LAG3) (fig. S1F). All clusters are present in all patients, although at varying frequencies (fig. S1G). As CD8 + T cells are critical players during tumor rejection, we focus our analysis on this compartment.

Memory-like precursors and terminally differentiated CD8 + TILs

CD8 + T cells (10,243 cells after selection) segregate mainly along a horizontal axis in the UMAP, with cells expressing a stemness signature (Fig. 1E), central memory-related genes (including IL7R, CXCR5, and TCF7; Fig. 2A), and oxidative phosphorylation signatures (also related to T cell memory; Fig. 2A) at the left end of the axis. At the right end of the axis, cells express higher levels of effectorrelated genes (Fig. 1E), negative ICPs (Fig. 1E), including HAVCR2/ TIM3 and ENTPD1 (Fig. 2A), TOX (a transcription factor involved in exhaustion, Fig. 2A), and a glycolysis signature (also associated with effector function, Fig. 2A). Consistent with these results, exhausted precursor (GMZK/ZNF683) and terminally exhausted (LAYN) signatures from Guo et al. (13) are also highly expressed on the left and right ends of this axis, respectively. Comparison of the expression of signatures from a study on chronic viral infection in mice (19) in our data showed progenitor and terminal exhausted signatures with enhanced expression on the left and right ends of this axis, respectively (fig. S2, A andB). These results suggest that CD8 + T cell clusters on the left of the UMAP are "memory-like precursors" (clusters CD8-GZMK, CD8-KLF2, and CD8-XCL1) and that those on the right part of the UMAP (CD8-GZMH and CD8-LAYN) are more differentiated and most likely related to late exhausted/dysfunctional cells.

Memory-like precursors include two main subclusters, CD8-GZMK/CD8-KLF2 and CD8-XCL1, expressing high levels of the GZMK and ZNF683 signatures from Guo et al. (13), respectively (Fig. 2A). Memory-like precursors are enriched in signatures of negative regulation of cell proliferation, whereas terminally differentiated cells are enriched for lymphocyte activation-related genes and adhesion integrins (fig. S2B). We conclude that early-stage NSCLC TILs include two main populations of memory-like precursors (CD8-GZMK/CD8-KFL2 and CD8-XCL1) and two populations of late, differentiated cells (CD8-GZMH and CD8-LAYN) that may correspond to terminally exhausted TILs.

Two distinct circulating precursor subsets

Previous studies have associated precursor CD8 + TILs to a favorable prognosis and better response to ICB-based immunotherapies (17). To investigate the possible relevance of the different clusters defined here to clinical outcomes, we compared our dataset to a previous single-cell study from Sade-Feldman et al. (17) and applied signatures from melanoma patients responding (good response) versus not responding (bad response) to ICB. CD8-KLF2-circulating and CD8-XCL1-resident precursor clusters express higher levels of the good response signature (Fig. 2, B andC). In contrast, CD8-LAYN and CD8-GZMH late clusters, and a part of CD8-GZMK early clusters, express higher levels of the poor response signature (Fig. 2, B andC).

CD8-KLF2 and CD8-GZMK show similar global expression profiles but are discordant for the ICB prognostic signatures; thus, we performed differential analysis of gene expression between these two populations. Cells in the CD8-GZMK cluster overexpress genes related to T cell activation [human leukocyte antigen (HLA)-DR/DP/DQ], effector-related genes [killer cell lectin like receptor G1 (KLRG1), GZMB, GZMH, and GZMA, in addition to GZMK], and chemokines (CXCR6 and CCL5) (Fig. 2, D andE). The CD8-KLF2 cluster overexpresses genes related to the control of CD8 + T cell differentiation, including NR4A1/2/3 and KLF2 (Fig. 2, D andE). These two closely related populations of CD8 + TILs are also present in other public datasets (13,16,17), including some generated on different technological platforms [switching mechanism at the end of the 5′-end of the RNA transcript sequencing (SMART-seq2)/massively parallel RNA single-cell sequencing (MARS-seq)] and in other tumor types, such as melanoma and breast carcinoma (fig. S2C). This supports the existence of subtle but important differences between these two populations and motivates further characterization of their relationship with respect to signaling and differentiation pathways.

Tissue-resident and transitional CD8 + populations

To better understand the different early and late TIL populations, and because one of the main driving genes among the three early clusters is ZNF683 [also known as Hobit, a transcription factor involved in the programming of tissue-resident memory CD8 + T cells (25)], we hypothesized that some clusters could correspond to tissue-resident cells and others to circulating cells. As shown in Fig. 3A and fig. S3A, single gene expression of ZNF683, as well as a core tissue residency gene signature, consisting of four main tissueresident markers (ITGAE, ITGA1, CXCR6, and ZNF683), are all expressed at higher levels in late-differentiated CD8-GZMH/ CD8-LAYN and in early CD8-XCL1 clusters, as compared with the CD8-GZMK/CD8-KLF2 clusters, which are characterized by enriched expression of KLRG1 (fig. S3A). These mRNA expression differences were validated at the protein level. Flow cytometry analysis of six additional patients with NSCLC showed that a higher proportion of CD103 + CD8 + TILs express markers related to T cell dysfunction, including PD-1, TIM3, and CD39, as compared with CD103 -TILs (Fig. 3, B to D and fig. S3B). The subset of CD103 + CD8 + TILs not expressing markers of T cell dysfunction, or expressing them at very low levels, likely corresponds to the CD8-XCL1resident precursors that we identified in the scRNA-seq analysis, because transcript levels for inhibitory ICPs are low in this cluster (fig. S1F). A large fraction of CD103 + CD8 + TILs coexpresses T cell dysfunction markers and GZMB while showing lower expression of KLRG1, as compared with CD103 -CD8 + TILs (Fig. 3C). These results are consistent with recent studies (26-28) that found subsets of tissue-resident (CD103 + ) CD8 + TILs enriched in checkpoint receptors and display features of enhanced cytotoxicity and tumor reactivity, indicative of a therapeutic potential of this population. Analysis of cluster membership as a function of clustering resolution showed that the split between CD8-XCL1 and CD8-KFL2/GZMK occurs early, indicating the importance and robustness of the differences (fig. S3C). On the basis of these results, we hypothesized that the two early, memory-like populations of CD8 + cells could correspond to tissue-resident cells (CD8-XCL1) and recent emigrants (most likely from blood, CD8-KFL2/CD8-GZMK; see below). We further speculated that the two types of precursors could both show converging differentiation into late-differentiated and dysfunctional/exhausted (CD8-GZMH/CD8-LAYN) subtypes.

CD3 + T cell selection Expression level Min Max C D 8 -F C G R 3 A C D 8 -K L F 2 C D 8 -G Z M K C D 8 -X C L 1 C D 8 -G Z M H C D 8 -L A Y N C D 8 -S L C 4 A 1 0 C D 3 -T R D C C D 4 /8 -T O P 2 A C D 4 /8 -M C M 5 C D 4 /8 -I S G 1 5 C D 4 -S E L L C D 4 -IL -7 R C D 4 -H S P H 1 C D 4 -C D 6 9 C D 4 -G Z M A C D 4 -S E S N 1 C D 4 -T N F R S F 1 8 C D 4 -I L 3 2 C D 4 -M A L C D 4 -C C R 8
C D 8 -K L F 2 C D 8 -G Z M K Expression level NR4A1 0 1 2 3 4 C D 8 -K L F 2 C D 8 -G Z M K NR4A2 0 1 2 3 C D 8 -K L F 2 C D 8 -G Z M K NR4A3 0 1 2 3 4 C D 8 -K L F 2 C D 8 -G Z M K KLF2 0 1 2 3 C D 8 -K L F 2 C D 8 -G Z M K SOCS3 0 1 2 3 4 C D 8 -K L F 2 C D 8 -G Z M K Expression level GZMA 0 1 2 3 4 5 C D 8 -K L F 2 C D 8 -G Z M K NKG7 0 2 4 C D 8 -K L F 2 C D 8 -G Z M K CCL5 0 1 2 3 4 C D 8 -K L F 2 C D 8 -G Z M K GZMB 0 1 2 3 C D 8 -K L F 2 C D 8 -G Z M K PRF1 Min Max Expression level -0.5 0.0 0.5 1.0 -0.2 -0.1 0.0 0.1 0.2 C D 8 -F C G R 3 A C D 8 -K L F 2 C D 8 -G Z M K C D 8 -X C L 1 C D 8 -G Z M H C D 8 -L A Y N Good response
To test this working model, we used two different unsupervised approaches to infer continuous transitions between clusters. Statistical analysis of connectivity in the k-nearest neighbor graph of cellcell expression similarities [partition-based graph abstraction (PAGA)] ( 13-15, 29, 30) shows that the CD8-GZMH cluster represents a transitional state between early and late clusters (Fig. 3E). Because discretizing elements in a continuum of differentiation might not be optimal, we also used pseudotime alignment to resolve relationships between continuous populations. Monocle3 analysis revealed a converging differentiation process, consistent with the original UMAP representation and with PAGA analysis (Fig. 3F). As anticipated earlier, pseudotime reconstitution distinguishes separate resident and circulating precursor differentiation branches (Fig. 3F). Examples of differentially expressed genes between branches are shown in fig. S3D. With both tools, the CD8-GZMH cluster localized at the intersection between three branches (CD8-XCL1, CD8-KFL2/CD8-GZMK, and CD8-LAYN), suggesting that the CD8-GZMH cluster corresponds to a transitional population of cells undergoing differentiation from early, memory-like (CD8-XCL1 and CD8-KLF2/CD8-GZMK) to late, terminally differentiated (CD8-LAYN) clusters.

If the CD8-GZMH cluster represents a transitional state, then it should present lower clustering robustness, because cells from dif-ferent clusters would enter this "state" and then differentiate into other states. We quantified robustness of clustering using the silhouette score (Fig. 3G and fig. S3E), which validates consistency by evaluating how close each cell inside a cluster is to its neighboring cells within the same cluster (high score), compared with its neighboring cells in other clusters (low score). As shown in Fig. 3G, the silhouette score is lower for cells in the CD8-GZMH cluster, as compared with cells in all other CD8 + clusters. Consistent with this analysis, label transfer, a method introduced for the integration of multiple single-cell datasets (18), of cells from the CD8-GZMH cluster shows even projection to the corresponding neighboring clusters, including CD8-GZMK (upper left part of the UMAP), CD8-XCL1 (lower left), and CD8-LAYN (to the right) (Fig. 3H and fig. S3F). Although our clustering algorithm identified one transitional GZMH cluster for statistical robustness, this is a heterogeneous population and subtle differences are likely to exist within this subpopulation. The likely dynamic nature of this transition process makes it challenging to track experimentally in human clinical samples. These results, as hypothesized earlier, are consistent with early precursor CD8 + populations (CD8-GZMK, CD8-KLF2, and CD8-XCL1) differentiating into terminally dysfunctional CD8-LAYN cells, via a transitional CD8-GZMH state.

Clonal sharing between clusters informs differentiation

To better understand the ontological relationships between the different TIL clusters, we sought to analyze their TCR repertoires. scTCR-seq for paired a and b chains from six patients was obtained, with matched 5′ RNA profiling, for about 80% of the cells (table S4), including cells harboring unique or shared TCRs, indicative of clonal expansion. TCR clonotypes were called using the 10X pipeline when the analysis was restricted to one sample and a custom algorithm when the samples were distinct but autologous (31,32). The clone size is defined as the total number of cellular barcodes associated with the same clonotype. Top clones cluster identities for each patient are shown in fig. S4A.

To investigate the extent of clonal expansion in different clusters, we first represented individual cells in colors according to the size of their TCR clones. As shown in Fig. 4 (A andB), larger clones are found within CD8 + T cell clusters as opposed to CD4 + and, more specifically, in late-differentiated (CD8-GZMH and CD8-LAYN) clusters, consistent with previous reports (13-15, 29, 30). This pattern was conserved across patients (Fig. 4B, bottom), as measured using either a published expansion index (12,14) or by the mean number of cells per TCR clonotype. By visualizing all cells from individual TCR clones in the UMAP, we found that TCR clones are confined to specific clusters in most cases (examples in Fig. 4C).

To quantify whether there is a link between individual TCRs and specific clusters, we calculated the conditional probability that if a particular TCR is observed in one cluster, then it will also be observed in the same or other clusters. Unsupervised clustering of these probabilities is shown in Fig. 4D. The probability that the same TCR is found in cells from the same clusters is consistently higher than the probability to find the same TCR in cells from different clusters (dark blue diagonal). Nevertheless, there is an observable degree of TCR sharing across clusters and the probability to find the same TCRs in two different clusters is not zero. As illustrated by the light blue cases in Fig. 4D, "sharing" of TCRs between clusters occurs and is preferential between clusters with similar functional abilities. For example, the subclusters of CD4-T regs share more TCRs among them than they share with CD4-Tfh or CD4-memory T cell clusters. Among the CD8 + T cell clusters, increased sharing occurs between circulating precursor memory-like clusters (CD8-GZMK and CD8-KLF2) and between differentiated clusters (CD8-LAYN and CD8-GZMH). The CD4/8-ISG15 (enriched for IFN-related genes) and the cycling clusters show increased sharing with early memorylike and late-differentiated clusters, respectively. We conclude that while high-intracluster TCR sharing highlights the relevance of our clustering strategy, intercluster TCR sharing can be used to infer transitions of cells between clusters: Higher TCR sharing between two clusters might indicate dynamic transitions of cells between the two transcriptional states. We conclude that intercluster TCR sharing can be used to infer a biologically relevant measure of proximity between clusters that can indicate either a dynamic transition between two transcriptional states or a static similarity not fully captured by transcriptome-based clustering methods.

To investigate potential transitions between clusters based on TCR sharing in more detail, we represented all the TCRs shared between two clusters in a "circos" representation (Fig. 4E), in which each line represents a single TCR clonotype. This representation illustrates that the two main circulating precursor memory-like clusters (CD8-GZMK and CD8-KLF2) and the two main late-differentiated (CD8-GZMH and CD8-LAYN) clusters share numerous TCRs. We also computed the transition index using the single T-cell analysis by Rna-seq and Tcr TRACking (STARTRAC) method (33), as the likeliness of a cluster to share clones (fig. S4B), and we see that CD8-LAYN and CD8-GZMH are the clusters that are most likely to have TCR sharing. The circos plot shows that GZMH displays higher clonal sharing with all the rest of the clusters (Fig. 4E and fig. S4C).

AB CD E F

Terminally differentiated Precursors This result was reproducible between patients. The CD8-GZMH cluster shares over half of its shared TCRs with the CD8-LAYN cluster, whereas the other half is shared with other CD8 + T cell clusters (fig. S4C). Supervised representation of the computed probabilities of shared TCRs among clusters (Fig. 4F) shows that the CD8-GZMH cluster could represent a central hub, sharing TCRs with both early, memory-like, and late-differentiated clusters. The high levels of TCR sharing between clusters were confirmed by analysis of the number of T cells from each cluster for shared clones. As shown in fig. S4D, "equilibrated" sharing between clusters, when the same TCR is present in three or more cells from each cluster, is only seeing between CD8-GZMH/CD8-LAYN and CD8-GZMK/CD8-KFL2 and, to a lower extent, between XCL1/GZMH. Equilibrated sharing suggests active transitions between clusters. As proposed below, TCR-sharing analysis is consistent with the GZMH cluster representing a transitional state between the two main populations of memory-like precursors (CD8-KLF2/CD8-GZMK and CD8-XCL1) and the late-differentiated cells present in the CD8-LAYN cluster.

CD8-LAYN CD8-GZMK CD8-GZMH CD8-XCL1 CD4-TNFRSF18-Tfh CD4-CCR8-T regs 3 2 1 0 -1 -2 Number of cells Expansion index C D 8 -K LF 2 C D 8 -G Z M K C D 8 - G Z M H C D 8 -L A Y N C D 8 -X C L 1 C D 4 / 8 -I S G 1 5 CD4/8-TOP2A C D 4/ 8-M C M 5 C D 4 / 8 -T R D C C D 4 / 8 -F C G R 3 A

Terminally differentiated, not memory-like, CD8 + T cells are in cell cycle

The question of the cycling activity of early memory-like versus late terminally differentiated TILs is still a matter of debate (34). As T cells require TCR signaling to enter the cell cycle, clones identified to be actively cycling intratumorally strongly suggest a local intratumor source of antigen. Cell cycle-related genes are expressed at high levels, and thus, cycling cells cluster independently in scRNA-seq analyses. These clusters, however, may include cells that also bear underlying signatures from other clusters, which are in some way masked by highly expressed cell cycle genes. We took two independent approaches to investigate whether the cells in the cycling cluster are related to other clusters.

We focused first on the transcriptomic dataset, and clustering of the infiltrating T cells in 11 patients with NSCLC shows two different cycling populations, one in G 2 -M phase (CD4/8-TOP2A) and the other one in S phase (CD4/8-MCM5) (fig. S5A). To better characterize these clusters, we used label transfer to interrogate the "second best" cluster to each cycling cell. The IFN-related cluster (CD4/8 ISG15) also includes CD4 + and CD8 + T cells, and ISGs are also highly expressed and can drive independent clustering of cells related to other clusters, so we included this cluster in the same label transfer analysis. Label transfer of CD8 + cells from the CD4/8 ISG15 cluster shows that they recluster mainly to CD8-GZMK and CD8-GZMH clusters (Fig. 5, A andB). Label transfer of CD8 + T cells from the cycling clusters results in almost exclusive reclustering to late-differentiated effector clusters CD8-GZMH and CD8-LAYN (Fig. 5, A andB). Almost no cells from the cycling clusters are reattributed to progenitor CD8 + T cell populations, suggesting that cycling cells are transcriptionally closer to late-differentiated T cells.

We validated these results using flow cytometry analysis in freshly isolated TILs from 10 additional patients with NSCLC. As shown in Fig. 5 (C andD), CD8 + cells coexpressing PD-1, TIM3, and CD39 [triple-positive (TP)] are also Ki67 + . CD8 + T cells negative for the three inhibitory ICP markers or expressing only PD-1 show very low Ki67 labeling (Fig. 5, C andD). t-distributed stochastic neighbor embedding (tSNE) visualization of cytometry data shows that a tissue-resident (CD103 + ) subset of these proliferative TP CD8 + cells also expresses PD-L1, suggesting that this subpopulation may play an important role in ICP response (fig. S5B).

The LAYN signature has previously been associated with exhausted/dysfunctional CD8 + T cells (13,33,35,36). Among cycling cells, the ones transcriptionally related to the CD8-LAYN cluster by label transfer also express higher levels of DNA damage repair genes (fig. S5C). These cells express higher levels of gene signatures for DNA damage repair and stress (Fig. 5E) and of the terminally differentiated LAYN signature (fig. S5D) (13), as compared with cycling cells reattributed to the CD8-GZMH cluster. Analysis of the dataset of Guo et al. (13) for the cycling signature shows that the LAYN and the cycling populations overlap (fig. S5E), suggesting that late differentiating cells, and not early progenitors, are the main cycling T cells in NSCLC TILs.

We reasoned that if late differentiating cells cycle, then clonally amplified TCRs in these clusters should also be found in cycling T cells. Circos representation shows that the TCRs expressed in cycling cells are found preferentially in cells from late-differentiated effector CD8-GZMH and CD8-LAYN clusters (Fig. 5F). As shown in fig. S5F, preferential TCR sharing of cycling cells with late-differentiated clusters is not due to the larger size of these late clones. CD8-GZMK/ CD8-KLF2 clones of more than 20 cells still share very few TCRs with the cycling cluster, as compared with CD8-GZMH/CD8-LAYN clones of similar sizes. Increased sharing can also be visualized in the sharing probability graph shown in Fig. 4D. Therefore, cycling cells share TCRs preferentially with late effectors, as compared with early memory-like precursors.

The results presented thus far are consistent with the possibility that transition between CD8-GZMH and LAYN clusters occurs while cells divide. To test this hypothesis, we measured the frequency of cycling cells among clones shared between the different clusters. As shown in Fig. 5G, most of the shared clones that include GZMH and LAYN cells also have TCRs in the cycling clusters (more red central dots in the left panel as compared with the right panel, Fig. 5G). Figure S5G uses a similar representation to show that the TCR clones found in the GZMH, LAYN, and cycling clusters include TCRs also present in the XCL1 cluster more frequently than TCRs also present in the GZMK cluster.

These results suggest that the preferential pathway toward terminal differentiation in tumors originates in the CD8-XCL1 cluster and transitions through CD8-GZMH to CD8-LAYN while cell divide and clonally expand. These data further support a model in which many expanded CD8 + TIL clones are actively cycling within the tumor, likely recognizing local sources of antigen in the tumor microenvironment. Furthermore, among CD8 subpopulations, the CD8-GZMH and CD8-LAYN clusters are potentially the most enriched for such local antigen specificity, as they show the strongest evidence for active cycling. It should be noted that this does not mean that tumor-specific TCRs are completely restricted to the GZMH and LAYN clusters nor does it guarantee that all local antigens are tumor specific. This model is consistent with these cells being chronically stimulated and with recent results in chronic human inflammation (37).

Resident versus circulating origin and tissue distribution of TILs

We analyzed scRNA-and TCR-seq in TILs from blood and juxtatumor tissue to further investigate the tissue versus blood origins of the different TIL populations described thus far. We obtained blood samples from four of the previously described patients and juxtatumor tissue samples from two patients. scRNA-and scTCR-seq in all four patients resulted in 54,247 validated cells, with ~80% RNAand TCR-coupled data. After removing contaminating cells and integrating the three different tissues, we used the clusters that have been identified before in the 11 different NSCLC tumor samples as a reference (as shown in Fig. 1B).

SCIENCE IMMUNOLOGY | RESEARCH ARTICLE 9 of 15 C D 8 - X C L 1 C D 8 -K L F 2 CD4 /8-M CM 5 C D 4 /8 -IS G 1 5 C D 4 /8 -T O P 2 A CD8-LAYN C D 3-T R D C CD8-F CGR3 A C D 8 -G Z M H C D 8 - G Z M K Cycling IFN Label transfer AB C IFN Cycling CD8-KLF2 CD8-GZMK CD8-XCL1 CD8-GZMH CD8-LAYN CD4-SELL CD4-IL-7R CD4-HSPH1 CD4-CD69 CD4-GZMA CD4-TNFRSF18 CD4-IL-32 CD4-CCR8 0
We quantified the proportions of the different populations in the different tissues by integrating cells from the three different tissues in these patients and used label transfer to map the populations from blood and juxtatumor to the tumor reference (Fig. 6,A to C). This allowed us to focus on conserved patterns between our tumor reference and the other tissues and minimize oversmoothing induced by integrating data from different tissues (fig. S6A). The CD4-SELL cluster (naïve T cells) contains mostly cells from blood, whereas the CD8-LAYN cluster (terminally differentiated) contains almost exclusively cells from tumors (Fig. 6, A to C, and fig. S6B). This result is consistent with flow cytometry analysis in four additional patients with NSCLC, showing that TIM3 + CD39 + CD8 + cells are enriched in the tumor samples, as compared with blood samples, which contain mainly single-positive PD-1 + or TN CD8 + cells (fig. S6C).

The CD8-XCL1 cluster (resident precursors) contains mostly cells from tumor and juxtatumor tissue (not from blood), indicating that it is a tissue-specific population, whereas cells in the CD8-GZMK and CD8-KLF2 clusters are also present in blood. These results are consistent with the proposed dual origin of memory-like progenitors: CD8-KLF2/CD8-GZMK cells originate from blood and could represent circulating precursors, whereas CD8-XCL1 cells absent from blood, but in juxtatumor tissue, may represent tissue-resident precursors (Fig. 6, B andC). The most abundant CD8 + T cell cluster in blood is CD8-FCGR3A, which are also present in juxtatumor samples but are rare in tumors (Fig. 6, B andC). Consistent with our model, CD8-GZMH (transitional) and CD8-LAYN (terminal) clusters are phenotypes acquired in tissue (Fig. 6B). The absence of CD8-LAYN cells in blood and juxtatumor is also consistent with CD8-GZMH cells becoming CD8-LAYN only in the tumor microenvironment.

We used scTCR information to further test this model. Most TCR expanded clones in the three tissues are present within CD8 + T cell populations (Fig. 6D). The clusters with the most expansion, however, are distinct in each tissue: the CD8-LAYN cluster in tumors (as shown before), the CD8-XCL1 cluster in juxtatumor, and the CD8-FCGR3A (together with the CD8-GZMK/CD8-KLF2 clusters) in blood (Fig. 6D). Analysis of shared TCR clonotypes present in tumor and blood (Fig. 6E) shows extensive TCR sharing between CD8-GZMK and CD8-KLF2 clusters in the two locations. Many expanded CD8-GZMK TCRs in blood are expressed in CD8-KLF2 cells in tumor, suggesting that CD8-KLF2 TILs may be derived from CD8-GZMK blood cells.

To further investigate possible ontological associations between clusters in blood and tumor, we analyzed the top 20 TCR clones per cluster in the blood (Fig. 6F, top) or in the tumor (Fig. 6F, bottom) for their transcriptional programming and TCR numbers in the other tissue (blood for tumor and tumor for blood). The top 20 TCR clonotypes from late-differentiated tumor clusters (CD8-GZMH and CD8-LAYN) were not observed in the blood, consistent with these clones expanding intratumorally. Among the top 20 CD8-XCL1 TCR clonotypes from the tumor, only one was also found in blood, where it has CD8-GZMK transcriptional programming. This result is consistent with most CD8-XCL1 clones not coming from blood but from tissue-resident origin. Among the top 20 tumor clonotypes in the CD8-KLF2 and CD8-GZMK clusters, all are found at substantial frequencies in blood. However, whereas top tumor CD8-GZMK clonotypes are also mapped to CD8-GZMK cells in blood, top tumor CD8-KLF2 clonotypes are also found in CD8-GZMK clustered cells in blood (as suggested by the circos analysis; Fig. 6E). Analysis of the top 20 clonotypes per cluster from the blood yields consistent results. First, most expanded clones in blood are in the CD8-GZMK cluster, and most of these clones are also found in the tumor, where they display transcriptomic reprogramming corresponding to both CD8-KLF2 and CD8-GZMH clusters, mainly. This result suggests that TCR expanded CD8-GZMK blood CD8 + T cells are the main blood precursors for tumor infiltration. This also suggests that after infiltration, blood TCR expanded CD8-GZMK precursors are reprogrammed to other phenotypes including CD8-KLF2 and CD8-GZMH.

DISCUSSION

We used scRNA-and scTCR-seq to analyze CD8 + TILs, juxtatumor tissue, and blood in early-stage resected patients with NSCLC. On the basis of integration of transcriptomic programming with TCR in the three tissues, we propose an integrative working model for TIL origin, differentiation, and functional organization in primary NS-CLC (Fig. 6G). Although we used multiple complementary approaches to derive this model, it should be noted that they are based on computational inference from static molecular snapshots and that directly assessing such inherently dynamic processes requires longitudinal sampling of individual patients (29).

The most abundant and clonally expanded population of CD8 + T cell effectors in blood (CD8-FCGR3A) is rare in tumors as are the clonal TCRs they express, indicating that these cells either do not infiltrate tumors effectively or die rapidly upon infiltration. If the cells entered tumors and changed phenotypes, then their TCRs would still be present in other clusters. This population also specifically expresses CX3CR1 (Fig. 1D), which has been used to identify a population of TILs with some markers of effector activity (13)(14)(15)33). The position of this population with respect to TIL differentiation pathways has been unclear, with some studies describing it as an intermediate state in a linear pathway to exhaustion (15), others as an alternative end state to exhaustion in a branched pathway (13,33), and yet others as an independent population not linked to intratumoral T cell differentiation (14). Our results, including trajectory analyses, TCR sharing, and cross-tissue comparisons, are most consistent with the final model. CD8-GZMK cells, in contrast, are abundant and clonally expanded in both blood and in tumors and share similar TCRs, suggesting that they represent the main blood precursor for TILs. A large proportion of the TCRs present in the CD8-GZMK cluster are also present in tumors in a second cluster of memory-like cells, CD8-KLF2, suggesting that these two populations are in dynamic equilibrium. Although trajectory reconstitution suggests that these two populations can differentiate into transitional CD8-GZMH cells, TCR sharing between the two is relatively low as compared with the sharing between CD8-XCL1 and CD8-GMZH. This result suggests that the transition between CD8-GZMK/CD8-KLF2 and CD8-GZMH is slow or inefficient. The other main population of early, memory-like, CD8 + T cells is in the CD8-XCL1 cluster. These cells express tissue-resident signatures and markers and are absent from blood, suggesting their tissue-resident origin. These cells share numerous expanded TCRs with CD8-GZMH cells, suggesting that they represent the main source of GZMH transitional cells. Previous studies have described TILs with markers shared by our CD8-XCL1 cluster not as a precursor population but rather as either an intermediate differentiation state (13) or as a final resident memory state in a branched differentiation pathway (33). One potential reason for these differences is that these previous studies used plate-based approaches for scRNA-seq, which generally profiles fewer cells per sample than the droplet microfluidics technology we used. This difference gives our dataset more power in terms of TCR repertoire analysis, both for quantifying expansion and intercluster sharing.

The effector-like, late CD8 + TILs are divided into two clusters: CD8-GZMH and CD8-LAYN. Cells in the CD8-GZMH cluster show limited robustness in clustering and share numerous clonal TCRs with all other CD8 + T cell clusters, indicating that early memory like cells (circulating or tissue resident) transit through this stage before entering late terminally differentiated, most likely dysfunctional states (CD8-LAYN). Determining whether consistent, subtle differences exist between cells within the GZMH cluster, for instance, depends on whether they were derived from tissue-resident or circulating precursors and is unclear with our current dataset. Cells in the CD8-LAYN cluster belong to the larger TCR clones, most of which are shared with the CD8-GZMH cluster. Consistent with the idea that transition between these two late clusters requires cell divisions, a large proportion (~40%) of these late-stage cells shows strong cell cycle signatures and are labeled by Ki67 antibodies. Intratumor expansion of these cells, together with expression of T cell activation markers, suggests that the cells recognize antigens intratumorally, whether tumor specific or not. Tissue-resident CD8-XCL1 cluster cells share most TCRs with the CD8-GZMH/CD8-LAYN clusters, so we hypothesize that tumor antigen-specific T cells are derived mostly from tissue-resident memory populations (rather than from clones recently stimulated by antigen in the tumor-draining lymph nodes). The highly expanded TCRs present in blood (mainly in CD8-GZMK cluster cells) remain in early/memory-like clusters and are rare in late clusters, suggesting inefficient intratumor progression from blood precursors to terminally differentiated cells.

We did not investigate whether cells in the CD8-LAYN cluster are active effectors or dysfunctional/terminally exhausted in this study. Several recent papers provide partial evidence in each direction, depending on the types of tumors and markers used (14,19,28,(36)(37)(38).

Here, we show that the CD8-LAYN cluster expresses high levels of a "bad response" signature to ICP in melanoma. We also found that cycling cells with LAYN-signature expression also highly express endoplasmic reticulum (ER) stress and DNA repair signatures, as compared with both noncycling cells from the same cluster and to cycling cells from the CD8-GZMH cluster (Fig. 5E). These results are consistent with cycling LAYN cells being exhausted/dysfunctional, rather than actively involved in tumor rejection.

In contrast to CD8-LAYN cells, very few memory-like CD8 + T cells were found to bear cell cycle signatures. In mice, several studies show that progenitor, memory-like TILs cycle (36)(37)(38), but the findings with human cells are more unclear. Recent papers suggest that TILs in early dysfunctional state in melanoma cycle (14), whereas other studies in NSCLC are consistent with our results, showing that more terminally differentiated cells that highly express inhibitory ICP molecules cycle more than early progenitors (12). The working model for TIL infiltration and differentiation proposed here would predict that TCR expansion in the CD8-LAYN cluster, most likely exhausted or dysfunctional T cells, will not indicate a good prognosis for clinical responses to ICP blockers. Clinical responses to ICP blockade should be distinguished from the ability of ex vivo expanded TILs to recognize autologous tumor. Our results suggest that the CD8-GZMH and CD8-LAYN clusters are likely to be enriched for tumor-specific TCRs, consistent with reports that these populations show the highest proportion of tumor specificity by such experiments (12,14). Among early, memory-like, and potentially reprogrammable CD8 clusters, our results suggest that TCR expansion in CD8-XCL1, rather than CD8-KLF2/CD8-GZMK, cluster cells could suggest a good prognosis for responses to ICB. Future scRNA-seq studies in patients with NSCLC responding or not to ICB will test this hypothesis.

MATERIALS AND METHODS

Study design

The overall objective of this study was to characterize the diversity and ontogeny of CD8 + T cells infiltrating untreated lung tumors. This was accomplished by performing scRNA-seq and TCR-seq and flow cytometry analysis on T cells isolated from tumor, juxtatumor, and blood samples from patients undergoing surgical resection for early-stage lung cancer. All samples were collected from the Institute Mutualiste Montsouris, under a dedicated protocol for lung cancer specimens approved by the French Ethics and Informatics Commission (EUdract 2017-A03081-52). All patients in this study provided written informed consent for sample collection and data analysis. Cohort size was selected to assess interpatient variability in subpopulation levels, and the number of cells per sample was selected to robustly estimate the frequency of subpopulations of at least 10% abundance.

Human specimens

Twenty-one patients, who were pathologically diagnosed with NSCLC, were enrolled in this study. Eleven were profiled by single-cell sequencing technologies, including 10 patients with adenocarcinoma and 1 patient with squamous cell carcinoma. The remaining 10 patients were profiled by flow cytometry analysis. All patients were untreated, with early-stage disease. Tumor tissue samples were obtained from all 21 resected patients with NSCLC. For four of the patients profiled by scRNA-seq, paired peripheral blood was collected and analyzed. Among those four, two patients had matched adjacent normal lung tissue collected and analyzed. For five of the patients profiled by flow cytometry analysis, paired peripheral blood was collected and analyzed. Among those five, two patients had matched adjacent normal lung tissue collected and analyzed. All samples were collected from the Institute Mutualiste Montsouris, under a dedicated protocol for lung cancer specimens approved by the French Ethics and Informatics Commission (EUdract 2017-A03081-52). All patients in this study provided written informed consent for sample collection and data analysis. All clinical information is summarized in table S1.

Tumor tissue and adjacent normal lung tissue were obtained from surgical specimens after macroscopic examination of the tissue by a pathologist. Tissue samples were stored in CO 2 -independent medium (Invitrogen) with 10% human serum and transferred within 1 hour after surgery to the research institute. For each specimen, a fragment was formalin-fixed and paraffin-embedded for histology and immunohistochemistry.

Tissue dissociation

Tumor and adjacent normal lung tissue samples were gently cut in about 1-mm 3 pieces. Tissues were digested enzymatically by a 20-to 40-min incubation, based on the size of the tissue, at 37°C under agitation, in CO 2 -independent medium (Invitrogen) with collagenase I (2 mg/ml; Sigma-Aldrich), hyaluronidase (2 mg/ml; Sigma-Aldrich), and deoxyribonuclease (25 mg/ml; Sigma-Aldrich).

The tissue pieces were gently dissociated with a 20-ml syringe plunger on a 40-mm cell strainer (BD) in 1× phosphate-buffered saline (PBS) (Invitrogen) with 1% fetal bovine serum (FBS) and 2 mM EDTA (Gibco) until uniform cell suspensions were obtained. The suspended cells were subsequently centrifuged for 10 min at 400g.

TILs isolation

TILs were isolated using Ficoll-Paque PLUS solution (Sigma-Aldrich). After tissue digestion, cells were resuspended in CO 2 -independent medium and layered onto Ficoll-Paque PLUS solution. Cells were centrifuged for 20 min at room temperature at 400 relative centrifugal force (RCF) without breaking. After centrifugation, TILs were carefully transferred to a new tube and washed with 1× PBS with 1% FBS and 2 mM EDTA (Gibco).

Peripheral blood mononuclear cell isolation

Peripheral blood mononuclear cells (PBMCs) were isolated using Ficoll-Paque PLUS solution (Sigma-Aldrich). Fresh peripheral blood was collected before surgery in EDTA anticoagulant tubes. Eight milliliters of fresh peripheral blood was resuspended in 1× PBS (Invitrogen) with 1% FBS and 2 mM EDTA (Gibco) and layered onto Ficoll-Paque PLUS solution. Cells were then centrifuged for 20 min at room temperature at 2000 rpm without breaks. After centrifugation, TILs were transferred to a new tube and washed with 1× PBS with 1% FBS and 2 mM EDTA (Gibco).

CD3 + T cell isolation and purification

CD3 + T cells from all samples were enriched using magnetic positive selection; CD3 microbeads and MACS separation (Miltenyi Biotec). Subsequently, dead cells and debris were removed, following the manufacturer's instructions (dead cell removal kit and debris removal kit, Miltenyi Biotec) resulting in ~80% of purity and CD3 + T cells were resuspended in 1× PBS with 0.04% bovine serum albumin. Cell numbers and viability were measured using a Countess II Automated Cell Counter (Thermo Fisher Scientific) as well as classical hemocytometer and trypan blue.

scRNA-seq and TCR profiling

Single-cell suspensions were loaded onto a Chromium Single Cell Chip (10X Genomics) according to the manufacturer's instructions for coencapsulation with barcoded gel beads at a target capture rate of 5000 to 10,000 individual cells per sample, based on the initial number of cells per sample. Referring to the blood samples, the target capture rate was 10,000 individual cells, whereas for tumor tissue and normal adjacent lung tissue samples, it was diverse between 5000 and 10,000 individual cells, due to the different cell number per tissue. For patients P34, P35, P42, P43, and P46, captured mRNA was barcoded during cDNA synthesis and converted into pooled scRNA-seq libraries for Illumina sequencing using the Chromium Single Cell 3′ Solution (10X Genomics) according to the manufacturer's instructions. For patients P47, P55, P57, P58, P60, and P61, RNA and TCR libraries were synthetized by following the Chromium Single Cell 5′ V(D)J Enrichment Kit, Human T Cell (10X Genomics).

Flow cytometry

Flow cytometry analysis was performed in a dataset of 10 additional untreated patients with NSCLC, in early disease stage. For the isolation of lymphocytes from tumor tissue and normal adjacent lung tissues and of PBMCs, the same steps predescribed above were followed. For negative selection of CD3 + T cells, we used the Pan T cell Isolation Kit (Miltenyi Biotec) that leads to about 80% purity. Gating strategy: Cells were gated avoiding doublets and debris. Dead cells were excluded by surface staining with a LIVE/DEAD fixable Zombie NIR viability kit (1:1000; BioLegend, lot number B262784) for all experiments. CD8 + T cells were stained with a combination of surface markers. Exclusive antibody list is included in Supplementary Materials and Methods. Gating was applied by using monocolors and fluorescence minus one (FMO) controls. Data analysis was performed using FlowJo v.10.6.1. Intracellular staining was performed using the eBioscience intracellular fixation and permeabilization buffer set (Thermo Fisher Scientific) according to the manufacturer's guidelines. All samples were analyzed using a Cytek Aurora spectral flow cytometer.

Software versions

Data were collected using Cell Ranger software (10X Genomics) v2.0.1/v3.0.2 and analyzed using R v.3.5.1, and the following packages and versions in R for analysis: Seurat v3.1.1, ENHANCE v1.0.0, DropletUtils v1.8, clustree v0.4.1, and cluster v2.1.0 two-dimensional gene expression maps, were generated using coordinates from the UMAP algorithm using the R package uwot v0.1.3 implementation. Figures were produced using the following packages and versions in R: RColorBrewer v1.1-2, pheatmap v1.0.12, ggplot v3.2.0, and ggsignif v0.6.0.

scRNA/TCR-seq data processing

Single-cell expression was analyzed using the Cell Ranger Single Cell Software Suite (v2.0.1 and v3.0.2 for P58, P60, and P61, 10X Genomics) to perform quality control, sample demultiplexing, barcode processing, and single-cell 3′ and 5′ gene counting. Sequencing reads were aligned to the GRCh38 human reference genome (Ensembl 84).

Pathway enrichment

Pathways enrichment tests were performed using Metascape (using https://metascape.org/ as of February 2020) (39) with default parameters using differentially expressed genes between early precursors (CD8-KLF2, CD8-GZMK, and CD8-XCL1) and terminally differentiated CD8 + T cells (CD8-GZMH and CD8-LAYN). Gene sets from MSigDB (v.7.0) were downloaded in GMT format from www.gsea-msigdb.org/gsea/msigdb/collections.jsp. These gene sets were used as modules for the AddModuleScore function in Seurat.

Data reprocessing

Data from Guo et al. (13) were downloaded from GSE99254. Data from Sade-Feldman et al. (17) were downloaded from GSE120575. Data from Azizi et al. (16) were downloaded from GSE114727 and reprocessed using Seurat v3 pipeline using the same preprocessing as described previously.

Trajectory analysis

To compute pseudotime alignment of our transcriptomes, we first used Monocle3 (v2.99.1) using the first 30 PCs of the integrated ARTICLE Human in vivo-generated monocyte-derived dendritic cells and macrophages cross-present antigens through a vacuolar pathway Tsing-Lee Tang-Huau 1,2 , Paul Gueguen 1 , Christel Goudot 1 , Mélanie Durand 1 , Mylène Bohec 3 , Sylvain Baulande 3 , Benoit Pasquier 2 , Sebastian Amigorena 1 & Elodie Segura 1 Presentation of exogenous antigens on MHC-I molecules, termed cross-presentation, is essential for cytotoxic CD8 + T cell responses. In mice, dendritic cells (DCs) that arise from monocytes (mo-DCs) during inflammation have a key function in these responses by crosspresenting antigens locally in peripheral tissues. Whether human naturally-occurring mo-DCs can cross-present is unknown. Here, we use human mo-DCs and macrophages directly purified from ascites to address this question. Single-cell RNA-seq data show that ascites CD1c + DCs contain exclusively monocyte-derived cells. Both ascites mo-DCs and monocytederived macrophages cross-present efficiently, but are inefficient for transferring exogenous proteins into their cytosol. Inhibition of cysteine proteases, but not of proteasome, abolishes cross-presentation in these cells. We conclude that human monocyte-derived cells crosspresent exclusively using a vacuolar pathway. Finally, only ascites mo-DCs provide costimulatory signals to induce effector cytotoxic CD8 + T cells. Our findings thus provide important insights on how to harness cross-presentation for therapeutic purposes. C ross-presentation is essential for the induction of cytotoxic CD8 + T cells and efficient immune responses against infections or cancer 1 . Numerous studies in mice have shown that cross-presentation is performed by dendritic cells (DCs). DCs can be classified into four subsets based on ontogeny 2 . "Classical" Batf3-dependent DC1 (cDC1), "classical" Batf3-independent DC2 (cDC2), and plasmacytoid DCs (pDCs) derive from pre-committed bone marrow precursors. Monocytederived DCs (mo-DCs) arise from monocytes recruited into tissues and become the most abundant DC population during inflammation. In mice, cross-presentation is mainly performed by cDC1 in lymphoid organs 1,3 , but mo-DCs have the unique ability to cross-present antigens to CD8 + T cells directly in peripheral tissues 4-6 . Cross-presentation by mo-DCs has a crucial role in the rapid activation of tissue-resident memory CD8 + T cells upon infection 4 and in the efficacy of anti-tumoral treatments based on immunostimulatory agents or chemotherapy 5,7 . Harnessing the cross-presentation capacity of mo-DCs for therapeutic intervention is therefore an attractive prospect. However, determining whether human mo-DCs that arise in tissues can cross-present, and the molecular mechanisms involved, will be a prerequisite.

We and others have shown that the functional specialization for cross-presentation is not conserved between mouse and human DC subsets. In contrast to mouse DCs, human cDC1, cDC2, and pDCs all have a similar ability to cross-present antigens 8-11 . Human mo-DCs generated in vitro from monocytes cultured with GM-CSF and IL-4 can cross-present, and have long been used as a model to understand the biology of cross-presentation, however this culture system gives rise to DCs that do not closely resemble naturally-occurring mo-DCs found in vivo in inflammatory fluids 12 . Therefore, the cross-presentation ability of human mo-DCs remains unclear.

Here, we address this question using human in vivo-generated mo-DCs, directly isolated from peritoneal ascites from cancer patients 12,13 .W efind that mo-DCs and monocyte-derived macrophages (mo-Mac) can both cross-present efficiently, using exclusively a vacuolar pathway. However, only mo-DCs are able to produce co-stimulatory signals for the induction of effector cytotoxic CD8 + T cells.

Results

Tumor ascites CD1c + DCs are monocyte-derived cells. Based on phenotype and gene expression analysis, we have identified the CD1c + DC population found in tumor ascites as naturallyoccurring mo-DCs 12,13 . Because of the sensitivity of the functional assay for cross-presentation, a minor population of cDC within ascites DCs could bias our results. Therefore, we first sought to address the heterogeneity of ascites DCs using singlecell RNA-seq analysis. We purified ascites DCs (gated as HLA-DR + CD11c + CD1c + CD16 -), ascites macrophages (gated as HLA-DR + CD11c + CD1c -CD16 + ) and, for comparison, tonsil cDCs (gated as HLA-DR + CD11c + CD14 -), and analyzed singlecell transcriptomes using a droplet-based method enabling 3′ mRNA counting 14 . To increase the power of the analysis, we combined this dataset with that of blood CD14 + monocytes that we had previously generated 12 . To evaluate the heterogeneity of these population, we performed unsupervised clustering using a graph-based approach with the Seurat package 15 . For visualization of the cell clusters, we used t-distributed stochastic neighbor embedding (t-SNE). Unsupervised clustering of the combined dataset identified 13 main clusters (Fig. 1a,b), and three minor clusters of contaminating cells that were removed from subsequent analysis (see Methods for details). Cluster 1 contained cells from monocytes; clusters 2-3: cells from the macrophage sample; cluster 4: cells from both ascites DCs and macrophages samples; clusters 5-7: cells from the ascites DCs sample; clusters 8-12: cells from the tonsil DCs sample; and cluster 13: cells from the ascites DCs, ascites macrophages, and tonsil DCs samples (Fig. 1a,b).

We then analyzed differentially expressed genes between clusters (Fig. 1c and Supplementary Fig. 1). Cluster 1 displayed high expression of monocyte genes such as CTSS, FCN1, S100A9, LYZ, VCAN. Clusters 2 and 3 shared high expression of macrophage genes such as LGMN, CTSB, CD14, APOE, C1QB, MARCO, CD163, FCGR3A. Cluster 4 expressed high levels of monocyte and macrophage-related genes such as FCN1, S100A9, VCAN, S100A8, MAFB. Clusters 5 and 6 expressed monocyterelated genes such as FCN1, S100A9, VCAN, FCGR1A, FCGR1B, as well as DC genes including CD1C, FCER1A, IFITM2, CLEC10A, FCGR2B. Clusters 7 and 8 showed high expression of DC activation genes CCR7, LAMP3, CCL19, MARCKSL1, CD83, IDO1. Cluster 9 expressed cDC2 hallmark genes CD1C, FCER1A, CLEC10A. Cluster 10 displayed high expression of macrophage genes LGMN, CTSB, CD14, APOE, MARCO, CD163, FCGR3A, S100A8, MAFB. Cluster 11 showed high expression of cDC1 genes such as RAB32, CLEC9A, IRF8, C1orf54, IDO1. Cluster 12 had high expression of genes expressed on progenitors or related to cDC development such as LTB, PRDM16, LST1, RUNX3, CD164. Finally, cluster 13 showed high levels of cell cycle genes including PCLAF, STMN1, MKI67, TOP2A, CDK1. Of note, cluster 13 contained cells from three different samples (ascites DCs, ascites macrophages, and tonsil DCs), showing that in this analysis, cells with similar transcriptional programs are grouped in the same cluster independently of their sample origin.

To confirm the identity of these clusters, we analyzed signature scores in individual cells for several sets of gene signatures (Fig. 2a). For each cell, we calculated the average expression of each signature, substracted by the aggregated expression of control gene sets 16 . We used published gene signatures for blood cDC1 17 , blood cDC2 17 , CD14 + monocytes 17 , and skin CD14 + cells 17 (Supplementary Data 1). We also designed signatures for tissue cDC2 by combining transcriptomic data from blood and spleen cDC2 18 , for in vitro-generated mo-Mac and mo-DCs (obtained with M-CSF, IL-4, and TNFa 12 ), for genes enriched in blood cDC2 compared to ascites DCs and ascites macrophages 13 , and for "activated DC" by selecting genes enriched in both blood cDC2 and in vitro-generated mo-DCs (obtained with GM-CSF and IL-4) exposed to the same stimulus (Menomune microbial vaccine) 19 . Complete lists of genes and strategy for each signature are shown in Supplementary Data 1. Only clusters 2 and 3 expressed the in vitro mo-Mac signature, confirming the identity of these cells as macrophages. As expected from its cellular origin, cluster 1 had the highest score for the CD14 + monocyte signature, but clusters 2, 3, 4, and 5 also displayed high scores for this signature. In addition, clusters 1-6 had high scores for the signature of skin CD14 + cells, which have been shown to derive from monocytes 20 . This is consistent with the notion that ascites DCs and macrophages are related to monocytes. Cluster 11 had the highest score for the blood cDC1 gene signature, confirming the results from differential gene expression. Clusters 6, 7, and 9 had the highest scores for the blood cDC2 signature. However, when using the signature for tissue cDC2, we found that clusters 2, 3, and 9 had the highest scores, with some cells from clusters 5 and 6 also displaying high scores. This suggests that markers for cDC2 may be less robust than for other cell types, possibly due to similar transcriptional programs between cDC2 and other antigen-presenting cells (APC). Consistent with this, when analyzing expression scores for genes enriched in cDC2 compared to ascites DCs and macrophages, we found that clusters 11 Cluster 7 (from ascites DCs) and cluster 8 (from tonsil DCs) shared a high number of markers genes (Fig. 1c) and displayed mixed gene signatures (Fig. 2a). Nevertheless, these cells are not grouped in the same cluster by the clustering algorithm, independently of the resolution used (Supplementary Fig. 2), suggesting that there are significant differences in their transcriptional profile driving their identification as distinct population. Similar profiles could be explained by the convergence of transcriptional programs of mo-DCs and cDC2, in particular for maturation genes, as observed for activated mouse DCs of distinct ontogeny 21 . Alternatively, the separation into distinct clusters could be driven by differences due to tissue origin (fluid versus lymphoid organ). To directly address whether ascites DCs contain a population of cDCs in cluster 7 or whether tonsil DCs contain a population of mo-DCs in cluster 8, we analyzed genes that are the most differentially expressed between cluster 7 (containing ascites DCs) and cluster 8 (containing tonsil DCs) (Supplementary Fig. 3A). Cluster 7 had higher expression for genes reported to be highly expressed in CD14 + monocytes, such as TYROBP (encoding DAP12) 22 , TNFSF13B (encoding BAFF) 23 , NMT1 (a gene essential for monocyte development) 24 , or genes upregulated when monocytes differentiate into DCs such as CST7 25 and CD1E 12 (Supplementary Fig. 3B). By contrast, cluster 8 had higher expression of genes preferentially detected in other clusters of tonsil cDCs such as RELB, FAM60A, IER2, TNFAIP2, SPI1, PTP4A2 (Supplementary Fig. 3C). These genes were found in an independent study to be expressed at similar levels in circulating cDCs from blood and resident cDCs from spleen (by both cDC1 and cDC2) (Supplementary Fig. 3D) 18 , indicating that their differential expression between clusters 7 and 8 is more likely related to distinct ontogeny rather than tissue type. This analysis suggests that cluster 7 corresponds to mo-DCs rather than cDCs.

Based on these results, we annotated cluster 1 as monocytes, clusters 2 and 3 as mo-Mac, cluster 4 as monocyte-derived cells at an early stage of differentiation, clusters 5 and 6 as mo-DCs, cluster 7 as end-stage mo-DCs, cluster 8 as activated cDC2, cluster 9 as cDC2, cluster 10 as contaminating tonsil macrophages, cluster 11 as cDC1, and clusters 12 and 13 as precursor cells (Fig. 2b).

Collectively, these results show that ascites CD1c + DCs do not contain a population of cDCs and support their identification as in vivo-generated mo-DCs.

Human mo-DCs and mo-Mac can both efficiently crosspresent. To address whether ascites mo-DCs can cross-present, we analyzed cross-presentation of a model antigen using a MelanA-specific CD8 + T cell clone (HLA-A2-restricted). Ascites mo-DCs, and mo-Mac for comparison, were incubated with a 34-aa long peptide (requiring processing for cross-presentation) or a pre-processed short peptide corresponding to the minimal epitope, as control for T cell activation ability (Fig. 3a). Ascites mo-DCs and mo-Mac could both cross-present the MelanA antigen, with ascites mo-Mac being more efficient than mo-DCs. mo-Mac were also more efficient for presentation of the short peptide, suggesting a better ability for T cell activation, possibly due to greater MHC class I molecules expression. We also compared the relative expression of genes involved in antigen processing and presentation using Gene Set Enrichment Analysis. Consistent with results from the cross-presentation assay, we did not find any enrichment for gene signatures of antigen presentation between ascites mo-DCs and mo-Mac, although both cell types were enriched for these signatures compared to blood monocytes (Supplementary Fig. 4A-B).

This finding was surprising because we have previously shown, using the same model antigen, that tonsil macrophages are poor cross-presenting cells 9 . To assess whether the ability to crosspresent was induced in macrophages by the ascites environment, we performed the same experiment using in vitro equivalents of these cells, obtained by culturing monocytes with M-CSF, IL-4, and TNFa 12 . In vitro-derived mo-DCs and mo-Mac could both

In vitro mo-Mac cross-present MelanA antigen (Fig. 3b). mo-Mac were again more efficient for presentation of the short peptide, due to higher expression of the MHC class I molecule HLA-A2 (Supplementary Fig. 4C). We confirmed this result using MelanA-coated beads as a model for particulate antigen (Fig. 3c). As an internal control for this assay, we used DCs obtained by culturing blood CD34 + precursors with GM-CSF, Flt3-L, and TNFa. CD1a + DCs could cross-present efficiently, in contrast to CD14 + DCs (Fig. 3d), as previously reported 9,26 . We conclude that mo-DCs and mo-Mac both have the intrinsic ability to cross-present antigens.

mo-DCs and mo-Mac are poor for endosome-to-cytosol transfer. Two main pathways have been described for crosspresentation 1,3 . In the "cytosolic pathway", exogenous antigens are transferred from endocytic compartments into the cytosol, where they are degraded by the proteasome. In the "vacuolar pathway", internalized antigens are degraded in endocytic compartments by lysosomal proteases. To address the intracellular pathway used by mo-DCs and mo-Mac for cross-presentation, we first analyzed their ability to transfer exogenous β-lactamase into their cytosol by measuring the cleavage of a cytosolic β-lactamase-sensitive FRET probe 9,27 . Ascites mo-DCs and mo-Mac were both poor at transferring exogenous β-lactamase into their cytosol (Fig. 4a,b). This was not due to the tumor ascites micro-environment, as the same was found for in vitrodifferentiated mo-DCs and mo-Mac (Fig. 4a,b), nor to deficient uptake of β-lactamase, as all population could efficiently internalize fluorescent β-lactamase (Fig. 4c). By contrast, CD1a + DCs could transfer exogenous β-lactamase into their cytosol more efficiently than CD14 + DCs, as previously reported 9 (Fig. 4a,b). These results suggest that human mo-DCs and mo-Mac do not use the cytosolic pathway for cross-presentation.

mo-DCs and mo-Mac use the vacuolar pathway. To confirm this finding, we analyzed cross-presentation by mo-DCs and mo-Mac in the presence of a proteasome inhibitor, lactacystin (Fig. 5a,b). Cross-presentation by mo-DCs or mo-Mac was not impaired in the presence of lactacystin. By contrast, lactacystin inhibited cross-presentation by CD1a + DCs, as previously reported 9 (Fig. 5c). To confirm that proteasome activity was inhibited by lactacystin in mo-DCs and mo-Mac at the concentration used in the cross-presentation assay, we performed a fluorometric assay for the chymotrypsin-like activity of the proteasome (Supplementary Fig. 5). Lactacystin significantly inhibited proteasome activity in both mo-DCs and mo-Mac. Collectively, these results show that proteasome activity was dispensable for crosspresentation by mo-DCs and mo-Mac.

To directly assess the role of the vacuolar pathway, we used a pan-cathepsin inhibitor to block the activity of lysosomal cysteine proteases. Transcriptomic analysis showed that mo-Mac express overall higher levels of lysosomal proteases than mo-DCs (Supplementary Fig. 6A). In the presence of the cathepsin inhibitor, cross-presentation by mo-DCs and mo-Mac was impaired compared to vehicle control (Fig. 5d,e). This was not due to toxicity of the inhibitor as cell viability was similar in all conditions (Supplementary Fig. 6B). In addition, cross-presentation by CD1a + DCs was not affected by the cathepsin inhibitor (Fig. 5f). These results show that, in mo-DCs and mo-Mac, antigens are degraded by lysosomal proteases for cross-presentation.

We conclude that human monocyte-derived cells use exclusively the vacuolar pathway for cross-presentation.

Only mo-DCs are efficient inducers of cytotoxic CD8 + T cells. To address the outcome of cross-presentation, we analyzed the ability of ascites mo-DCs and mo-Mac to induce the differentiation of cytotoxic effectors from naïve CD8 + T cells. For this assay, we turned to an allogeneic culture system. We co-cultured purified mo-DCs or mo-Mac with allogeneic naïve CD8 + T cells, and assessed T cell proliferation and expression of effector molecules (Granzyme A, Perforin, and IFN-γ). Help from CD4 + T cells is necessary for the differentiation of effector cytotoxic CD8 + T cells 28-30 . In the setting we used, CD4 + T cells have been reported to be essential for CD8 + T cell proliferation 10 , which we confirmed (Fig. 6a,b). When cultured with both naïve CD4 + and CD8 + T cells, only mo-DCs could induce significant proliferation of CD8 + T cells and expression of Granzyme A, Perforin, and IFN-γ (Fig. 6c,d). To evaluate the efficiency of effector differentiation induced by mo-DCs, we performed the same experiment with tonsil DC subsets. We purified cDC1, cDC2, pDCs, and tonsil macrophages 9 . cDC1 and cDC2 were the most efficient inducers of CD8 + T cell proliferation and effector differentiation, while macrophages and pDCs were poor at it (Supplementary Fig. 7). Overall, proliferation and induction of effector molecules were comparable between cDC1, cDC2, and ascites mo-DCs, suggesting that mo-DCs are indeed efficient activators of cytotoxic CD8 + T cells.

Finally, to address the mechanisms underlying the superior ability of ascites mo-DCs to prime effector CD8 + T cells, we compared the capacity of ascites mo-DCs and mo-Mac to provide co-stimulatory signals. Transcriptome analysis showed that ascites mo-DCs express higher levels of co-stimulatory m o l e c u l e st h a na s c i t e sm o -M a c( F i g .7a). Consistent with this, in our allogeneic culture system, mo-DCs were better stimulators of naïve CD4 + T cell proliferation (Fig. 7b), potentially inducing more CD4 + T cell help. To address the ability of ascites mo-DCs and mo-Mac to secrete cytokines involved in the acquisition of CD8 + T cell effector functions 31 , we measured the production of IL-12p70 after ex vivo restimulation. Only ascites mo-DCs were able to secrete IL-12p70 (Fig. 7c

Discussion

Here, we show that human mo-DCs and mo-Mac, both naturally occurring in vivo in peritoneal ascites and generated in vitro from monocytes cultured with M-CSF, IL-4, and TNFα, cross-present exclusively using a vacuolar pathway. However, only ascites mo-DCs induce the differentiation of cytotoxic CD8 + T cells. Numerous studies have shown that DCs are the most efficient cross-presenting cells 1,3 . However, mouse macrophages can cross-present in vitro 32-34 , and in vivo in some settings 35,36 . Human in vitro-generated macrophages can also cross-present antigens 37,38 , in contrast to macrophages isolated from tonsils 9 . Here, we found that human macrophages from tumor ascites can cross-present as efficiently as mo-DCs from the same samples. These discrepancies could be explained by functional differences related to the adaptation of macrophages to their tissue environment 39,40 .

Contradictory data exists on the pathway used for crosspresentation by in vitro-differentiated DCs derived from monocytes using GM-CSF and IL-4, with some studies showing for the cross-presentation of soluble antigens a vacuolar pathway 41,42 and others a cytosolic pathway 38,43-45 , while cross-presentation of cell-associated antigen was reported to be proteasomedependent 46 . The reasons for these discrepancies are not clear. Using a culture model that yields mo-DCs closely resembling in vivo-generated mo-DCs 12 , we found that mo-DCs were inefficient for transferring exogenous proteins into their cytosol, and use a vacuolar pathway for cross-presentation.

It has been proposed that the choice of cytosolic versus vacuolar pathway for cross-presentation is dictated by the form of antigen. Consistent with this, human pDCs cross-present soluble and cell-associated antigens using a cytosolic pathway 9,47 , but cross-presentation of viral antigens has been shown to be proteasome-independent and to use a vacuolar pathway 48 . Furthermore, human cDC1 cross-present soluble antigens using a cytosolic pathway 9,42,49 , but cross-presentation of immune complexes is inhibited by both proteasome and lysosomal proteases inhibitors 50 . The possibility to use either pathway for crosspresentation may also be subset-specific. In support of this, the molecular requirements for cross-presentation are different between mouse cDC1 and mo-DCs both in vivo and in vitro 51,52 .

Our results extend these observations to human DC subsets. Mouse cDC1 are recognized as the main cross-presenting cells 1,3 . mo-DCs can also cross-present efficiently in various inflammatory settings 4-7,51,53,54 . Accumulating evidence indicate that mo-DCs play a key role in the induction and regulation of cytotoxic T cell responses, complementary to that of cDC1. In contrast to cDC1 that interact with CD8 + T cells in lymphoid organs, mo-DCs are able to cross-present antigens and to stimulate effector CD8 + T cells directly in situ, in inflammed tissues 4-6 . In addition, mouse mo-DCs express high levels of co-stimulatory signals involved in the differentiation of cytotoxic CD8 + T cells, including CD70 54,55 .

Our results suggest that human in vivo-generated mo-DCs are equipped for playing a similar role.

Enhancing cross-presentation represents a way of improving vaccination efficiency or anti-tumor immune responses. Deciphering cross-presentation in human DCs that are present in vivo is a pre-requisite for its manipulating for therapeutic purposes. By providing a better understanding of cross-presentation mechanisms in human DC subsets, our results should have important implications for the design of DC-targeted therapies.

Methods

Human samples. Buffy coats from healthy donors (both male and female donors) were obtained from Etablissement Français du Sang (Paris, France) in accordance Cell isolation. Tonsil samples were digested as described previously 56 . In brief, samples were cut into small fragments, digested with 0.1 mg mL -1 Liberase TL (Roche) in the presence of 0.1 mg mL -1 DNAse (Roche) for 40 min at room temperature before addition of 10 mM EDTA. Cells were filtered on a 40 μm cell strainer (BD Falcon) and washed. Light density cells were isolated by centrifugation on a Ficoll gradient (Lymphoprep, Greiner Bio-One). DCs were enriched by depletion of cells expressing CD3, CD15, CD19, CD56, and CD235a using antibody-coated magnetic beads (Miltenyi). Cell subsets were further isolated by cell sorting on a FACSAria instrument after staining for CD11c, HLA-DR, CD14, CD304, CD1c, and CD141 (BD Biosciences). Peripheral blood mononuclear cells (PBMC) were prepared by centrifugation on a Ficoll gradient. Blood CD14 + monocytes were isolated from healthy donors' PBMC by positive selection using anti-CD14-coated magnetic beads according to manufacturer's instructions (Miltenyi). DCs and macrophage population from ascites were isolated after centrifugation of total ascites cells on a Ficoll gradient, enrichment by depletion of cells expressing CD3, CD15, CD19, CD56, and CD235a using antibody-coated magnetic beads (Miltenyi), and cell sorting on a FACSAria instrument. Ascites DCs were gated as HLA-DR + CD11c + CD1c + CD16 -and ascites macrophages as HLA-DR + CD11c + CD1c -CD16 + .

Cell culture. Blood CD34 + cells were isolated from PBMC by positive selection using anti-CD34-coated magnetic beads and magnetic columns according to manufacturer's instructions (Miltenyi). CD34 + cells were cultured for 9-10 days in Yssel medium supplemented with 10% fetal calf serum (FCS), penicillin/streptomycin, 50 ng mL -1 GM-CSF (Miltenyi), 100 ng mL -1 Flt3-L (Miltenyi), and 10 ng mL -1 TNF-α (Miltenyi). DC subsets were isolated by cell sorting on a FACSAria instrument (BD Biosciences) after staining for CD1a and CD14. Monocytes (1×10 6 cells mL -1 ) were cultured for 5 days in RPMI-Glutamax medium (Gibco) supplemented with antibiotics (penicillin and streptomycin) and 10% FCS in the presence or absence of 100 ng mL -1 M-CSF (Miltenyi), 40 ng mL -1 IL-4 (Miltenyi), and 5 ng mL -1 TNF-α (Miltenyi). Cell population was isolated by cell sorting on a FACSAria instrument after staining for CD1a and CD16.

Flow cytometry. Non-specific binding was blocked using TruStain (Biolegend). Unless otherwise stated, cell viability was assessed using DAPI (Sigma Cytosolic translocation assay. Cells (10 × 10 6 cells/condition for ascites cells and 5×10 6 cells/condition for in vitro-generated cultures) were incubated with 0.5 μgmL -1 CCF4-AM (Invitrogen) for 30 min at room temperature at 10 × 10 6 cells mL -1 in loading buffer (120 mM NaCl, 7 mM KCl, 1.8 mM CaCl 2 , 0.8 mM MgCl 2 , 5 mM glucose, 25 mM Hepes, pH 7.3) containing solution B (dilution 1/20, LiveBLAzer FRET-B/G loading kit, Invitrogen) and 1 mM probenecid (Invitrogen). After washing, cells were incubated at 10 × 10 6 cells mL -1 in loading buffer containing 1 mM probenecid in the presence or absence of 2 mg mL -1 of β-lactamase (Penicillinase from Bacillus cereus, Sigma) for 3 h at 4 or 37 °C. Cell viability was assessed using Fixable Viability Dye eFluor780 (eBioscience). After washing, cells were stained for surface markers (CD11c-PeCy7 and CD1c-APC for ascites cells; CD1a-APC-Vio770 and CD16-APC for in vitro-generated monocyte-derived cells; CD1a-APC and CD14-APC-Vio770 for CD34 + cells-derived cultures). Cells were analyzed on a FACSVerse Instrument (BD Biosciences). CCF4-AM fluorescence was assessed by measuring the 450 and 520 nm channels.

Internalization assay. β-Lactamase was conjugated to Atto dye 633 according to manufacturer's instructions (Sigma). Cells were incubated with fluorescent β- Cross-presentation assay. Antigen sources were MelanA short peptide (EAA-GIGILTV), MelanA long peptide (KGHGHSYTTAEEAAGIGILTVILGVL), or beads coated with 750 μM of MelanA long peptide. In brief, Polybead 3-micron polystyrene microspheres (Polypeptide) were washed 3 times in PBS, then incubated at 4 °C overnight in PBS containing 750 μM of MelanA long peptide (in a volume 4 times that of the initial volume of the beads). Beads were washed 3 times in PBS and resuspended in PBS to their initial volume. Cross-presentation assay was performed as described 57 . In brief, purified HLA-A2 + APCs were incubated (10 4 cells per well) for 3-4 h in Yssel medium in V-bottom 96-well plates (Corning) with different concentrations of MelanA long peptide, MelanA short peptide, or MelanA-coated beads, in the presence or absence of 2.5 μgmL -1 lactacystin (clasto-Lactacystin β-Lactone; Merck/Millipore), 10 μM Cathepsin Inhibitor I (inhibiting cathepsin B, cathepsin L, cathepsin S, and papain; Calbiochem) or the corresponding concentration of DMSO (Sigma). After extensive washing, APCs were cultured for 24 h with CD8 T cell LT12 clones 58 (2×10 4 cells per well) in Yssel medium supplemented with 10% FCS. Supernatants were collected and kept at -20 °C until measurement of IFN-γ concentration by ELISA (BD Biosciences). Background levels (APC cultured with LT12 cells without peptide) was substracted for each cell type.

Proteasome activity assay. Proteasome activity was assessed using a proteasome 20S activity kit (Sigma). For proteasome activity assay, mo-DCs and mo-Mac from in vitro cultures were purified using anti-CD16-or anti-CD1a-coated magnetic beads and magnetic columns according to manufacturer's instructions (Miltenyi). Cells (10 5 cells per well, triplicate wells) were incubated with or without 2.5 μgmL -1 lactacystin for 30 min at 37 °C in RPMI supplemented with 10% FCS. Cells were then incubated with the assay loading solution for 3 h at 37 °C. Fluorescence ratio between 490 nm (excitation) and 525 nm (emission) was measured using a FLUOstar Omega instrument (BMG Labtech). Background fluorescence of blank wells (medium without cells) was substracted.

Cytotoxic T lymphocyte polarization assay. Naive CD8 + T cells and CD4 + T cells were isolated from healthy donors' PBMC using EasySep human Naïve CD8 or CD4 isolation kit according to manufacturer's instructions (StemCell Technologies). APC (1×10 4 cells per well) were cultured with naive CD8 + T cells (5×10 4 cells per well) with or without naive CD4 + T cells from the same T cell donor (5×10 4 cells per well) for 7 days in Yssel medium supplemented with 10% FCS. To analyze T cell proliferation, CD8 + T cells were stained with Cell Trace Violet (CTV, Thermo Fisher) prior to culture. To assess the expression of intracellular effector molecules, T cells were stimulated with PMA (50 ng mL -1 ) and ionomycin (1 µg mL -1 ) for 6 h in the presence of BFA (4 µg mL -1 ) for 6 h (all from Sigma). After washing, cells were stained for surface CD4 for 30 min at 4 °C, washed and stained with Live/dead eFluor780 (Thermo Fisher Scientific) for 20 min at 4 °C. Then the cells were fixed and permeabilized (Intracellular Fixation & Permeabilization Buffer Set, eBioscience) and stained for intracellular proteins (Granzyme A, Perforin, and IFN-γ) at room temperature for 45 min in a buffer containing 2% of normal mouse serum. The samples were acquired on a FACSVerse instrument (BD Biosciences).

Cytokine secretion. Sorted cell population (2.5×10 4 cells per well) were incubated during 24 h in Yssel medium in the absence or presence of 1 μgmL -1 dimerized CD40-ligand (Alexis), 1000 IU mL -1 IFNγ (Miltenyi), and 1 µg mL -1 R848 (Invivogen). Supernatants were collected and kept at -20 °C. IL12p70 secretion was assessed by CBA (BD Biosciences).

Gene expression analysis.

For differential gene expression analysis of ascites mo-DCs and mo-Mac and blood monocytes, we used transcriptomic data from Affymetrix microarrays that we had previously generated (GSE40484) 13 . Analysis was performed using R (v3.3.3). Raw data was preprocessed using the Robust Multiarray Average (RMA) method from oligo package 59 . Gene expression levels were analyzed on a base-2 logarithmic scale. Moderated t-tests were performed using the limma package 60 and the p-values were corrected for multiple testing with the Benjamini Hochberg method. Heatmaps were produced using R package gplots.

For gene expression analysis of blood and spleen cDC1 and cDC2, we used published datasets (GSE77671) 18 .

Gene Set Enrichment Analysis. GSEA 61 was performed using the GSEA software (v6.3) 62 and gene signatures from MSigDB (v6.1) 63 . GSEA has been performed with the default parameters except for the number of permutations that we fixed at n = 1000 and the number of min gene at n = 15.

Single-cell RNA-seq library preparation. Cellular suspension (3500 cells) of cellsorted ascites DCs (gated as HLA-DR + CD11c + CD1c + CD16 -), ascites macrophages (gated as HLA-DR + CD11c + CD1c -CD16 + ), or tonsil DCs (gated as HLA-DR + CD11c + CD14 -) was loaded on a 10× Chromium instrument (10× Genomics)

according to manufacturer's protocol based on the 10× GEMCode proprietary technology. Single-cell RNA-Seq libraries were prepared using Chromium Single Cell 3′ v2 Reagent Kit (10× Genomics) according to manufacturer's protocol. Briefly, the initial step consisted in performing an emulsion where individual cells were isolated into droplets together with gel beads coated with unique primers bearing 10× cell barcodes, unique molecular identifiers (UMI), and poly(dT) sequences. Reverse transcription reactions were engaged to generate barcoded fulllength cDNA followed by the disruption of emulsions using the recovery agent and cDNA clean up with DynaBeads MyOne Silane Beads (Thermo Fisher Scientific). Bulk cDNA was amplified using a GeneAmp PCR System 9700 with [START_REF] Lun | EmptyDrops: distinguishing cells from empty droplets in droplet-based single-cell RNA sequencing data[END_REF] Single-cell RNA-seq data analysis. Single-cell expression was analyzed using t h eC e l lR a n g e rS i n g l eC e l lS o f t w a r eS u i te (v2.0.1) to perform quality control, sample de-multiplexing, barcode processing, and single-cell 3′ gene counting 14 . Sequencing reads were aligned to the UCSC hg38 transcriptome using the Cell Ranger suite with default parameters. Samples were merged using Cellranger aggregate function with default parameters. A total of 8404 single cells were analyzed. Mean raw reads per cell were 59,333. Further analysis was performed in R (v3.4) using the Seurat package (v2.2.1) 15 . The gene-cell-barcode matrix of the samples was log-transformed and filtered based on the number of genes detected per cell (any cell with less than 400 genes or more than 5000 genes per cell was filtered out). Any cell with more than 6% of mitochondrial UMI counts and more than 50% of ribosomal UMI was filtered out. Regression in gene expression was performed based on the number of UMI and the percentage of mitochondrial genes. Only genes detected in at least three cells were included. Cells were then scaled to a total of 1 e4 molecules. Altogether, 6964 cells were kept for statistical analysis. To reduce data dimensionality, 5789 variable genes were selected based on their expression and dispersion (expression cut-off = 0, and dispersion cut-off = 0.5). PCA was run on the normalized gene-barcode matrix. Barnes-hut approximation to t-SNE 64 was then performed on the first 19 principal components to visualize cells in a two-dimensional space. The first 19 principal components were used for the t-SNE projection and clustering analysis using the Elbow Plot approach. Clusters were identified using the "Find_Clusters" function in Seurat with a resolution parameter of 0.8. This graph-based clustering method relies on a clustering algorithm based on shared nearest neighbor (SNN) modularity optimization. Unique cluster-specific genes were identified by running the Seurat "Find_All_Markers" function using the MAST framework 65 . Three clusters containing contaminating cells were removed from the analysis: a cluster of 65 cells from the tonsil DCs sample corresponding to NK T cells (top genes: CTSW, KLRB1, CD7, TRDC, XCL2, XCL1, AC092580.4, GNLY, IL2RB, TRBC1, KLRC1, CD3E), a cluster of 58 cells from the ascites DCs sample corresponding to inflammatory CD11c + B lymphocytes (top genes: IGKC,CD79A,JCHAIN,IGLC2,CPNE5,ISG20,CD79B,MZB1,MS4A1,IGHA1, IGHG3, AL928768.3), a cluster of 13 cells from both tonsil DCs and ascites DCs samples corresponding to epithelial cells (top genes: CCDC80, KRT18, TM4SF1, KRT8, CALD1, SLPI, PRG4, NNMT, PLA2G2A, KRT19, DSRN, C3). Heatmaps and violin plots were plotted using Seurat. Data is available at GEO (accession numbers GSE115007 and GSE115006). Scripts used to perform this analysis are available on GitHub (https://github.com/p-gueguen/tang_et_al_2018).

Analysis of gene signatures at the single-cell level. Signature scores were computed using the Seurat function "AddModuleScore" using the gene signature of interest. This function calculates for each individual cell the average expression of each gene signature, subtracted by the aggregated expression of control gene sets 16 . All analyzed genes are binned into 25 bins based on averaged expression, and for each gene of the gene signature, 100 control genes are randomly selected from the same bin as the gene. Featureplots were plotted using minimum and maximum cutoff values for each feature were respectively quantile 3 and quantile 97. We used published gene signatures for skin CD14 + cells,bloodcDC1,bloodcDC2,andbloodCD14 + monocytes 17 . To design genes signature, we used the GeneSign module of BubbleGUM software 66 with our transcriptomic data (GSE40484 and GSE102046) 13,12 . To extract genes enriched in blood cDC2 compared to ascites DCs and ascites macrophages, we used the Mean(Test)/Mean(Ref) method and cut-off of 1% for the adjusted p-value and 1.5 for fold change. To extract genes enriched in in vitro mo-DCs compared to in vitro mo-Mac and monocytes or in vitro mo-Mac compared to in vitro mo-DCs and monocytes, we used the Minimal Pairwise (Mean(Test)/Mean(Ref)) method and cut-off of 1% for the adjusted p-value and 2 for fold change. To design a gene signature for tissue cDC2, we used the GeneSign module of BubbleGUM software 66 with published transcriptomic data (GSE77671) 18 ,using the Minimal Pairwise (Mean(Test)/Mean(Ref)) method and cut-off of 1% for the adjusted p-value and 2 for fold change. To design a signature of genes upregulated upon DC activation, we use published gene expression data of blood cDC2 and in vitro-generated mo-DCs cultured with GM-CSF and IL-4, exposed to the same stimulus (GSE56744 and GSE44721) 19 .Weidentified the genes that are concomitantly (i) up-regulated genes in mo-DCs activated with Menomune (Neisseria meningitidis vaccine, MGL) compared to unstimulated, and (ii) upregulated genes in cDC2 activated with MGL compared to unstimulated. Gene expression levels were analyzed on a base-2 logarithmic scale. Moderated t-tests were performed using the limma package 60 and the p-values were corrected for multiple testing with the Benjamini Hochberg method. We used a cut-off of 1% for the adjusted p-value and 2 for fold change.

Software and statistical analysis. Flow cytometry data were analyzed using FlowJo software v9.9 or v10 (Tree Star). Statistical analyses were performed using the Prism software v7 (GraphPad). Wilcoxon non-parametric test was used. Variance was similar between the groups being compared.

12

NATURE

Titre : Voies de différenciation des lymphocytes T CD8 + infiltrant la tumeur à l'échelle de la cellule unique

Résumé : Il est désormais possible de sonder l'identité cellulaire au niveau de la cellule unique. C'est une excellente opportunité afin d'étudier l'hétérogénéité des cellules T CD8+ afin de mieux comprendre le mécanisme des immunothérapies. Nous avons premièrement examiné l'effet de l'enzyme épigénétique Suv39h1 chez des souris infectées par Listeria monocytogenes (Pace et al. 2018, Science). Nous avons montré que l'expression des gènes mémoires était déréprimée principalement dans les cellules T effectrices Suv39h1-KO. En conséquence, les cellules T Suv39h1 KO CD8+ présentent une survie soutenue et une capacité accrue de reprogrammation de la mémoire à long terme. Sur la base de ces résultats, nous avons examiné l'effet de la déplétion de cette même enzyme dans les lymphocytes infiltrant les tumeurs (LIT) dans un modèle de mélanome de souris, combiné à l'immunothérapie par blocage des points de contrôle immunitaires (Niborski et al. en révision). Le scRNA-seq a révélé qu'en l'absence de Suv39h1, l'anti-PD-1 induit des voies d'activation alternatives pour la survie et la différenciation des cellules effectrices produisant de l'interféron et du GZMB qui n'atteignent pas l'épuisement final, et lie clairement l'émergence de sous-populations de cellules T CD8+ effectrices à l'efficacité du blocage des points de contrôle immunitaires. Enfin, en utilisant le séquençage combiné de l'ARN et du récepteur des cellules T (TCR) au niveau unicellulaire, nous avons étudié l'organisation fonctionnelle des populations des LITs chez des patients atteints de CPNPC primaire [START_REF] Gueguen | Contribution of resident and circulating precursors to tumor-infiltrating CD8 + T cell populations in non-small cell lung cancer patients[END_REF], Science Immunology). Nous avons identifié deux sous-populations de LIT CD8+ exprimant des modules de gènes de type mémoire: l'un est également présent dans le sang (précurseurs circulants), et l'autre dans les tissus juxta-tumoraux (précurseurs résidents dans les tissus). Dans les tumeurs, ces deux populations de précurseurs convergent par un état transitoire unique en cellules différenciées en phase terminale. Ces résultats fournissent un modèle cohérent expliquant l'origine des LITs, leur filiation et leur organisation fonctionnelle dans le CPNPC primaire.

Mots clefs : scRNA-seq, scTCR-seq, Immunologie, Biologie des tumeurs, Biologie computationnelle Title: Tumor-infiltrating CD8+ T cell differentiation paths at the single-cell level Abstract: It is now possible to probe cellular identity at the single-cell level. It is a formidable avenue to study CD8+ T cell heterogeneity, a pressing topic in immunology that has a direct role in better understanding immunotherapies. We started by looking at the effect of the epigenetic enzyme Suv39h1 in Listeria monocytogenes infected mice (Pace et al. 2018, Science). We showed that the expression of stem/memory-related genes was derepressed primarily in Suv39h1-KO effector T cells. As a result, Suv39h1 KO CD8+ T cells exhibit sustained survival and increased long-term memory reprogramming capacity. Building on these results, we looked at the effect of depletion of this enzyme in tumor-infiltrating lymphocytes (TILs) in a mouse melanoma model, combined with immune checkpoint blockade immunotherapy (Niborski et al. in revision). scRNA-seq revealed that in the absence of Suv39h1, anti-PD-1 induces alternative activation pathways for the survival and differentiation of effector cells producing interferon and GZMB that do not reach final exhaustion and links CD8+ effector T cells to the efficacy of blocking immune control points. Using single-cell RNA and T cell receptor (TCR) sequencing, we finally investigated TIL populations' functional organization in mice human primary NSCLC [START_REF] Gueguen | Contribution of resident and circulating precursors to tumor-infiltrating CD8 + T cell populations in non-small cell lung cancer patients[END_REF], Science Immunology). We identified two CD8+ TIL subpopulations expressing memory-like gene modules: one is also present in the blood (circulating precursors), and the other one in juxta-tumor tissue (tissue-resident precursors). In tumors, these two precursor populations converge through a unique transitional state into terminally differentiated cells. These results provide a working model for TIL origin, filiation, and functional organization in primary NSCLC.

Keywords: scRNA-seq, scTCR-seq, Immunology, Cancer Biology, Computational Biology Thank you to my supervisors Sebastian Amigorena and Josh Waterfall for their continuous and complete support, trust, and freedom I had during this thesis. I can say with confidence that you both taught me how to do good science. Thank you Seb for initially welcoming me into the lab, trusting me to do a purely computational Ph.D., and teaching me the inner workings of science. Thank you Josh for your precious guidance that helped me navigate when I was lost in my single-cell datasets.

Thank you to Elodie Segura, Christel Goudot, and Luigia Pace for initially taking me in my M2 internship, training me, and showing me their visions of science, this was fundamental and shaped me for the rest of my thesis.

I also want to thank Thomas Dupic and Aleksandra Walczak, and Thierry Mora at ENS for our discussions, which were really key on my main project. I was humbled to work with such great scientists. I want to especially thank Thomas for giving me an inspiring model on how to do good bioinformatics! I want to thank Myriam for initially taking on the NSCLC project, setting up connections with IMM surgeons, and smooth out the rough edges in these challenging clinical samples' processing.

A massive thank you to Christina for taking on the NSCLC project, for working overtime in the lab many times, for many fruitful discussions, and for successfully going through the publication process with me! Finally, the NSCLC project wouldn't have been possible without patients agreeing on their sample being retrieved and staff and surgeons at IMM being there to take care of patients. I thank them deeply.

Thank you to Eleonora for such a smooth collaboration on myeloids in TNBC. While not in the direct scope of my thesis project, I am certain scientific projects can feed on each other, which was definitely the case here.

Thank you to U932 PIs that I've worked with, giving me different flavors of how to do science: Olivier Lantz, Eliane Piaggio, Emanuela Romano. I truly benefited from Curie's rich work environment.

I want to thank All Curie U932, where lab meetings and lab retreats were great places of scientific growth. I want to thank all Waterfall lab members, namely Julien Vibert and Wilfrid Richer, for always having interesting discussions regarding single-cell methods.

Thank you to Jason Cosgrove and the single-cell group meeting that is always a good place to exchange ideas. I hope to continue exchanges in this group after my defense.

I would also like to thank my advisory committee member Leila Perié for her advice and encouragement.

I want to thank people that make our research possible: Curie IT department and Core sequencing platform, including Sylvain Baulande, Mylène Bohec, Sonia Lameiras. This work wouldn't have been possible without funding. I want to thank la Ligue contre le cancer for funding my Ph.D. and Curie / U932 for their generous COVID extension program and funding me until my defense.

I also want to thank people that make science openly available via preprints, code sharing, tutorials via Twitter, Github, bioRxiv, etc. Keeping myself actively updated in the fast-moving single-cell field was crucial and gave me new ideas to analyze data. This wouldn't have been possible without the rest of the single-cell community.

Finally, I want to thank family and friends for being there, especially at the end, during some tough COVID times.

Table of contents

Acknowledgments

Table of contents

List of appendix

List of figures

Lists of abbreviations

Résumé en Français

Preface

Published papers Papers sent for revisions This thesis has been highly collaborative. In a fast-moving field such as single-cell analyses, the diversity of projects and people made me question methods, find new ways to ask questions and answer them. Below are some other publications that I have contributed to during my candidature but are not presented in this thesis.

Chapter 1. Introduction

Cancer is a dynamic disease. Malignant tumors have diverse phenotypic and molecular profiles at both inter-patient and intra-tumor levels, and cancers grow more heterogeneous over the course of the illness. As a consequence of this heterogeneity, the tumor mass can contain a complex set of cells carrying distinct molecular signatures with differential degrees of treatment sensitivity (cellular heterogeneity), stromal cells, immune cells, tumor cells, among many others. This heterogeneity can result in a non-uniform distribution of genetically distinct cell subpopulations around and within disease sites (spatial heterogeneity) or temporal differences in cells' composition as the tumor progresses, regresses, or evades treatments (temporal heterogeneity). This cellular, spatial and temporal heterogeneity contribute to treatment resistance and metastasis. Therefore, an accurate assessment of tumor heterogeneity, preferably at the single-cell level, is essential for developing effective therapies.

Immunotherapies rely on harnessing the immune system to fight tumors. Recent advances such as CAR-T cell therapy and Immune Checkpoint Blockade (ICB) to fight tumors have proven to be highly successful in certain contexts. However, while having achieved tremendous clinical results, immunotherapy only works for a small fraction of patients. Among the pressing challenges that the field faces, we currently lack good biomarkers of response to immunotherapy, and we are also missing fundamental concepts explaining response to immunotherapies. Among all the immune populations present in the tumor microenvironment, we specifically lack a dynamic view of an important player during immunotherapy: T cells. A precise model of T cell pathways at work during the response to immunotherapy remains to be defined. Specifically, CD8+ T cells, a subpopulation of lymphocytes, play an important role in tumor immunity. Naïve CD8+ T cells are activated by recognition by the T Cell Receptor (TCR) of specific peptides presented by Major Histocompatibility Complex (MHC) class I on antigen-presenting cells (APCs) in peripheral lymphatic organs. After this activation, CD8+ T cells undergo differentiation and a massive expansion to generate large numbers of cytotoxic T cells, effector cells that can migrate into the periphery. The process of CD8+ T-cell differentiation is dictated by antigen strength, co-stimulatory molecules, and cytokines. These external cues induce transcription factors, further specifying CD8+ T-cell decision. More recently, tissue-resident CD8+ T cells have emerged as critical mediators of health and disease. In tumors, the dynamics of infiltrating or resident T cell subsets are largely unknown.

Single-cell genomics is a set of techniques used to probe cellular identity at the single-cell level.

The field has rapidly developed over the last few years, both in terms of throughput and the molecular levels researchers can explore. One of the most commonly used assays, single-cell RNA-seq (scRNA-seq), generates a snapshot of the dynamic gene expression processes within a cell and can identify cell types and explore differences between them. scRNA-seq has been rapidly adopted by researchers in areas where cellular heterogeneity is of particular interest, such as developmental biology, neuroscience, and immunology. scRNA-seq is of particular interest in tumor immunology, where tumor-infiltrating lymphocyte numbers are often too low to be processed using bulk RNA-seq and where cellular heterogeneity and differentiation paths were previously masked by pooling cells together. Coinciding with the rapid development of protocols for producing data, there has been a simultaneous need to analyze and interpret these high dimensional and often noisy datasets and integrating modalities. Indeed, in the last decade, scRNA-seq datasets have grown from a few single cells to millions of cells. These datasets need specialized tools to be accurately stored, processed, and interpreted. My thesis projects' overarching goal is to better understand the different states and differentiation paths of CD8+ T cells, which are one of the main cellular players during immunotherapy. I present here three projects tackling cellular heterogeneity of CD8+ T cells at the single-cell level. First, I analyzed CD8+ T cells from mice infected with Listeria monocytogenes, specifically analyzing the effect of depleting Suv39h1, an epigenetic enzyme implicated in CD8+ T cell differentiation (Pace et al., 2018). In the second project, I looked at the effect of depleting this enzyme in tumor-infiltrating lymphocytes in a mouse melanoma model (B16-OVA), combined with checkpoint blockade immunotherapy (Niborski et al., in revision). Lastly, I focused on human CD8+ T cell differentiation in NSCLC patients. We used a combination of scRNA and scTCR-seq in the tumor, juxta-tumor, and blood to understand ontogeny and recirculation among the different tumor-infiltrating CD8+ subsets [START_REF] Gueguen | Contribution of resident and circulating precursors to tumor-infiltrating CD8 + T cell populations in non-small cell lung cancer patients[END_REF]. Taken together, my thesis covers a wide set of tools and techniques for single-cell RNA and single-cell TCR analysis. I present interpretations of single-cell data and show how tools can be used to extract meaningful biological insights from these big and noisy datasets, with the ultimate goal to improve our understanding of fundamental immunological processes.

Chapter 2. State of the art

Cancer biology

Origin of cancers

Cancer is the second leading cause of mortality worldwide [START_REF] Wang | Global, regional, and national life expectancy, all-cause mortality, and cause-specific mortality for 249 causes of death, 1980-2015: a systematic analysis for the Global Burden of Disease Study 2015[END_REF]. Cancer origin has been attributed to intrinsic DNA replication errors and exogenous exposures, like chemical compounds or smoking, with much debate over the relative contribution of each to cancer risk.

Viruses, bacteria, and radiation rays are other carcinogenesis factors, comprising about 7% of all cancers. Understanding cancer causation is mandatory for effective cancer prevention and treatment.

Cancer occurs by a series of successive mutations in genes so that these mutations change cell functions. In general, cancer disrupts cellular relations and results in the dysfunction of vital genes. This disturbance leads to abnormal proliferation. Historically, this was linked to the discovery of mutations that produce oncogenes with a gain of function and tumor suppressor genes with a loss of function [START_REF] Loeb | Significance of multiple mutations in cancer[END_REF].

One of the classical ways of viewing cancer is to summarize it by a small number of underlying principles [START_REF] Hanahan | The Hallmarks of Cancer[END_REF]. This simplification derives from classical cell biology that postulates that all mammalian cells carry similar molecular machinery regulating their proliferation, differentiation, and death. Six hallmarks of cancer were defined at first, among which: sustaining proliferative signaling, evading growth suppressors, activating invasion and metastasis, enabling replicative immortality, inducing angiogenesis, and resisting cell death.

Progress in the last ten years has added two hallmarks to this list: reprogramming of energy metabolism and avoiding immune destruction. (Fig. 1. Hanahan and Weinberg, 2011), as well as two enabling factors: genome instability, which generates the genetic diversity that expedites their acquisition, and inflammation, which promotes multiple hallmarks. The development of targeted therapy has accompanied these new advances in cancer mechanisms. With the understanding of these ten hallmarks of cancer mechanisms and the main changes in the normal cell growth pathways, the treatment of cancer cells with specific therapies has become possible in many cases. For example, this is the case for NSCLC. The most significant paradigm change in the last 10 years for NSCLC management was EGFR tyrosine kinase inhibitors (TKI) as first-line therapy for patients with a targetable EGFR driver mutation.

However, prolonged treatment with TKIs often results in the development of acquired drug resistance that limits their clinical benefit duration. A potential explanation for acquired TKI resistance may be the acquisition of Receptor Tyrosine Kinases (RTK) mutations, epithelial-mesenchymal transition (EMT), and the activation of alternative signaling pathways, potentially in different tumor clones.

Tumor clonal heterogeneity

Models of carcinogenesis are typically based on the Darwinian principle that evolution requires genetic and/or epigenetic changes that generate new phenotypes. As a tumor grows, it is It is now well established that the tumor immune microenvironment (TIME) plays an important role in tumor clearance and resistance. Currently, three broad classes of TIME can be described according to recent human and mouse data [START_REF] Binnewies | Understanding the tumor immune microenvironment (TIME) for effective therapy[END_REF]. First, TIMEs that are broadly populated with immune cells but are relatively void of cytotoxic lymphocytes (CTLs) in the tumor core are termed infiltrated-excluded (I-E) or "cold" tumors, supposed to be poorly immunogenic. Infiltrated-inflamed TIMEs are considered to be immunologically 'hot' tumors and are characterized by high infiltration of CTLs expressing Programmed cell death 1 (PD-1) and leukocytes and tumor cells expressing the immune-dampening PD-1 ligand PD-L1. Finally, a subclass of TIMEs, here termed TLS-TIMEs, display histological evidence of tertiary lymphoid structures (TLSs), lymphoid aggregates whose cellular composition is similar to that in lymph nodes. TLSs are often correlated with a positive prognosis (Dieu-Nosjean et al., 2014). Definition of the fundamental immunological processes at play during tumor development is still lacking. In such heterogeneous systems, we currently need models and concepts describing the function, localization, circulation, and differentiation of immune cells.

The immune system

Immunity overview

The immune system function is to recognize the self from non-self when the system encounters pathogens and cancer cells. In mammals, the host uses both innate and adaptive mechanisms to detect and eliminate tumor cells. Both of these mechanisms include self-nonself discrimination. These two different kinds of immunity are represented by a growing compendium of cell types, with sometimes overlapping functions.

Innate immunity

The innate immune system consists of cells and proteins that are always present and ready to mobilize and fight microbes at the site of infection. It is an evolutionarily conserved host defense system with key features shared between plants, invertebrates, and mammals [START_REF] Buchmann | Evolution of Innate Immunity: Clues from Invertebrates via Fish to Mammals[END_REF]. Innate immune defenses in mammals encompass virtually all tissues, particularly barrier surfaces such as the skin or the mucosal surfaces of the respiratory and gastrointestinal tract. The innate immune system can sense such situations through germline-encoded receptors (such as toll-like receptors). Innate immune responses can be mediated through cell-dependent mechanisms (e.g., phagocytosis and cytotoxicity) or secreted factors, including antimicrobial peptides, complement factors, alarmins, cytokines/chemokines, and other less-categorized molecules.

Adaptive Immunity

On the other end, the adaptive immune system is called into action against pathogens that can evade or overcome innate immune defenses. Adaptive immunity involves specialized immune cells and antibodies that attack and destroy foreign invaders and can prevent disease in the future by remembering what those substances look like and mounting a new immune response.

Adaptive immunity may last for a few weeks or months or for a long time, sometimes for a person's entire life. Unlike the innate immune system, the adaptive immune system is highly specific to each pathogen the body has encountered and relies on memory to fight infections. In that regard, T and B cells are considered part of the adaptive immunity. An important aspect of the adaptive immune system is to recognize host cells that are infected by viruses, intracellular bacteria, other intracellular parasites, and tumor cells. These cells have evolved an elegant mechanism that recognizes foreign antigens together with self-antigens as a molecular complex [START_REF] Chaplin | Overview of the Immune Response[END_REF]. Fig. 5. Elements of the innate and adaptive immune system. The innate immune system consists of mast cells, macrophages, natural killer (NK) cells, dendritic cells (DCs) and granulocytes. γδ T cells and natural killer T (NKT) cells build the bridge between the adaptive and the innate immune system. Antibody-producing B cells, CD8 + T cells and CD4 + T cells are part of the slow responding but specific adaptive immune system. The figure is taken from Ritter. 2016.

Classical cancer treatments like surgery, chemotherapy, and radiation therapy, have demonstrated limited efficacy for some patients, especially with late-stage disease. Therefore, original and effective cancer treatments are needed for cancer patients. In that regard, harnessing the tumor-infiltrating immune system could represent a great avenue for cancer therapy.

Immunotherapy

Therapeutic manipulation of either the innate or adaptive arms of the immune system, or both, to optimize recognition and elimination of tumors is broadly known as cancer immunotherapy.

The past decade has seen a revolution in cancer treatments by moving away from drugs that target tumors broadly (for example, chemotherapy and radiation) and toward the use of antibody-based immunotherapies or cell-based therapy that modulate immune responses against tumors. a-PD1 checkpoint blockade was the driving force and is now the standard of care for many advanced tumors like melanoma or NSCLC [START_REF] Waldman | A guide to cancer immunotherapy: from T cell basic science to clinical practice[END_REF].

Anticancer immunotherapy can be divided into three general categories: cytokine-targeting therapies (such as interleukins, interferons, and colony-stimulating factors), Immune Checkpoint Inhibitors (ICPI), and cell-based therapies (such as chimeric antigen receptor (CAR) technologies).

Cytokine-targeting therapies

Cytokines are soluble proteins that mediate cell-to-cell communication. Based on the discovery of the potent antitumor activities of several pro-inflammatory cytokines in animal models, clinical research led to the approval of recombinant interferon-alpha and interleukin-2 for the treatment of several malignancies. However, clinical trials' results failed to meet the high expectations from preclinical models and highlighted the limitations of approaches based on unmodified recombinant proteins. These limitations include the short half-life of most cytokines and narrow therapeutic windows with only modest antitumor efficacy, at least as monotherapies. Only two cytokines, IL-2 and IFN-α, demonstrated mild clinical benefit and consequently received The Food and Drug Administration (FDA) approval to treat several malignant diseases. The synergistic combination of anti-CTLA-4 and anti-PD-1 monoclonal antibodies has been approved to treat advanced melanoma, microsatellite instability-high, or mismatch repair-deficient metastatic colorectal cancer by the FDA. It has improved the overall response rate of patients. In this context, cytokines are incorporated into combination clinical trials, mainly in conjunction with anti-PD-1 and anti-PD-L1 monoclonal antibodies.

Antibodies-based

Immune checkpoint blockade

Immune checkpoints are mechanisms of immune cell inhibition that restrains activation. These pathways are crucial for self-tolerance, tolerance to sustained inflammation, reducing tissue-damage, and prevents the immune system from attacking cells indiscriminately. ICPIs are one of the most rapidly growing categories of immune-related agents developed to treat cancer.

They are monoclonal antibodies that block key molecules in immune checkpoint pathways, such as programmed cell death protein 1 (PD-1) and cytotoxic T-lymphocyte-associated protein 4 (CTLA-4). In one seminal clinical report, Pembrolizumab (a humanized monoclonal antibody against PD-1) showed antitumor activity in advanced NSCLC, melanoma, or renal-cell cancer [START_REF] Topalian | Safety, Activity, and Immune Correlates of Anti-PD-1 Antibody in Cancer[END_REF]. This paper, among others [START_REF] Garon | Pembrolizumab for the Treatment of Non-Small-Cell Lung Cancer[END_REF][START_REF] Robert | Pembrolizumab versus Ipilimumab in Advanced Melanoma[END_REF], paved the way for immune checkpoint blockers to be included in the standard of care for many cancer types, and specifically for late-stage and metastatic patients. PD-L1, which is expressed on both tumor cells and immune cells and binds PD-1 to attenuate T-cell activity, is also an immunotherapy target. Several receptor-ligand immune checkpoints have been identified, such as PD-1 interacting with its ligands programmed death-ligand 1,2 (PD-L1,2). Nivolumab and pembrolizumab are highly specific antibodies that bind PD-1, while atezolizumab, avelumab, and durvalumab target PD-L1, and Ipilimumab targets CTLA-4. Treating solid tumors with immune checkpoint inhibitors can be highly effective and can lead to durable responses in approximately around 20-30% of patients with NSCLC, melanoma, and renal-cell cancer [START_REF] Sharma | Immune Checkpoint Targeting in Cancer Therapy: Toward Combination Strategies with Curative Potential[END_REF]. Further studies confirmed the utility of anti-PD-1 antibodies in treating cases of advanced NSCLC, with an objective response rate (ORR) of 19.4% and a median duration of response of 12.5 months [START_REF] Garon | Pembrolizumab for the Treatment of Non-Small-Cell Lung Cancer[END_REF].

Multispecific antibodies

The term bispecific antibody (bsAb) is used to describe a large family of molecules designed to recognize two different epitopes or antigens. The classical bsAb application which is T cell redirection and engagement was first described in the mid-1980 (Staerz and Bevan, 1986) but did not reach patients until 2009 with the European Union approval of catumaxomab for the intraperitoneal treatment of malignant ascites. Recent results showed impressive clinical results of another approved T cell-engaging bsAb (bsTCE), blinatumomab (CD3 × B lymphocyte antigen CD19) [START_REF] Goebeler | Blinatumomab: a CD19/CD3 bispecific T cell engager (BiTE) with unique anti-tumor efficacy[END_REF], provoked renewed interest and investment in this concept.

More recently, results showed that trispecific antibodies could be efficiently used in patients. A trispecific antibody that interacts with CD38, CD3 and CD28 was developed to enhance both T cell activation and tumor targeting (CD3 and CD28 for T cell stimulation and the anti-CD38 domain directs T cells to myeloma cells) (Wu et al., 2020).

TIL-therapy

Adoptive cell therapy with tumor-infiltrating lymphocytes (TIL) could be one of these additional treatment options. In TIL therapy, T cells residing in patient-specific tumor material are isolated and expanded ex vivo in a dedicated production facility and given back to the patient as a single intravenous infusion after a lymphodepleting non-myeloablative preparative regimen and subsequent treatment with interleukin-2 (IL-2). TIL treatment was introduced in small clinical trials in the '80s [START_REF] Rosenberg | Use of tumor-infiltrating lymphocytes and interleukin-2 in the immunotherapy of patients with metastatic melanoma. A preliminary report[END_REF]

Cancer vaccines

In the case of diseases caused by viruses and bacteria, vaccines work by exposing people to a weakened or inactivated version of the threat. In the case of cancer, however, the situation is more complicated for several reasons. In particular, unlike bacteria and viruses, which appear foreign to our immune system, cancer cells more closely resemble our normal, healthy cells. For this reason, the choice of antigen is the single most important component of cancer vaccine design. Ideally, the antigen should be expressed specifically by cancer cells (and not in normal cells), present on all cancer cells, necessary for cancer cell survival (such that the cancer cannot escape immune attack by downregulating the antigen), and highly immunogenic. Few if any antigens meet all of these criteria, yet there are several classes of antigens that have been employed in cancer vaccines. To date, most cancer vaccines have targeted tumor-associated antigens, which are self-proteins that are abnormally expressed by cancer cells. TAAs include: cancer/germline antigens (also known as cancer/testis antigens), which are normally expressed only in immune-privileged germline cells (e.g., MAGE-A1, MAGE-A3, and NY-ESO-1). DC vaccines, in which patient-derived, autologous DCs are either loaded with peptide antigen or transfected with antigen genes, have also been studied extensively. The first US Food and Drug Administration-approved cancer vaccine, sipuleucel-T (Provenge), is in use for metastatic castration-resistant prostate cancer (mCRPC). This vaccine is generated by enriching DCs derived from each patient by leukapheresis and activation ex vivo with a chimeric protein, GM-CSF fused to the antigen PAP. Several other DC vaccines are now being developed and tested (Perez and De Palma, 2019). However, despite this major oncologic advance, all these immunotherapy strategies are not universally successful across patients and tumor types. The response rate is indeed highly variable, ranging from complete to no response. Only a minority of patients achieve long-term survival. Thus, there is a great need to find predictive biomarkers for patient stratification and selection.

Biomarkers of response

Tumor mutational burden Tumor mutational burden (TMB)-the number of somatic mutations per DNA megabase has emerged as a proxy for neoantigen burden that is a biomarker associated with ICPI outcomes [START_REF] Klempner | Tumor Mutational Burden as a Predictive Biomarker for Response to Immune Checkpoint Inhibitors: A Review of Current Evidence[END_REF]. It is believed that tumor cells with high TMB may have more neoantigens. These neoantigens can be recognized by T cells, inciting an anti-tumor response.

TMB can be reliably estimated using algorithms from next-generation sequencing assays that interrogate a sufficiently large subset of the exome. Biological processes contributing to elevated TMB can result from exposure to cigarette smoke and ultraviolet radiation, from deleterious mutations in mismatch repair leading to microsatellite instability, or mutations in the DNA repair machinery. Clinical studies have shown that patients with higher TMB experience longer survival and greater response rates following treatment with ICPIs than those with lower TMB levels.

Expression of PD-L1

Expression of PD-L1 is another biomarker of response. PD-L1 is one of the first and remains the most widely investigated biomarker for predicting response to ICB. PD-L1 is present on various cell types, including cancer cells, dendritic cells, macrophages, and activated T and B lymphocytes. PD-L1 normally plays a role in protecting tissues from excessive inflammation and autoimmune reactions. However, tumors are also able to express PD-L1, leading to their escape from the immune system. Escape from immune response occurs when PD-L1 produced by malignant cells binds to PD-1 on T cells. The binding of PD-L1 to T cells leads to attenuation or inhibition of their activity, allowing tumors to avoid immune surveillance. Based on these findings, pembrolizumab (aPD1) monotherapy was initially approved by the US Food and Drug Administration (FDA) for the first-line treatment of patients with metastatic NSCLC whose tumors contain high concentrations of PD-L1 (≥50% of tumor cell staining) [START_REF] Garon | Pembrolizumab for the Treatment of Non-Small-Cell Lung Cancer[END_REF].

MSI status

Microsatellites (MSs) are short sequences of DNA (usually 1-6 nucleotides long) tandemly repeated throughout the genome. These sequences are located in both gene and intergenic regions. MS instability (MSI) occurs when the genome gains or loses ≥1 repeats. The DNA repair mechanism responsible for correcting these errors is known as the mismatch repair (MMR) system. Germline or somatic mutations in any of these genes or hypermethylation in the promoter result in defective MMR and, thus, an inability to repair errors during DNA replication.

As these errors tend to occur predominantly in MS regions, tumors with such errors have high MSI. As a result of defective MMR, both the number of mutations and predicted number of neoantigens are higher in tumors with this defect. The increased number of neoantigens might render tumors more immunogenic and, thus, more likely to respond to immunotherapy.

The recent surge of immunotherapy strategies and drugs are backed by objective clinical responses. These therapies made manifest that we can harness the immune system to fight cancer. Nevertheless, we still lack precise pathways that are at play during tumor recognition by the immune system. One attractive target to understand immunotherapy is one of the major players during immune response: CD8+ T cells, known for recognizing and killing tumor cells effectively.

CD8+ T cell CD8+ T cell origin

T cells with multiple specificities are generated in the thymus before antigen exposure. After maturation into single-positive (CD4+CD8-or CD4-CD8+) thymocytes, they are then released from the thymus to peripheral tissues. Naive T cells, which have not yet recognized and responded to antigens, circulate throughout the body and reside in peripheral lymphoid organs in a resting state, and acquire powerful functional properties only after activated by their T Cell Receptor (TCR). The goal of T cell activation is to generate, from a small pool of naive cells with predetermined receptors for any antigen, a large number of functional effector cells that can eliminate that antigen (called clonal expansion) and a population of memory cells that remain for long periods to react against the antigen in case it is reintroduced rapidly. (Abbas, Cellular and Molecular Immunology, 2011). The use of CD8+ T cells has been a focus of clinical cancer therapy for over 20 years. T effector (TEFF) cells have long been recognized as important mediators of tumor protection based on their ability to clonally expand and exert cytotoxic function. However, studies in mice and humans have also revealed that the presence of tumor-infiltrating TEFF cells is not sufficient for tumor rejection because of apoptosis induced clonal contraction. Cancer is a chronic disease and, unlike immunity to acute infections, curative immunity to cancer is thought to require long-lived T cell immunity. Indeed, the persistence of adoptively transferred T cells in peripheral blood has been shown to correlate with tumor regression in cancer patients. In contrast to TEFF cells, memory CD8+ T cells are capable of durably persisting and functioning throughout host tissues and tumors, making them a gold standard for anti-tumor immunity.

CD8+ T cell differentiation

After antigen recognition, naive cells progressively differentiate into effector T cells. CD8+ effectors kill the infected or cancer cells, yet, immune responses to an antigen usually also result in the generation of memory T cells specific for that antigen, which persist for many years. When looking at CD8+ T cell differentiation process, we can distinguish multiple differentiation subsets. Naive cells can differentiate into many subsets such as Stem Cell Memory, Central Memory, Effector Memory, and Effectors. However, these populations are still ill-defined.

There exist currently 2 main differentiation models: the linear and the branching model.

According to the linear model of memory T cell differentiation, most effector cells die, and some survivors develop into the memory population. According to the branched differentiation model, some naive cells can directly differentiate into memory cells. These models are usually coming from acute infection models, like Lymphocytic Choriomeningitis Virus (LCMV) infection, and do not necessarily recapitulate well what happens in human tumors.

CD8+ T cell dysfunction

Upon activation driven by antigen recognition, naive CD8 T cells differentiate into memory precursor effector cells (MPECs) and eventually fully differentiate into terminal effector (TEFF) T cells (in acute infection/vaccination) to control peak antigen load or exhausted (TEX) precursor cells (in chronic/cancer). During acute infection or vaccination, most TEFF cells. The remaining cells differentiate into memory T cells (TMEM) that patrol different sites in the body or resident memory T cells (TRM) cells retained in tissues. There are notable differences between acute infection/vaccination and chronic infection/cancer shown in Fig. 10.

Namely, after TCR chronic stimulation, such as in chronic infection or cancer, it has been documented that T cells become hyporesponsive (Wherry and Kurachi, 2015). This is believed to be a physiological adaptation, effectively reducing tissue damage and inflammation, taken advantage of by the tumor cells. T cell exhaustion is the dysfunction or physical elimination of T cells in response to specific antigens during chronic diseases, including viral infection and cancer [START_REF] Mclane | CD8 T Cell Exhaustion During Chronic Viral Infection and Cancer[END_REF]. It is believed that T cells undergoing exhaustion lose their effector functions, and their memory recall is compromised (Wherry and Kurachi, 2015), but there is still debate about which T cells can be defined as exhausted or not [START_REF] Blank | Defining 'T cell exhaustion[END_REF].

Nevertheless, exhausted T cells are characterized by increased expression of co-inhibitory receptors, including PD-1, CTLA-4, LAG3, CD223, TIM-3, CD39, TIGIT, CD244. It is now believed that during chronic infection or cancer, Tex precursors are unable to fully clear antigen and gradually develop into early Tex cells that retain some proliferative potential. In the setting of persisting antigen, Tex progenitor cells differentiate into terminal Tex cells and, in some cases, can be eventually deleted entirely. Early exhausted cells form during chronic infections and cancer, where antigen persists. Early Tex cells are characterized by intermediate expression of PD-1 and low Eomes, with high expression Tcf-1 (also known as TCF7 in humans) in the progenitor population. Early Tex cells give rise to terminally exhausted T cells characterized by high expression of immune checkpoints and loss of Tcf-1, and loss of the ability to proliferate further upon additional antigen stimulation. This phenotype is also believed to have specific transcriptome and epigenetic features distinct from senescence [START_REF] Zhao | Exhaustion and senescence: two crucial dysfunctional states of T cells in the tumor microenvironment[END_REF] and anergy [START_REF] Crespo | T cell anergy, exhaustion, senescence, and stemness in the tumor microenvironment[END_REF]. It was also described that other cell types from the surrounding TME drive dysfunction through soluble and cellular factors, such as Tregs [START_REF] Li | Regulatory T cells in tumor microenvironment: new mechanisms, potential therapeutic strategies and future prospects[END_REF], tumor-associated macrophages [START_REF] Nixon | IRF8 Governs Tumor-Associated Macrophage Control of T Cell Exhaustion[END_REF], and myeloid-derived suppressor cells (MSDCs) [START_REF] Bird | MDSC metabolite stuns T cells[END_REF]. Since the exhausted phenotype is linked with chronic TCR activation, notably in tumor tissue, this phenotype is tightly linked to the tissue-resident memory subset.

CD8+ T cell tissue residency

Pioneering work in the early 2000s discovered memory CD8+ T cells in peripheral tissues, which were initially thought to be TEM in recirculation from the blood [START_REF] Bromley | Chemokine receptor CCR7 guides T cell exit from peripheral tissues and entry into afferent lymphatics[END_REF]. Early findings in models of vesicular stomatitis virus and Listeria monocytogenes infections [START_REF] Sheridan | Oral Infection Drives a Distinct Population of Intestinal Resident Memory CD8+ T Cells with Enhanced Protective Function[END_REF] revealed that these T cells were actually durably resident in tissues, and such resident memory responses have since been documented in response to multiple pathogens [START_REF] Schenkel | Tissue-Resident Memory T Cells[END_REF]. TRM cells are known to express CD69, a marker of T cell activation, which blocks T cell expression of S1PR1 [START_REF] Skon | Transcriptional downregulation of S1pr1 is required for the establishment of resident memory CD8 + T cells[END_REF], and lack the expression of CCR7, which cooperates with S1PR1 for tissue egress, thereby promoting tissue retention and residency. Loss of KLF2, which regulates S1PR1, has also been shown to enhance tissue residency. Taken together, these findings have revealed TRM cells as a distinct memory lineage.

Epigenetics of CD8+ T cell

DNA condensation is required for DNA to fit in the eukaryotic nucleus. This compaction is achieved thanks to structures called nucleosomes. However, nucleosomal DNA can be inherently restrictive to access of transcriptional machinery. Within a nucleosome, part of the DNA faces the histones' globular domains, causing the DNA sequence on that side to be hidden sterically. Recent advances in molecular and computational technologies have enabled the query of higher-order chromatin architecture and accessibility [START_REF] Buenrostro | ATAC-seq: A Method for Assaying Chromatin Accessibility Genome-Wide[END_REF][START_REF] Lieberman-Aiden | Comprehensive mapping of long range interactions reveals folding principles of the human genome[END_REF].

Post-translational modification of histones can further facilitate chromatin compaction. For example, the histone H3 tail can be covalently methylated at lysine 9 (H3K9me2 or H3K9me3) or lysine 27 (H3K27me3), and these modifications are linked to additional chromatin compaction.

Nucleosomes carrying these modifications flank accessible regulatory elements acting as enhancers and promoters. In T lymphocytes, similar to other developmental programs, the hierarchical genome organization is shaped by a highly coordinated labor division among different classes of sequence-specific transcription factors. The first lineage decision task is carried out by a small number of lineage-determining transcription factors often referred to as "master regulators". Although many studies have characterized the transcription factors that control the differentiation of T cells, the corresponding epigenetic states and associated chromatin dynamics involved in the establishment and maintenance of CD8+ T cell memory and effector identities are still incompletely understood. Chapters 3 and 4 of this thesis deal with CD8+ differentiation in regard to Suv39h1, an enzyme responsible for the trimethylation of histone H3K9.

As we've seen, CD8+ T cell functionality, heterogeneity, tissue specialization, and epigenetic commitment is yet to be fully evaluated. As technological progress advances, CD8+ T cell definitions are being refined. All these different facets of CD8+ T cell phenotype (transcriptomics, epigenetics, TCR clonality) are necessary to fully understand the range of CD8+ T cell heterogeneity. Since population-averaged measurements might mask differences among these closely related CD8+ T cell subsets, all these different facets should be studied, if possible, at the single-cell level.

Single-cell technologies

From fruit smoothie to fruit salad

The adult human body comprises ~37 trillion cells [START_REF] Bianconi | An estimation of the number of cells in the human body[END_REF], which are the functional units of organismal systems. The majority of experimental and clinical results are based on the assumption that all of the cells in a culture or tissue are homogeneous. The growing "omics" fields of study (genomics, proteomics, transcriptomics, etc.) are mainly based on the bulk of cells or tissue samples. However, this averaging can miss some critical information because the sample's heterogeneity is ignored, while the nature of biology is diverse. Even samples considered as homogeneous, such after cell-sorting, can show substantial heterogeneity (Wagner et al. 2016), because of cell subtypes, functional variation, or stochastic events. Thus, it is necessary to conduct omics studies at the single-cell level.

Flow cytometry analysis has been used for single-cell profiling for the past several decades, with some limitations. First, it is a targeted analysis method for only a preselected set of molecules.

Second, due to the spectral limitation of fluorescent proteins, this method can profile up to 17 proteins simultaneously, which can be extended to ~40 proteins by mass cytometry.

Single-cell technologies aim to define cell types in terms of distinctive molecular profiles (e.g., gene expression, protein expression, chromatin accessibility). The scRNA-seq approach seeks to determine, in a data-driven way, categories and spectrums of cellular identity, rather than relying on ones imposed by prior conceptions. The remaining challenge is to connect this information with classical cellular descriptions (e.g., location and morphology). A popular analogy for single-cell technologies, as compared to bulk technologies, was coined by Aviv Regev and Alex Shalek. If a different fruit represents each cell, bulk RNA-seq is the equivalent of a fruit smoothie, scRNA-seq can be seen as a fruit salad, up to a fruit tart if we are able to keep the orginal spatial localization. Although droplet-based techniques are currently the most commonly used cell capture technologies, other approaches have been proposed to capture even more cells at a lower cost per cell, like the Seq-Well protocol [START_REF] Gierahn | Seq-Well: portable, low-cost RNA sequencing of single cells at high throughput[END_REF]. Some datasets are generated using combinatorial indexing technologies, such as sci-RNA-seq [START_REF] Cao | Comprehensive single-cell transcriptional profiling of a multicellular organism[END_REF], or in combination with droplet microfluidics, such as scifi-RNA-seq [START_REF] Datlinger | Ultra-high throughput single-cell RNA sequencing by combinatorial fluidic indexing[END_REF], which increases the throughput of the original droplet-based single-cell RNA-seq up to 15-fold. This technology, sometimes called 'split and pool barcoding', relies on multiple barcoding rounds, where cells are split into wells, barcoded, and then pooled to be sequenced together.

Single-cell multi-omics techniques

Many protocols have been developed to enable multiple measurements from the same individual cells. Some approaches that allow multiple measurements from the same individual cells include G&T-seq [START_REF] Macaulay | G&T-seq: parallel sequencing of single-cell genomes and transcriptomes[END_REF], scMT-seq, which measures RNA and DNA methylation [START_REF] Hu | Simultaneous profiling of transcriptome and DNA methylome from a single cell[END_REF], scTrio-seq which can measure the genome, transcriptome and methylome simultaneously [START_REF] Hou | Single-cell triple omics sequencing reveals genetic, epigenetic, and transcriptomic heterogeneity in hepatocellular carcinomas[END_REF] and ECCITE-Seq [START_REF] Mimitou | Multiplexed detection of proteins, transcriptomes, clonotypes and CRISPR perturbations in single cells[END_REF] which allows to profile up to four modalities per single cell (proteins, transcriptomes, clonotypes, and Chosen technology: 10x Genomics Chromium system For my three projects, we used the newly developed 10x Genomics Chromium system. This machine had the advantage of segregating and performing scRNA-seq preparation for up to 10 000 cells.

A droplet-based system 10x Genomics launched their Chromium single-cell technology on February 11th, 2016. This droplet-based system can generate up to thousands of single-cell libraries per sample. This machine relies on the encapsulation of cells in oil droplets and the 3' or 5' capture of the mRNA.

First, cells are dissociated if they come from solid tissue, using either mechanical dissociation or enzymes. Then, cells are isolated using microfluidics and barcoded that allow the pooling of everything together. Reverse Transcriptase (RT) is performed on the mRNA to obtain cDNA. cDNA is then amplified via PCR or IVT. Sequencing adaptors are added to the reads to create the sequencing library. This library is then sequenced, and cell-specific expression profiles are obtained.

Cell and RNA Barcoding

Each cell and RNA molecule is barcoded with a unique randomer (a random sequence of nucleotides).

-Cell barcoding: Cells are tagged using a 14 bp barcode drawn from ∼750,000 designed sequences to index GEMs (Gel bead in EMulsion).

-These sequences will constitute the whitelist in which the alignment algorithm looks into when mapping barcodes and correcting barcode sequencing errors.

-RNA barcoding: UMIs are unique randomers that are used to tag mRNA molecules. UMIs are used to differentiate between PCR duplicates and real transcripts abundance [START_REF] Islam | Quantitative single-cell RNA-seq with unique molecular identifiers[END_REF].Fig. 18). UMIs are used to reduce amplification noise. UMIs are tagged on individual RNA molecules before amplification. Here 10x Genomics uses 10bp randomer to index RNA molecules. In the 10X Genomics (Chromium controller), cells are encapsulated in oil droplets, each containing a gel bead (Fig. 20). mRNA is captured on the barcoded gel beads. mRNA is then retro-transcribed into cDNA and amplified using PCR. Finally, barcoded oil droplets are lysed, and cDNAs are pooled together to be sequenced. scRNA-seq and scTCR-seq Though critical in any understanding of human health, the adaptive immune response has proven challenging to study. The immune response is largely governed by two factors-lymphocyte diversity and antigen specificity. Healthy individuals typically have hugely diverse populations of T cells, each with an antigen-specific T-cell Receptor (TCR) formed from diverse alpha and beta chain genes. This TCR diversity is generated through V(D)J Recombination, the somatic recombination of Variable (V), Diversity (D), and Joining (J) gene sequences. This process creates sequence variability at the V/D and D/J recombination junctions, providing T cells with a layer of genetic diversity that does not exist in other cell types and can be measured only with a system that can profile the exact V, D, and J sequences for both alpha and beta genes within each cell. The V(D)J system combines microfluidics and 5' molecular barcoding to assemble full-length V(D)J sequences on a cell-by-cell basis. Then, the accompanying bioinformatics software, the Cell Ranger pipeline, profiles, and groups the V(D)J sequences to build clonotypes, groups of T cells sharing the same TCRs. Bulk TCR sequencing cannot assign paired TCR subunits to individual T cells; therefore, this information was historically obtained by isolating T cells using in vitro culture methods, then subcloning cells through a low-throughput, labor-intensive process of creating homogeneous populations of T cells that could then be sequenced to determine individual TCR αβ subunit combinations.

Paired αβ sequences can now be obtained at higher throughput using 10X Genomics technology. Combining scRNA-seq with single-cell TCR sequencing (scTCR-seq) links T cell phenotypes with the specific TCR clonotypes of individual T cells. In 10X technology, RT enzyme and poly(dT) primer are delivered to all GEMs as part of the master mix. Barcoded 5' template switch oligo delivered to GEMs from Gel Beads. The RT reaction generates unbiased cDNA with a sequencing adapter, a cell barcode, and a UMI on the 5' end. This is followed by enrichment PCR with a universal primer for the 5' adapter and successive nested primers for the TCR constant regions. Then fragmentation and sequencing optimized for assembly of the full V(D)J sequence (5' UTR to constant regions) from short paired-end reads on a cell-by-cell basis.

By complementing single-cell phenotypic profiling of T cells (such as via scRNA-seq) to identify their respective TCR clonotypes, these integrated approaches add information on T cell antigen specificity to immune cells' analyses and allow to explore of the roles of antigen-specific T cells. To align transcripts to the genome, pseudo alignment methods can be applied to scRNA-seq datasets but have required adaptations such as the alevin method for UMI-based datasets [START_REF] Srivastava | Alevin efficiently estimates accurate gene abundances from dscRNA-seq data[END_REF]in the Salmon package or Kallisto via the Bustools workflow [START_REF] Melsted | Modular and efficient pre-processing of single-cell RNA-seq[END_REF]. When using conventional aligners, UMI samples need extra processing with tools like UMI-tools [START_REF] Smith | UMI-tools: modeling sequencing errors in Unique Molecular Identifiers to improve quantification accuracy[END_REF] or zUMIs [START_REF] Parekh | zUMIs -A fast and flexible pipeline to process RNA sequencing data with UMIs[END_REF] to assign cell barcodes and deduplicate UMIs. In our analyses, sed Cellranger, a complete pre-processing pipeline that includes an automated downstream analysis based on the STAR aligner [START_REF] Dobin | STAR: ultrafast universal RNA-seq aligner[END_REF].

Cellranger processing workflow is based on some well-established tools from bulk RNA-seq, but some other tools are designed specifically for scRNA-seq data [START_REF] Poirion | Single-Cell Transcriptomics Bioinformatics and Computational Challenges[END_REF]. 10x

Genomics develops the CellRanger Workflow (Fig. 22) as a way to demultiplex, align and count transcripts. CellRanger is composed of 2 main functions: mkfastq, count, and an optional function: aggr.

-cellranger mkfastq: performs a barcode-aware demultiplexing from BCL to FASTQ files. This function is built as a wrapper around Illumina's BCL2FASTQ tool.

-cellranger count: performs read-level analysis of a single library, transcriptome alignment with STAR, barcode processing, gene counting, produces gene/cell matrix, produces expression analysis and static visualizations.

-cellranger aggr: aggregate `cellranger count`results. It combines data from multiple libraries and normalizes depth, resulting in a combined gene/cell matrix. This function will remove reads from the original matrix to match the sample with the least sequencing depth. This function can be replaced by more recent approaches designed to integrate samples.

Quality Control

This step is usually similar to RNA-seq. First, adapters are removed. Reads are trimmed according to their quality.

Alignment

Cell Ranger uses an aligner called STAR, which performs splicing-aware alignment of reads to the genome. Cell Ranger then uses the transcript annotation GTF to bucket the reads into exonic, intronic, and intergenic, and determine whether the reads align (confidently) to the genome. A read is exonic if at least 50% of it intersects an exon, intronic if it is non-exonic and intersects an intron, and intergenic otherwise.

UMI Counting

Before counting UMIs, Cell Ranger attempts to correct for sequencing errors in the UMI sequences. Reads that were confidently mapped to the transcriptome are placed into groups that share the same barcode, UMI, and gene annotation. If two groups of reads have the same barcode and gene, but their UMIs differ by a single base (i.e., are Hamming distance 1 apart), then one of the UMIs was likely introduced by a substitution error in sequencing. In this case, the UMI of the less-supported read group is corrected to the UMI with higher support.

Cell Ranger again groups the reads by barcode, UMI (possibly corrected), and gene annotation.

If two or more groups of reads have the same barcode and UMI but different gene annotations, the gene annotation with the most supporting reads is kept for UMI counting, and the other read groups are discarded. In case of a tie for maximal read support, all read groups are discarded, as the gene cannot be confidently assigned.

After these two filtering steps, each observed barcode, UMI, gene combination is recorded as a UMI count in the unfiltered feature-barcode matrix. The number of reads supporting each counted UMI is also recorded in the molecule info file. multinomial model to infer cell-free droplets' transcriptome profile [START_REF] Lun | EmptyDrops: distinguishing cells from empty droplets in droplet-based single-cell RNA sequencing data[END_REF]. By estimating the deviations from this profile, EmptyDrops assigns barcodes with significant deviations as cell-containing droplets. This method is now directly applied in Cellranger (v3 and after). Current multi-omic pipelines like Cellranger v5 can leverage multi-omic datasets (scRNA-seq/scTCR-seq or scRNA-seq/scATAC-seq) to do cell calling on both modalities and increase the number and quality of cell recovered.

After empty droplets have been removed, the second step of filtering is to remove low-quality cells. Cells with less than a defined number of UMIs/genes are discarded. This approach is already done directly by CellRanger but can be done later with more stringent parameters. Specific tools to remove bad quality cells are also available, like Celloline [START_REF] Ilicic | Classification of low quality cells from single-cell RNA-seq data[END_REF]. Automated methods were also developed to perform droplet calling. One example is the application of neural networks to filter barcodes called CellBender [START_REF] Fleming | CellBender remove-background: a deep generative model for unsupervised removal of background noise from scRNA-seq datasets[END_REF] [START_REF] Xi | Benchmarking Computational Doublet-Detection Methods for Single-Cell RNA Sequencing Data[END_REF]. As heterotypic doublets (different cell types/states in the same doublet) can be identified, homotypic doublets are still hard to distinguish.

Gene filtering

Filtering and selection of features is also an important step. Genes or transcripts that are lowly expressed are typically removed from datasets to reduce computational load and to reduce the noise provided by these lowly expressed genes. These features are often selected based on how expressed or variable they are across the dataset. Alternative selection methods have been proposed, such as M3Drop which identifies genes with more zeros than would be expected based on their mean expression (Andrews and Hemberg, 2019a). To remove housekeeping genes, genes with above-average expression and low coefficient of variance can also be removed. It is also possible to directly remove cycling genes from the dataset, as these might confound interesting signals. However, since many genes will co-vary with cycling genes during the cell cycle, other methods to remove the cell-cycle were developed.

Cell cycle removal

The cell cycle often has an enormous impact on cellular gene expression and can, in turn, be affected by them. The cell cycle sums up the challenges of dissecting multiple co-occurring processes in the same single cell, especially when these processes are not fully independent of each other. Recent methods such as Revelio propose to order and remove this effect [START_REF] Schwabe | The transcriptome dynamics of single cells during the cell cycle[END_REF]. The large variation induced by the cell cycle in some scRNA-seq experiments may hide other important sources of biological variation. On the other hand, removing its impact altogether (e.g., regressing it out) may eliminate important facets of biological variation. For example, cellular differentiation processes are strongly linked to the cell cycle. However, some of their aspects may not be apparent in a cell's molecular profile due to the cell cycle's large effect. This kind of reasoning can also be applied to mitochondrial and ribosomal genes, where mitochondrial genes correlate to cell stress [START_REF] Islam | Quantitative single-cell RNA-seq with unique molecular identifiers[END_REF] but still might retain interesting biology.

Imputation

On many occasions, zero expression for a specific gene will be recorded. It might be that these zero counts represent the true biological state of a cell. However, zeros can also result from confounding biological factors such as stage in the cell cycle, transcriptional bursting. Besides, some effects are purely technical factors. The small amounts of starting material and low sequencing depth mean that there are many occasions where no counts are recorded for a specific gene. Sampling effects can result in "dropout" events where a transcript is truly expressed in a sample but is not observed in the sequencing data. In bulk experiments these effects are limited because we average across the pool of cells in a sample. However, for single-cell experiments they can present a challenge for analysis as methods must account for the missing information.

One way to tackle sparsity in single-cell datasets is to use zero-inflated versions of common distributions. However, it is still debatable whether scRNA-seq datasets, particularly those from droplet-based capture protocols, are truly zero-inflated or if the additional zeros are better modeled with standard distributions with lower means [START_REF] Svensson | Droplet scRNA-seq is not zero-inflated[END_REF]. Another approach to take into account dropouts is to impute some of the zeros, replacing them with estimations of how expressed those genes should be based on their expression in closely related cells using Since it is difficult to analyze the data accurately until it is normalized, a wide variety of approaches have been developed for normalization-these range from simple scaling methods to the analysis of generative models of the data.

On the other side of the coin, biological variation such as the cell size should be considered since big cells will have more transcripts. We can also cite transcriptional bursts and mRNA variable degradation rate.

Effective normalization is just as crucial for single-cell experiments as it is for bulk RNA-seq datasets. Fragments Per Kilobase Million (FPKM) or Transcripts Per Kilobase Million (TPM) transformations can be used for datasets without UMIs. However, for UMI datasets, the gene length correction is not needed as reads only come from the ends of transcripts [START_REF] Islam | Quantitative single-cell RNA-seq with unique molecular identifiers[END_REF]. The widely used log-normalization approach, implemented in the Seurat package, applies a scale factor to diminish the importance of extreme values, and highlight robust data trends.

One limitation of this method is that in this case, all cells have the same underlying size. This might not be the case for highly heterogeneous samples. The approach proposed by SCONE is to compare multiple normalization techniques for one dataset while selecting the most appropriate [START_REF] Cole | Performance Assessment and Selection of Normalization Procedures for Single-Cell RNA-Seq[END_REF]. A recent benchmark pointed-out that normalization should be handled with care since this step was having among the most dramatic effect on all the rest of the downstream results [START_REF] Germain | pipeComp, a general framework for the evaluation of computational pipelines, reveals performant single cell RNA-seq preprocessing tools[END_REF].

Cell-free RNA removal

Droplet-based single-cell RNA sequence analyses assume that all acquired RNAs are endogenous to cells. However, any cell-free RNAs contained within the input solution are also captured by these assays. This sequencing of cell-free RNA constitutes background contamination that confounds the biological interpretation of single-cell transcriptomic data.

Approaches like SoupX (Young and Behjati, 2020) quantify the extent of the contamination and estimating "background-corrected" cell expression profiles that seamlessly integrate with existing downstream analysis tools. While it is not clear to which extent datasets or technology are most contaminated by cell-free RNA, applying this kind of method can improve biological interpretation.

Datasets Integration

Early scRNA-seq studies often used only a single sample. However, as technologies have become cheaper and more widely available, it is frequent to see studies with multiple batches or leveraging publicly available data, or even comparing cell populations between species. A range of computational approaches performing sample integration has been developed, including bbknn [START_REF] Polański | BBKNN: fast batch alignment of single cell transcriptomes[END_REF], LIGER [START_REF] Welch | Single-Cell Multi-omic Integration Compares and Contrasts Features of Brain Cell Identity[END_REF], Scanorama [START_REF] Hie | Efficient integration of heterogeneous single-cell transcriptomes using Scanorama[END_REF] and Harmony [START_REF] Korsunsky | Fast, sensitive and accurate integration of single-cell data with Harmony[END_REF]. The alignment approach in the Seurat package (Butler et al., 2018) uses Canonical Correlation Analysis (CCA) to identify a multi-dimensional subspace that is consistent between datasets. Dynamic Time Warping is then used to stretch and align these dimensions so that the datasets are similarly spread along them. A recent update to the Seurat v2 method combines these approaches by applying CCA before using the MNN approach to identify "anchor" cells with common features in the different datasets [START_REF] Stuart | Comprehensive Integration of Single-Cell Data[END_REF]. Batch correction vectors can then be calculated and subtracted using these anchors. Clustering and dimension reduction is then performed using these aligned dimensions, but as the original expression matrix is unchanged, the original matrix is used for tasks such as differential expression testing.

Feature selection

A standard clustering pipeline typically involves a step of highly variable gene selection, complicated by the mean-variance relationship of scRNAseq. Seurat's early approaches involved using dispersion estimates standardized for the mean expression levels, while more recent versions (≥ 3.0) rely on a variance-stabilizing transformation. While adjusting for the mean-variance relationship removes much of the bias towards highly expressed genes, it is plausible that this relationship may sometimes reflect biological relevance and would help classify cell types.

Dimension reduction

Although scRNA-seq datasets have high dimensionality, their intrinsic dimensionality is typically low. Many genes are co-expressed/co-regulated, and a few variables, such as cell type, a gene program, or the number of detected transcripts, can explain a substantial portion of the variation in a dataset. As a result, dimensionality reduction, followed by visualization or downstream analyses, has become a classical strategy for exploratory data analysis in single-cell genomics. Proper dimensionality reduction can allow for effective noise removal and facilitate many downstream analyses that include cell clustering and lineage reconstruction. robust clusters. Some of these methods include SINgle CEll RNA-seq profiling Analysis (SINCERA) [START_REF] Guo | SINCERA: A Pipeline for Single-Cell RNA-Seq Profiling Analysis[END_REF], Single-Cell Consensus Clustering (SC3) [START_REF] Kiselev | SC3: consensus clustering of single-cell RNA-seq data[END_REF]. Methods specifically designed to detect rare populations and tools have been developed, such as CellSIUS [START_REF] Wegmann | CellSIUS provides sensitive and specific detection of rare cell populations from complex single-cell RNA-seq data[END_REF] or GiniClust [START_REF] Jiang | GiniClust: detecting rare cell types from single-cell gene expression data with Gini index[END_REF]. These tools attempt to cluster similar cells together based on their expression profiles, forming groups of cells that are similar to each other. Among all the techniques, graph-based techniques have been more efficient at processing big and noisy datasets than k-means or hierarchical clustering.

The clustering method in the Seurat package has been shown to perform well. A Shared Nearest Neighbours (SNN) graph is first constructed by considering the distance (Euclidean distance by default) between cells in this multidimensional space and the overlap between shared neighborhoods. To define clusters, a community detection algorithm such as Louvain (Blondel et al., 2008) or Leiden [START_REF] Traag | From Louvain to Leiden: guaranteeing well-connected communities[END_REF] detection is then run on the graph, maximizing modularity. Modularity is a metric used to define communities. In Louvain or Leiden, the modularity is maximized using a greedy heuristic. Fine-tuning the algorithm is possible thanks to a resolution parameter that controls the number of produced clusters.

Selecting the resolution parameter is difficult but important as the number of clusters selected can affect the interpretation of results. The tool Clustree can give a visual interpretation of these resolutions by showing clusters relationships on a tree representation as the resolution increases (Zappia and Oshlack, 2018). This tool can notably help determine when further splitting the data does not improve interpretability and leads to a noisy picture, usually because clusters are not robust anymore. This phenomenon can be referred to as over clustering.

Assessing cluster relevance is then a work that necessitates domain knowledge, but some tools still have proven useful. The silhouette score is a metric for cluster robustness. A silhouette measure is assigned to each cell representing how close a cell is to its own cluster compared to other clusters. In this case, a silhouette score can give more confidence to clusters with a high average score. Other tools like ROGUE (Liu et al., 2020) assess cluster purity. With this metric, if the cluster purity is low, it means that we might be able to uncover additional cell subtypes by increasing clustering resolution.

Cell type annotation

Annotating cell types is a critical step in single-cell RNA sequencing (scRNA-seq) data analysis.

As a complement to unsupervised clustering, some supervised or semi-supervised classification Seurat [START_REF] Stuart | Comprehensive Integration of Single-Cell Data[END_REF]. Some of them use user-defined markers like Garnett, however, reliable marker genes are not always known for every cell-type of interest. [START_REF] Huang | Evaluation of Cell Type Annotation R Packages on Single-cell RNA-seq Data[END_REF].

As cellular consortium emerges, this kind of tool will be of growing interest to map small samples into a big reference, leveraging the power of big databases and growing consortium like the Human Cell Atlas [START_REF] Regev | The Human Cell Atlas White Paper[END_REF]. Classification has the advantage of making use of existing knowledge and avoids manual annotation of clusters, which can often be difficult and time consuming. However, it is biased by what is present in the reference datasets used and cannot reveal previously unknown cell types or states.

Pseudotime alignment

It is not always clear whether considering cells as discrete or as a continuous spectrum is the most appropriate approach. Some facets of a cell's identity can be better represented by a continuous phenotypic spectrum, for instance, as cells gradually differentiate from one phenotype into another. Even when a cell population partitions into discrete types, hybrid cells may be observed. In this case, it can make sense to order cells along a continuous trajectory from one cell type to another instead of assigning them to distinct groups. Trajectory analysis was initially introduced by the Monocle package, which used dimensionality reduction and computation of a minimum spanning tree to explore a model of skeletal muscle differentiation [START_REF] Trapnell | The dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells[END_REF]. Since then, a range of others tools have been developed, including TSCAN [START_REF] Ji | TSCAN: Pseudo-time reconstruction and evaluation in single-cell RNA-seq analysis[END_REF], Slingshot [START_REF] Street | Slingshot: cell lineage and pseudotime inference for single-cell transcriptomics[END_REF], ElPiGraph (Chen et al., 2019).

The process of learning trajectories can generally be broken down into two steps. In the first step, dimensionality reduction techniques such as PCA, t-SNE or UMAP are used to project cells into lower dimensions where a graph is constructed between them. The trajectory is then created by finding a path through the cells and ordering the cells along the trajectory.

Both grouping and ordering can be informative, and it is often useful to attempt both on a dataset and see how they compare. The recently developed PAGA tool tries to reconcile clustering with trajectory inference through a topology preserving map of single cells [START_REF] Wolf | PAGA: graph abstraction reconciles clustering with trajectory inference through a topology preserving map of single cells[END_REF]. This tool has compared favorably to the other trajectory tools in computation time and preserving the original dataset's topology [START_REF] Saelens | A comparison of single-cell trajectory inference methods[END_REF]. In that regard, PAGA seems to be a good first approach to try.

Marker gene detection

Once cells are annotated by clustering or reference-based annotation, the next task is to interpret what these groups represent. For clustered datasets, this is usually done by identifying genes that are differentially expressed across the groups or marker genes expressed in a single cluster. Many methods have been suggested for testing differential expression, some of which take into account the unique features of scRNA-seq data, for example, D3E (Delmans and Hemberg, 2016), DEsingle [START_REF] Miao | DEsingle for detecting three types of differential expression in single-cell RNA-seq data[END_REF], MAST (Finak et al., 2015), SCDE [START_REF] Kharchenko | Bayesian approach to single-cell differential expression analysis[END_REF]. The large number of cells in scRNA-seq datasets means that some of the problems that made standard statistical tests unsuitable for bulk RNA-seq experiments do not apply.

Simple methods like the unpaired Wilcoxon rank-sum test (or Mann-Whitney U test), Student's t-test or logistic regression may give reasonable results in this setting. A comprehensive evaluation of 36 differential expression testing tools found that methods developed for bulk RNA-seq as well as scRNA-seq specific methods [START_REF] Soneson | Bias, robustness and scalability in single-cell differential expression analysis[END_REF]. The authors also observed that some methods were biased in the types of features they tended to detect as differentially expressed.

When cells have been ordered along a continuous trajectory, the task is slightly different.

Instead of testing for a difference in means between two groups, the goal is to find genes with a relationship between expression and pseudotime. Fitting splines can accomplish this to the relationship between pseudotime and expression and testing the fitting coefficients. Genes that are associated with a trajectory are important as they describe the biology along a path.

Gene set analysis / Gene signatures

Interpreting the meaning of detected marker genes is a difficult task and is likely to remain so.

Gene set testing to identify related categories such as Gene Ontology, KEGG, WikiPathways,

Alternative analyses

The analyses presented here are not directly in the standard single-cell workflow. While these analyses might not be amenable or relevant to every single-cell experiment, these represent interesting propositions and highlight the diversity of approaches and creativity of the single-cell community.

Cell-cycle assignment

As we have previously seen, Cell cycle variation is a common source of variation in single-cell RNA-seq data. If we did not remove this signal in the preprocessing step, it is possible to assign cells to their cell-cycling phase. We can cite cyclone from the scran package (L. [START_REF] Lun | Pooling across cells to normalize single-cell RNA sequencing data with many zero counts[END_REF] or CellCycleScore, a function directly embedded into Seurat [START_REF] Stuart | Comprehensive Integration of Single-Cell Data[END_REF]) that assign to each cell a score based on its expression of G2/M and S phase markers. The traditional approach to classify and sort cells into distinct cell cycle states relies on a few known markers and quite arbitrary gating cutoffs. Most cells of any given non synchronized culture do not, in fact, show an unambiguous signature of being in one of the standard discrete cell cycle phases.

We know that the gene expression of cell state transitions does not occur in abrupt steps but rather in a continuous process. To address this specificity, a predictor of a continuous cycling state (Hsiao et Analyzing scTCR data 10x Genomics 5' V(D)J kit allows clonotypes and transcriptomes to be co-assayed in the same individual cells. This kind of data has specific considerations and should be considered on its own. In the case of single-cell TCR (scTCR data), Cellranger is also used to preprocess scTCR-seq datasets. This pipeline takes FASTQ files as input and proceeds as follows.

Preprocessing

The main function of the cellranger pipeline is cellranger vdj, which is a wrapper around multiple processing steps. Cellranger vdj takes FASTQ files from cellranger mkfastq or bcl2fastq for V(D)J libraries and performs sequence assembly and paired clonotype calling. It uses the Chromium cellular barcodes and UMIs to assemble V(D)J transcripts per cell. Here are the sequential steps the algorithm goes through.

-Assembly Algorithm -Assembles read-pairs into one contig for each transcript sequence.

-VDJ Cell Calling Algorithm -Identifies barcodes that contain T or B cells.

-Annotation Algorithm -Annotates contigs with V(D)J segment labels and locates CDR3 regions.

-Clonotype Grouping -Groups cells into clonotypes and in doing so filters out some cells. At the end of the pipeline, cellranger will output a folder with two main files of interest. The "Filtered contigs" file with each contig defined for each single-cell, and the "clonotypes" file will contain merged clonotypes definition, with their frequencies, and TCR alpha and beta sequences. The "clonotype_id" column name can be used to match the two main files.

Single-cell TCR analysis is a relatively new field. As such, the number of tools to process, compare, and analyze this kind of data is quite low as compared to scRNA-seq tools, for instance. Still, some automated pipelines have been recently developed, like Scrypy, as an extension of Scanpy [START_REF] Sturm | Scirpy: a Scanpy extension for analyzing single-cell T-cell receptor-sequencing data[END_REF], or scRepertoire as an extension of Seurat [START_REF] Borcherding | scRepertoire: An R-based toolkit for single-cell immune receptor analysis[END_REF], which are streamlining some of the processes. Here are some of the most classical types of analyses for this kind of data.

scTCR-seq analyses

Clonality One of the most straightforward analyses is to look for clonotype abundance by cluster and/or condition. In this case, the total number of different clonotypes are measured, as well as the clonal expansion sizes. One of the ways to estimate the diversity of samples is to evaluate clonality. It measures the amount of the most or the least frequent clonotypes. There are several methods to assess clonality, either in absolute numbers, or considering the relative abundance of clones among the total population.

Length distribution of CDR3

The TCR α and β chains possess three hypervariable regions termed complementarity determining regions (CDR1, 2 and 3). CDR3 is responsible for recognizing processed antigen peptides. As the sequence and length of CDR3 differs according to the type of T cell clone, the sequence of CDR3 determines the structure and specificity of the TCR, where one type of CDR3 sequence represents a specific T cell clonotype. With scTCR-seq data, it is possible to analyze CDR3 polymorphisms and sequence length diversity, in order to investigate the pattern of TCR utilization.

Public clonotypes shared between repertoire

Repertoire overlap is the most common approach to measure repertoire similarity. It is achieved by computation of specific statistics on clonotypes shared between given repertoires, also called "public" clonotypes. Overlap can be quantified using many metrics, such as cosine similarity, Tversky index, Jaccard index, Morisita's overlap index. These data can be compared with previously identified TCR sequences matching specific epitopes. A growing number of public TCR studies are now available, such as the databse TCRdb that can be queried TCR sequences with known specificities [START_REF] Chen | TCRdb: a comprehensive database for T-cell receptor sequences with powerful search function[END_REF].

Gene usage computation

To compute the distributions of genes, immunarch includes the geneUsage function. It receives a repertoire or a list of repertoires as input and genes and species for which you want to get the statistics. E.g., if you plan to use TRBV genes of Homo Sapiens Diversity of repertoires Several approaches to the estimation of repertoire diversity are implemented, sometimes derived from ecology and information theory. These indexes comprises :

-Chao1, a nonparametric asymptotic estimator of species richness, -Gini-Simpson index, the probability of interspecific encounter, i.e., probability that two entities represent different types -Inverse Simpson index, is the effective number of types that is obtained when the weighted arithmetic mean is used to quantify average proportional abundance of types in the dataset of interest.

-Gini coefficient, which measures the inequality among values of a frequency distribution (for example levels of income). K-mer analysis K-mer analysis is a way to analyse recurring sequence motif analysis of length K among CDR3 regions (2-mer, 3-mer etc.). Many repertoire-based approaches rely on k-mer encoding as this enables the establishment of a common comparative framework across individuals, which is otherwise difficult when considering entire sequences, which are mostly unique to an individual.

For the analysis of TCR sequences between tissues and clusters, the tool STARTRAC (Zhang et al., 2018) is able to quantitatively analyse the dynamic relationships between clusters. In distribution, expansion, migration and transition, to quantitatively describe tissue distribution, clonal expansion, migration and developmental transition or differentiation, respectively. It is now also possible with Tcellmatch to predict antigen specificity of single T cells based on TCR CDR3 regions [START_REF] Fischer | Predicting antigen specificity of single T cells based on TCR CDR3 regions[END_REF]. This approach is using deep learning models to predict the binding of T-cell receptors (TCR) to peptide MHC complexes (pMHC) from defined antigen panels. It is expected that the field of single-cell TCR analyses will grow as interest in immunotherapies continues to expand.

Thesis overview and aims

The full spectrum CD8+ T cell heterogeneity is still being explored. As these cells patrol blood and tissues, as antigens activate them or as external factors halt them, these cells display a wide range of functionality. Understanding CD8+ T cell dynamics is critical, notably to understand the processes at play during immunotherapy.

The tools for analyzing single-cell datasets have developed rapidly along with the technology. A continuous task during my Ph.D. was to keep myself updated on the single-cell literature, evaluate tools' relevance to our biological system, and finally try these new methods. One good example of this is the emergence of sample integration methods. In 2017, when my Ph.D. started, integration methods were not available. As the paper "The epigenetic control of stemness in CD8 + T cell fate commitment" contains two different conditions (Wild Type (WT) and Suv39h1-KO), these samples are analyzed separately. Implementing integration in my pipeline (first with Seurat v2, then Seurat v3) allowed me to compare conditions on the same embedding and the same clustering. Integration methods also led to an improvement in cluster robustness because of the increased number of cells overall. This was also the case with droplet calling. Before it was included in the Cellranger pipeline by default, Emptydrops increased the number of cells analyzed up to 30% (cf Appendix 4). Combined with integration methods, this significant change allowed me to better resolve T cell subpopulations, especially small populations like γδ T cells.

Ultimately, I applied those methods to scRNA-seq and scTCR datasets to understand better the dynamics and differentiation paths of CD8+ T cells, especially when they infiltrate tumors. I have completed these aims in the following ways:

Introduction

The paper presented in this chapter is "The epigenetic control of stemness in CD8 + T cell fate commitment", published in Science in 2018.

When I began my Ph.D. in mid-2017, single-cell RNA sequencing was only ever used in another project in Curie [START_REF] Goudot | Aryl Hydrocarbon Receptor Controls Monocyte Differentiation into Dendritic Cells versus Macrophages[END_REF]. I was able to leverage my previous experience with scRNA-seq data. Still, much of the preprocessing, quality control steps, visualization and analysis tools needed to be tested and refined. I'm presenting here an example of the classical single-cell pipeline, including clustering, differential and compositional analyses, and gene sets analyses.

Memory T lymphocytes provide lifelong protection against pathogens and cancer. In contrast to naïve and effector T cells, memory cells possess unique properties of "stemness," enabling long-term survival and plasticity to replenish effector pools after renewed antigen challenges.

Understanding the lineage relationships among naïve, effector, and memory T cells and the molecular pathways that regulate gene expression during the transitions from one to another of these distinct states, is essential for the rational design of vaccines and the development of new immune-therapeutic protocols.

Although many studies have characterized the transcription factors that control the differentiation of T cells, the corresponding epigenetic states and associated chromatin dynamics involved in the establishment and maintenance of CD8+ T cell memory and effector identities is still incompletely understood. Several epigenetic pathways, including trimethylated histone H3 Lys9 (H3K9me3)/HP-1/Suv39h1 and Polycomb repressive complexes, can contribute to establishing or maintaining transcriptional silencing. The H3K9me3 modification is considered to be a repressive mark, a hallmark of both constitutive and facultative heterochromatin, most often associated with silent gene loci. Mouse Suv39h1 and Suv39h2, two H3K9 site-specific histone methyltransferases, are critical heterochromatin regulators. Suv39h1 is involved in heterochromatin organization, gene silencing, and lineage stability. It also limits somatic reprogramming of differentiated cells into induced pluripotent stem cells. In B cells and CD4+ T cells, Suv39h1 is involved in gene silencing and lineage plasticity. Although the mechanisms underlying the induction of genes critically involved in effector and memory T cell generation have been extensively analyzed, the impact of heterochromatin-dependent gene expression silencing on the fates of T lymphocytes during differentiation has not been addressed. Here, we explore the role of Suv39h1-dependent gene silencing in the establishment and maintenance of memory CD8+ T cell stemness, plasticity, and transition to terminally differentiated effectors.

Conclusion

In this study, we argue that after priming, cycling CD8+ T lymphocytes reprogram both self-renewing and effector gene expression profiles. These cycling cells may represent bipotent intermediates, which would then repress either the effector or stem cell/memory programs while they differentiate to memory precursors or effectors, respectively. The silencing of the In this project, cell selection with Emptydrops led to a dramatic increase in the number of cells retrieved. This led to a much clearer picture of our closely related CD8+ T cell subsets. Assessing sample merging strategies was also key. I compared Seurat v2 and v3 with TCR "ground truth". I noticed that version 3 was indeed much better at conserving heterogeneity. Small subsets like B cells, Interferon high cells appeared clearly with Seurat v3, even at low clustering resolutions. I could also see the breaking point of the algorithm, where merging too diverging datasets led to overcorrection by the algorithm. This emphasises the importance of carefully considering decisions during analysis, especially preprocessing. Manually checking, visualising, and getting a feel of the data was critical here. The code used to generate the figures of this paper is available here: https://github.com/p-gueguen/single_cell_NSCLC

Conclusion

This study leveraged careful sample merging to define high resolution clusters. TCR sharing was then to inform T cell differentiation and circulation between compartments. We describe how CD8+ T cells in lung tumors differentiate through convergent pathways to an expanded and actively cycling terminal state. This work led to a comprehensive annotation of tumor-infiltrating T cells, and how they relate to each other in primary NSCLC.

Taking a step back, this study could be useful to design tumor liquid biopsy as a diagnostic tool. Indeed, recent evidence shows that circulating T cells identified by TCRs shared with T cells in tumors have unique transcriptional expression patterns [START_REF] Lucca | Circulating Clonally Expanded T Cells Reflect Functions of Tumor Infiltrating T Cells[END_REF]. In this case, one could imagine tracking CD8+ precursors signature as a readout of the ongoing immune response. This work also sheds light on the process of tissue-related exhaustion. We show that the resident cells seem to be more likely to differentiate into terminal stage than their circulating counterpart.

One could argue that datasets (especially big ones like scRNA-seq) are rarely fully analyzed, and we only chose to stop analyzing it. More work could have been done on comparing CD4+ vs CD8+ exhaustion, since we see similar ICP patterns in both Tfh and Tregs specific clusters.

Another remarkable feature is that γδT cells also express exhaustion features. This wasn't really explored to the best of my knowledge. Another question in the field is to characterize CD4+ tissue residency features, especially in the context of tumors. TLS status on matched patients slides are currently ongoing and could shed light on relationships between TLS status and specific subpopulations. Finally, to understand the full picture, work remains to understand T cell dynamics on other lymphoid compartments, like draining lymph nodes or ascites, ideally before and after immunotherapy treatment.

Since no consensus currently exists on CD8+ T cell subsets definitions, cluster annotation required manual annotation and curation, which was both time-consuming and prone to errors. This initial laborious step happened to be a blessing in disguise and highlighted the need for a common framework to understand and annotate T cell function. Thanks to careful patient and tissue merging, I believe I made an important contribution to this question during my thesis.

Comprehensive annotation of tumor-infiltrating T cells is currently being investigated by other groups [START_REF] Andreatta | Projecting single-cell transcriptomics data onto a reference T cell atlas to interpret immune responses[END_REF][START_REF] Nieto | A Single-Cell Tumor Immune Atlas for Precision Oncology[END_REF].

My thesis contains three projects that explored different facets of CD8+ T cell heterogeneity in acute viral infection or cancer, in both mice and humans. Chapter 3 analyses revealed a defect in the silencing of stem/memory genes selectively in Suv39h1-defective T cell effectors. Moreover, when analyzing cycling subpopulations, scRNA-seq suggests that cycling cells may represent bipotent differentiation intermediates expressing both effector and stem/memory potential.

Chapter 4 described an IFN-a response signature and a higher cytolytic potential upon PD-1 blockade. The anti-PD-1 response signature of Suv39h1-deficient CD8+ T cells is discriminating response to anti-PD-1 in melanoma patients. Chapter 5 described two CD8+ precursor populations in primary NSCLC: one circulating through the blood, and one resident in tissues.

We show how these two precursor populations differentiate through convergent pathways to an expanded and actively cycling terminal state. Overall I believe I have made a significant contribution to the field of T cell heterogeneity during my Ph.D.

Future of single-cell technologies in immunology

While powerful, scRNA-seq remains one technique among many others, with its flaws. Discoveries will always be fundamentally limited by what can be measured. Notably, although various protocols exist for scRNA-seq analysis, nearly all involve poly-A selection, thereby limiting the ability to examine non-polyadenylated transcripts, such as small nucleolar RNAs, pre-mRNAs, long noncoding RNAs, and transposable elements which may play diverse regulatory roles in cancer. Even within poly-A selection scRNA-seq protocols, droplet-based protocols that are restricted to only 3′ or 5′ ends will inherently be more limiting for allele-based CNV inference, mutation calling, and fusion gene detection compared to full transcript single-cell RNA-seq protocols. In that regard, orthogonal validation, and wet-lab experiments can not be replaced.

The numbers of analyzed cells and samples should continue to increase exponentially, increasing results' robustness. This is largely due to international and collaborative efforts such as the Human Cell Atlas [START_REF] Regev | The Human Cell Atlas White Paper[END_REF], Allen Brain Atlas [START_REF] Sunkin | Allen Brain Atlas: an integrated spatio-temporal portal for exploring the central nervous system[END_REF], Pediatric Cell Atlas [START_REF] Taylor | The Pediatric Cell Atlas: Defining the Growth Phase of Human Development at Single-Cell Resolution[END_REF], HuBMAP [START_REF] Snyder | Mapping the Human Body at Cellular Resolution --The NIH Common Fund Human BioMolecular Atlas Program[END_REF], Human Tumor Atlas Network (Rozenblatt-Rosen et al., 2020), LifeTime EU Flagship [START_REF] Rajewsky | LifeTime and improving European healthcare through cell-based interceptive medicine[END_REF]. As such, there will be a growing need to improve the scalability of computational methods through implementation improvement through cloud computing and parallelization and algorithmic optimization.

On the computational side, benchmarking and comparison studies remain vital to understand which tools are effective, but they can be difficult and time consuming to perform. For such complex systems as cells, it can be hard to construct a gold-standard dataset. The lack of reference datasets can be mitigated by performing simulation experiments, like with the splatter package [START_REF] Zappia | Splatter: simulation of single-cell RNA sequencing data[END_REF]. Benchmarking efforts will continue to drive progress in the single-cell field.

New platforms and protocols are being developed and used to generate even more modalities at the same time. Experiments that combine scRNA-seq with other measurements such as chromatin accessibility or protein abundance will become more common. Importantly, spatial transcriptomics technologies can retain a cell's position in tissue while measuring its transcriptional profile [START_REF] Wang | Spatial organization of the transcriptome in individual neurons[END_REF]. Such spatially resolved single-cell transcriptomic data generated from these different technologies will require new computational pipelines and methods for proper processing (like cell segmentation or deconvolution for Visium datasets).

For solid tumors, spatially resolved technologies have the ability to add phenotypic dimensions to our understanding of how cells interact locally in the tumor immune microenvironment and will help elucidate the mechanisms behind immunotherapies.

Single-cell technologies present a great avenue to lay the foundation for new therapeutic innovations. Breakthroughs were already made, like the discovery of a rare cell type in the lung: the pulmonary ionocyte (Plasschaert et al., 2018). This cell type is the most expressing CFTR, which is the mutated gene in cystic fibrosis. It is also possible to look back and unlock previous knowledge. This is the case in genome-wide association study (GWAS) studies, in which we didn't previously know the cell types implicated in genetic variants associated with a trait. On a broader scale, single-cell technologies pave the way to conceptualize what are cell type and cell states (and eventually if this distinction is relevant). The field might be ripe for such classification, as we currently don't have the equivalent of the periodic table for cell types. The high number of cells probed allows us to have a noisy but broad view of cellular heterogeneity, leading to data-driven classifications. As we've seen in other fields, the human mind is likely to be stuck in local minima [START_REF] Bengio | Evolving Culture Versus Local Minima[END_REF], and algorithms might be better tailored at picking up on non-linear trends, like in radiology images [START_REF] Mei | Artificial intelligence-enabled rapid diagnosis of patients with COVID-19[END_REF]. The hope is to have a data-driven way to understand the big picture: how cells organize themselves in high dimensional space. Algorithms can actually pick up on covariation patterns, and this is why the number of cells was able to mitigate the relatively low RNA coverage (only ~5-10% of the total RNA of a cell is sequenced). Using single-cell technologies, we can go both deep in scope (like with molecular biology, albeit at a lower resolution for now) and broad (like with genomics), while previously we needed to choose. This represents a new way to investigate biology. It is now possible to address biological problems with a computational mindset, and actually devise experiments based on computational ideas. In that regard, single-cell biology can be used to construct a novel language of health and disease.

Although immunotherapies have transformed the management of many advanced-stage malignancies, the majority of patients either do not respond to these agents or will ultimately acquire resistance. Thousands of clinical trials exploring combination immunotherapy approaches, with other immunotherapy agents, with non-immunotherapy agents or both, have either been completed or are currently underway. Nonetheless, the design of rational immunotherapy-containing combinations designed to overcome primary and/or acquired resistance will require the comprehensive characterization of the tumor immune microenvironment of each cancer type, and a detailed understanding of the determinants of effective anti-tumor immunity in the context of therapy.

Single-cell technologies enable comprehensive profiling of the TME and are therefore especially well-suited to studying tumor immune cell heterogeneity and how such differences in tumor biology contribute to therapeutic resistance. In order to effectively leverage single-cell technologies to improve our understanding of mechanisms of immunotherapy response and resistance and to aid the development of novel therapeutics, clinical trials will have to be carefully designed to incorporate single-cell approaches.

If clinical studies are carefully designed in collaboration with translational scientists, incorporating single-cell analysis in a subset of patients could yield informative results that can be applied to larger cohorts. Careful consideration must be given to what types (among others, tumor tissue, non-malignant adjacent tissues and/or peripheral blood) and formats (fresh or frozen) of samples should be collected, the optimal timing of collection (pretreatment, on-treatment and/or post-treatment). While challenging, the insights gleaned from an in-depth investigation of a relatively small number of patients can inform the analysis of a larger sample set using a more cost-effective and logistically more straightforward method.

An emerging approach is to exploit published scRNA-seq data rather than to generate these data de novo. Mining high-dimensional single-cell data to define low-dimensional panels of genes and/or proteins that show promise as predictive biomarkers provide a clear path towards the clinics. Another avenue for exploiting data-rich single-cell information is provided by the use of bioinformatics algorithms to map bulk RNA-seq data onto the cell types defined by scRNA-seq and related techniques. The use of these in silico dissection approaches is less subject to the compromises inherent in using low-dimensional panels alone. A good example of this approach is provided by CIBERSORTx [START_REF] Newman | Determining cell type abundance and expression from bulk tissues with digital cytometry[END_REF], a machine learning-based method that takes a signature matrix generated from scRNA-seq data and applies it to bulk tissue RNA-seq profiles in order to estimate the composition and abundance of cell types in each tissue sample.

Ultimately, the design of rational immunotherapy-containing combinations designed to overcome primary and/or acquired resistance will require the comprehensive characterization of the tumor immune microenvironment of each cancer type, and a detailed understanding of the determinants of effective anti-tumor immunity in the context of therapy. full Rapid flow cell, coverage was around 100M reads per sample corresponding to 100000 reads per cell.

Single-cell RNA-seq analysis

Single-cell expression was analyzed using the Cell Ranger Single Cell Software Suite (v2.0.1) to perform quality control, sample de-multiplexing, barcode processing, and single-cell 3′ gene counting. Sequencing reads were aligned to the mm10 transcriptome using the Cell Ranger suite with default parameters. Samples were merged using Seurat CCA algorithm (Butler et al. 2018). A total of 15121 single cells were analyzed. Further analysis was performed in R (v3.4) using the Seurat package (v2.3) (Satija et al., 2015). The gene-cell-barcode matrix of the samples was log-transformed and filtered based on the number of genes detected per cell (any cell with less than 400 genes or more than 6000 genes per cell was filtered out). Any cell with more than 6% of mitochondrial UMI (Unique Molecular Identifiers) counts and more than 50% of ribosomal UMI was filtered out. Genes detected in at least three cells were included. Cells were then scaled to a total of 1e4 molecules.

CCA merging was performed using Seurat CCA (alignment score 0.87) and accounting for the 1000 most variable genes, and taking into account the 15 first CCs for downstream analysis. These canonical components were used for the UMAP projection and clustering analysis. Clusters were identified using the "Find_Clusters" function in Seurat with a resolution parameter of 0.7. Clustree analysis was performed using Clustree R package (Version 0.2.2) (Zappia and Oshlack, 2018). 16 clusters were initially defined, and are visualized using UMAP in Fig. S4A). To simplify the analysis, some clusters which shared defined signatures and which did not evidently change among the different conditions were grouped as follow: original clusters 1, 3, 6, 7 and 11, corresponded to cycling cells in the G1/S and G2/M phase (Fig. S4B-C), were merged into one cluster (cluster 8 "cycling cells", Fig. 3B); and original clusters 0, 8, 12, 13 and 14 (in Fig. S4A) which shared a "late exhaustion" signature, were merged together (as cluster 7 "late exhausted", Fig. 3B).

Single-cell RNA/TCR seq data processing

We used emptyDrops (39) function from R package dropletUtils with default parameters when analyzing the samples processed with version Cellranger v2.0.1. The emptyDrops function is implemented starting from Cellranger v3.0. Further analysis was performed in R (v3.5.1) using the Seurat package (v3.1.1) (18). Cells were then filtered out when expressing less than 500 genes or more than 5000 genes, or when expressing more than 10% mitochondrial genes.

Altogether, among the 17 samples (11 from tumors, 6 from blood and 2 from juxta-tumor) 54247 cells were kept for statistical analysis. For each sample, the gene-cell-barcode matrix was then normalized to a total of 10^4 molecules. Specifically, UMI counts for each gene in a given cell are divided by the total UMI for that cell and multiplied by a scale factor of 10^4. This value is then log transformed using log1p. The top 2000 variable features were identified using the "vst" method from Seurat where both lowly and highly expressed genes are transformed onto a common scale. To reduce the technical effect associated with 5' versus 3' technologies (Sup Fig. S1) and to correct for other batch effects between samples, we used the integration method from Seurat v3 on all 11 tumor samples simultaneously. This method is leveraging closely-related cells (termed anchors) between datasets to compute a batch-corrected matrix (18). Default parameter values as well as the top 30 CCA components were used to find transfer anchors between datasets and used to generated the integrated matrix for the 11 tumors.

Dimension reduction and unsupervised clustering and contaminant removal

Top 30 Principal Components were computed and UMAP was performed using the top 30 PCs of the integrated matrix. Clusters were identified using the FindNeighbors and FindClusters function in Seurat v3.1.1 with a resolution parameter of 1.2 and using the first 30 principal components. To choose the optimal number of clusters and prevent overclustering, clustree analysis was performed using the clustree package v0.4.1 (40). Unique cluster-specific genes were identified by running the Seurat FindAllMarkers function using Wilcoxon test on the uncorrected matrix. Potential contaminant clusters were identified based on lack of CD3 expression. Investigation of marker genes in these clusters identified myeloid cells expressing APOE and S100A8, B cells expressing CD19 and MS4A1, epithelial cells expressing KRT1 and KRT10, and alveolar cells expressing SFTPA1 and SFTPB. All cells from these clusters were removed from further analysis. Signature scores (eg Figure 1E, 2A, 2B) were computed using the Seurat function AddModuleScore using the gene signature of interest and the integrated matrix and setting the number of control genes from the same bin of expression at 5. This function calculates for each individual cell the average expression of each gene signature, subtracted by the aggregated expression of control gene sets matched for individual gene expression level. All the signatures used to compute signature scores, as well as corresponding references, are listed in Table S2. The average heatmap was generated using the Average function in Seurat which is averaging gene expression across clusters. Cell cycle scoring was also performed using Seurat CellCycleScore function using cell cycle genes.

Label transfer using a reference

Fig. 1 .

 1 Fig. 1. CD8 + Tcell-mediated host protection is impaired in Suv39h1defective mice. (A) Experimental design. (B) Littermates and Suv39h1-KO mice were intravenously (i.v.) injected with LM-OVA (primary infection) and, on the indicated days post-infection (p.i.), the number of bacterial colonyforming units (CFU) was determined. (C) Littermates and Suv39h1-KO mice previously i.v. immunized with LM-OVA were rechallenged 48 days later. Three days after LM-OVA secondary infection, protection was assessed by counting CFU in spleen and liver. (D) Littermate and Suv39h1-KO mice were immunized with LM-OVA; 7 days later, primary CD8 + T cell responses were evaluated in the peripheral blood, using K b -SIINFEKL (K b -OVA + )m u l t i m e r s

Fig. 2 .Fig. 3 .

 23 Fig. 2. Gene expression patterns and differentiation programs of K b -OVA + CD8 + T lymphocytes are enriched in stem cell-like gene signatures in Suv39h1-KO mice. (A) Experimental design. Naïve littermate and Suv39h1-KO CD8 + T cells and day 7 p.i. K b -OVA + CD8 + Tcells (from LM-OVA-infected mice) were isolated by FACS. RNA was isolated and analyzed using Affymetrix microarrays (three mice per condition were analyzed). (B) Volcano plots of K b -OVA + CD8 + T cells versus naïve CD8 + Tcells from littermate and Suv39h1-KO mice show the adjusted P value (-log 10 ) versus fold change (log 2 ). Up-regulated and down-regulated mRNAsareshowninredandblue, respectively. (C) Venn diagrams summarize the overlap between differentially expressed genes that are up-regulated (left) or downregulated (right) in K b -OVA + CD8 + Tcells versus naïve CD8 + Tcells from wild-type and Suv39h1-KO mice. Total common gene numbers for each group are indicated within the areas. (D) Expression pattern of mRNAs down-regulated in wildtype K b -OVA + CD8 + T cells versus naïve cells, shown for wild-type and Suv39h1-KO CD8 + T cells. Representative genes are shown with red lines. (E)GSEAwasperformedto determine the specific enrichment in gene signatures (GeneSet) in wild-type and Suv39h1-KO K b -OVA + CD8 + T cells isolated on day 7 after LM-OVA infection. (F)GSEAfor stem/memory signature (table S1) in wild-type and Suv39h1-KO K b -OVA + CD8 + T cells. The top highly expressed mRNAs in Suv39h1-KO or in wild-type K b -OVA + CD8 + Tcells are shown. (G) BubbleGUM analysis of homologies between different CD8 + T subsets by high-throughput GSEA. The right panel summarizes the normalized enrichment score (NES) and false discovery rate (FDR) parameters obtained by GSEA of wild-type and Suv39h1-KO K b -OVA + CD8 + T cells, as in (F). The different gene signatures were analyzed in wild-type K b -OVA + CD8 + T cell subsets, isolated on day 7 after LM-OVA infection. Left: Sorted subsets; right: all possible pairwise comparisons.

PaceFig. 4 .

 4 Fig. 4. Suv39h1 is a critical regulator of peripheral effector versus memory CD8 + T differentiation. (A) Littermate and Suv39h1-KO mice were infected with LM-OVA. Seven days later, dump -CD44 hi K b -OVA + CD8 + T cells in the spleen were analyzed according to CD127 and KLRG1 staining (CD127 hi KLRG1 -memory precursors and CD127 lo/-KLRG1 + differentiated effector cells). Representative dot plots are shown; numbers represent the percentages. (B) Numbers of memory precursors and effectors in the spleen were measured. (C) Wild-type and Suv39h1-KO mixed bone marrow chimeras were infected with LM-OVA, and 7 days later the ratio between memory precursors and short-lived effector CD8 + T cells was evaluated in peripheral blood using K b -OVA + multimers. (D to F) Mice were infected with LM-OVA on day 0 (primary infection) and day 48 (secondary infection). On days 7, 18, 47, and 51 p.i. (with respect to the primary infection), numbers of short-lived effectors (D) and central memory cells (E) on gated dump -CD44 hi K b -OVA + CD8 + Tcells per ml of peripheral blood were measured. (F) Longitudinal analysis of the dynamics of central memory K b -OVA + CD8 + Tcells. (G) Representative dot plots of gated blood CD8 + Tl y mphocytes; numbers represent percentages. (H) Percentage of blood endogenous polyclonal central memory CD8 + Tcells. (I) Peripheral blood cells from LM-OVA-infected littermate and Suv39h1-KO mice were stimulated ex vivo 7 days after infection with the OVA peptide SIINFEKL (OVA 257-264 ) and analyzed for intracellular T-bet expression and IFN-g production. Representative dot plots are shown; numbers represent percentages. Data are shown as geometric means in (B), (D), and (E) or as m e a n si n( C ) ,( F ) ,a n d( H ) .* P < 0.05, **P < 0.01, ***P < 0.001 (Wilcoxon Mann-Whitney test).

  and fig.S9, left). The proportion of central memory (CD44 hi CD62L + CD127 + K b -OVA + )C D 8 + T cells(27)i ncreased over time in Suv39h1-K Om i c ea sc o mpared to wild-type mice (Fig.4F). Notably, the percentage of endogenous subsets of memory CD8 + TcellsinSuv39h1-KO mice was also higher, including polyclonal CD62L + central memory CD8 + T cells in the blood and secondary lymphoid organs (Fig.4, G and H, and fig.S10, Ato C), CD62L -effector memory CD8 + T cells in the bone marrow and blood (fig.S10D), and tissueresident CD8 + T cells in the liver (fig.S10, F to H). The endogenous Suv39h1-KO central memory CD8 + T cells also expressed increased levels of memory and stem cell-like memory markers SCA-1 and CD95 (fig.S10E) (22, 27, 28) Thus, Suv39h1-KO mice show higher levels of central and effector memory CD8 + T cells both before and after LM-OVA challenge.

Fig. 5 .

 5 Fig. 5. Suv39h1-mediated defective silencing during CD8 + Tcell lineage commitment results in the accumulation of stem cell-like central memory cells. (A) Experimental design. Dump -K b -OVA + CD8 + T cells or CD44 hi CD62L + CD127 + K b -OVA -central memory CD8 + Tcells, isolated 7 days after LM-OVA infection from littermate or Suv39h1-KO mice, were adoptively transferred to naïve congenic CD45.1 + recipient mice. The recipients were challenged with LM-OVA 40 days later. (B) Donor CD45.2 + K b -OVA + CD8 + T cells were analyzed in the spleen 4 days after infection. Representative dot plots and a phenotypic analysis of subsets are shown; numbers represent percentages. (C) Numbers of total donor CD45.2 + CD8 + T cells. (D)Numbers of donor CD127 + KLRG1 -memory precursors and CD127 lo/-KLRG1 + effectors from donor K b -OVA + CD8 + T cells. (E) Donor CD45.2 + central memory CD8 + T cells were analyzed in the spleen 4 days after LM-OVA challenge. Representative dot plots with subset phenotypic analyses are shown; numbers represent percentages. (F) Representative histograms of stem cell-like and memory markers on total gated donor central memory CD8 + T cells analyzed 4 days after LM-OVA challenge. (G) Total numbers of adoptively transferred donor central memory cells, before and after LM-OVA challenge. (H)Numbers of donor CD127 + KLRG1 -memory precursor and CD127 lo/-KLRG1 + effector CD8 + Tcells differentiated from donor central memory CD8 + Tcells. Graphs are representative of two experiments with three mice per group. All graphs show geometric means. *P < 0.05 (Wilcoxon Mann-Whitney test).

Fig. 6 .

 6 Fig. 6. Stem cell-like/memory genes are expressed by memory precursors andcyclingCD8 + T intermediates, and silenced by Suv39h1 in terminal effectors. (A to H) scRNA-seq analysis of naïve and K b -OVA + CD8 + T cells isolated from wild-type and Suv39h1-KO mice 7 days after LM-OVA infection. (A) Graph-based clustering by t-SNE projection of naïve and K b -OVA + CD8 + Tcells isolated from wild-type and Suv39h1-KO mice. The colors indicate sorted CD8 + T cell subsets. Each dot represents an individual cell. (B) Heat maps of six representative genes for each subset category identified. Columns represent cells; rows represent genes. The color scale is based on z-score distribution. (C) t-SNE projection of single cells, colored according to four major subset categories defined by semi-supervised clustering based on specific and distinct gene expression profiles. (D) Numbers and percentages of cells for each category. Violin plots in (E) and (F) show expression distribution of memory, effector, cycling, and stem cell/memory representative markers. (G) t-SNE projection of single cells, showing the number of stem cell-like memory genes expressed in each cell. (H) Density scatterplots represent numbers of effector versus stem cell/memory genes expressed in individual cells for each subset category. Color scale indicates gene density; marginal distributions show memory/stem (x axis) and effector ( y axis) gene number distributions. Pearson linear regression is displayed for each plot. (I) Working model depicting the pivotal role of Suv39h1 during CD8 + T cell lineage differentiation and commitment.

  ). It is most likely that the possibly reversible silencing of effector gene expression in memory cells occurs through other mechanisms, as memory cells do effectively reprogram into effectors upon rechallenge. These results open new perspectives for the manipulation of epigenetic programming of T lymphocyte identityi nt h ec o n t e x to fv a ccination, checkpoint-based immunotherapies, and adoptive T cell therapies. the manuscript; the flow cytometry platform, the genomic and animal facilities, the Nikon Imaging Center of Curie Institute; and the PICT-IBISA@Pasteur Imaging facility. We thank T. Jenuwein for providing the SUV39H1 gene and the Suv39h1-KO mice. LM-OVA was kindly provided by H. Shen. Supported by the Institute Curie, Institut National de la Santé et de la Recherche Médicale, Centre National de la Recherche Scientifique, ANR "ChromaTin" grant ANR-10-BLAN-1326-03, and "EPICURE" grant ANR-14-CE16-0009 (S.A. and G.A.); la Ligue Contre le Cancer (Equipe labelisée Ligue, EL2014.LNCC/SA) and Association de RechercheContreleCancergrantERC(2013-AdGN°340046DCBIOX) (S.A.); la Ligue Nationale contre le Cancer (Equipe labelisée Ligue), the European Commission Network of Excellence EpiGeneSys (HEALTH-F4-2010-257082), ERC Advanced Grant 2009-AdG_20090506 "Eccentric, " the European Commission large-scale integrating project FP7_HEALTH-2010-259743 "MODHEP, " ANR-11-LABX-0044_DEEP, ANR-10-IDEX-0001-02 PSL, ANR "CHAPINHIB" ANR-12-BSV5-0022-02, ANR "CELLECTCHIP" ANR-14-CE10-0013, and the Aviesan-ITMO cancer project "Epigenomics of breast cancer" (G.A.); and SiRIC/INCa grant INCa-DGOS-4654 (J.J.W.). High-throughput sequencing was performed by the ICGex NGS platform of the Institut Curie supported by the grants SESAME (région Ile-de-France), ANR-10-EQPX-03, and ANR-10-INBS-09-08. The bioinformatics data are available in the Gene Expression Omnibus (GEO) database under accession numbers GSE105163 (microarrays) and SuperSeries GSE106268 for ChIP-seq and scRNA-seq data (linked to subseries GSE106264, GSE106265, and GSE106267). The authors have no conflicting financial interests. The data are tabulated in the main paper and in the supplementary materials. Author contributions: L.P. conceived and designed the project, carried out experimental work, and wrote the manuscript; C.G. contributed to bioinformatics design, carried out bioinformatics work, and interpreted data; E.Z. contributed to project design, carried out experimental work, and interpreted data; P.G. carried out bioinformatics work and participated in data interpretation; N.B. carried out experimental work; J.J.W. participated in bioinformatics data analysis and interpretation; J.-P.Q. participated in project design and data interpretation and wrote the manuscript; G.A. designed the project and wrote the manuscript; and S.A. conceived and designed the project and wrote the manuscript. L.P., C.G., and S.A. are inventors on patent application EP17305757.1 submitted by Institut Curie that covers the targeting of Suv39h1 in the context of T cell adoptive transfer.

Figure 2 .

 2 Figure 2. PD-1 blockade in Suv39h1-KO CD8+ TILs modifies their exhaustion program and enhances their effector capacity (A) Representative histogram and frequency (%) of inhibitory receptors (PD-1 + , TIM-3 + , LAG-3 + and 2B4 + ) on CD8 + TILs from B16F10-OVA tumors. (B) Pie chart of co-expression (PD-1 + , TIM-3 + , LAG-3 + and 2B44 + ) on CD8 + TILs from B16F10-OVA tumors. (C) Representative histogram and frequency (%) of GZMb + cells among CD8 + TILs. Cells were re-stimulated in vitro with PMA and ionomycin for 4 hr. (GZMb: Granzyme b). (D) Blood cells from mice were re-stimulated ex vivo with OVA-I SIINFEKL peptide. The number of OVA-specific T cells producing IFN-γ per 2,5x10 5 blood cells was determined by ELISPOT analysis. (E) Representative histogram and frequency (%) of ISG15 + among CD8 + TILs from B16F10-OVA tumors. (F) Representative histogram and frequency (%) of Ki67 + among CD8 + TILs from B16F10-OVA tumors. (G) Representative histogram and frequency (%) of BCL-2 + among CD8 + TILs from B16F10-OVA tumors. (H) Frequency (%) of Annexin V + 7ADD -among CD8 + TILs from B16F10-OVA tumors. (I) Representative contour plots and frequency (%) of TBET + EOMES + among CD8 + PD-1 + TILs from B16F10-OVA tumors. (J) Representative contour plots and frequency (%) of progenitor exhausted (TCF1 + TIM-3 -) and late exhausted (TCF1 - TIM-3 + ) among CD8 + PD-1 + TILs from B16F10-OVA tumors. (K) Representative histogram and frequency (%) of transitory (CD101 -) and exhausted (CD101 + ) among CD8 + PD-1 + TCF1 -TIM-3 + TILs from B16F10-OVA tumors. A representative experiment out of two is shown. p values were calculated using Mann-Whitney test. *p<0.05; **p<0.01; ***p<0.001; ****p<0.0001. See also Figures S3.

Figure

  Figure 6. Pharmacological inhibition of Suv39h1 potentiates tumor rejection by anti-PD-1 Ab

6 .

 6 Figure 6. Pharmacological inhibition of Suv39h1 potentiates tumor rejection by anti-PD-1 Ab (A) Graphical representation of model system of experimental groups, including C57BL/6 mice receiving B16F10-OVA melanoma cells followed by ETP-69 oral treatment, PBS or anti-PD-1 Ab injection. (B) Tumor growth kinetics represented as means of one representative experiment out of 2; with n=X o X mice per group. Black arrows indicate time of initial ETP-69 or vehicle control dose administration and anti-PD-1 Ab injection. (C) Tumor volumes in cm 3 (on the day of sacrifice, day 19). (D) Frequency (%) and quantification (number) of CD8 + TILs (CD45 + TCRb + CD4 -). (E) Frequency (%) of memory and effector CD8 + TILs from B16F10-OVA tumors. (F) Representative histogram and frequency (%) of inhibitory receptors (PD-1 + , TIM-3 + , LAG-3 + and 2B4 + ) on CD8 + TILs from B16F10-OVA tumors. (G) Pie chart of co-expression (PD-1 + , TIM-3 + , LAG-3 + and 2B44 + ) on CD8 + TILs from B16F10-OVA tumors. (H) Representative histogram and frequency (%) of GZMb + and IFNγ + cells among CD8 + TILs. Cells were re-stimulated in vitro with PMA and ionomycin for 4 hr. (I) Representative histogram and frequency (%) of Ki67 + among CD8 + TILs from B16F10-OVA tumors. (J) Frequency (%) of TBET + EOMES + among CD8 + PD-1 + TILs from B16F10-OVA tumors. A representative experiment out of two is shown. p values were calculated using Mann-Whitney test. *p<0.05; **p<0.01; ***p<0.001; ****p<0.0001. See also Figure S8. See also Figures S6.

Fig. 1 . 10 P

 110 Fig. 1. Characterization of CD3 + TILs in NSCLC. (A) Graphical overview of the study design. scRNA-and TCR-seq were applied to tumor tissues, normal adjacent lung tissues, and blood samples derived from 11 primary, untreated patients with NSCLC. Spectral flow cytometry analysis was applied in an additional dataset of patients with the same clinical characteristics. (B) UMAP of 28,936 single CD3 + tumor-infiltrating T cells from 11 patients with NSCLC, showing the formation of 21 main clusters, including 7 for CD8 + cells, 11 for CD4 + cells, 1 for T cells highly expressing IFN-related genes, and 2 for cycling T cells. (C) Summary of the distribution of the number of cells contributing to each cluster. (D) Heatmap of expression values for the top genes with enriched expression in CD3 + T cells, found by k-nearest neighbors' subclustering. Expression values are zero-centered and scaled for each gene. (E) Projection of a selected set of marker genes and gene signatures, identifying T cell state. Each cell is colored on the basis of the normalized expression.

Fig. 2 .

 2 Fig. 2. Stem-like precursors and late dysfunctional CD8 + cells. (A) Expression of single marker genes and gene signatures for memory-like, precursors CD8 + cells (left) and terminally differentiated CD8 + cells (right), respectively. (B) UMAP projection of tumor-infiltrating CD8 + cells, with each cell colored based on the gene signature score of a published study in melanoma patients. (Left) The expression of genes associated with therapeutic response and improved outcome to anti-PD-1 treatment (good response). (Right) The expression of genes associated with bad response to anti-PD-1 treatment. (C) Violin plots of the distribution of the expression scores of the two gene signatures among the different CD8 + clusters. (D) Violin plots of differentially expressed genes between CD8-KLF2 and CD8-GZMK clusters. (E) Volcano plot of differentially expressed genes between CD8-KLF2 and CD8-GZMK clusters. Each red dot denotes an individual gene passing P value and fold difference thresholds.

Fig. 3 .

 3 Fig. 3. Tissue-resident and transitional CD8 + subsets. (A) UMAP projection of tumor-infiltrating CD8 + cells, with each cell colored on the basis of the relative normalized expression of a gene signature for tissue residency, consisting of the following genes: ITGAE, ITGA1, XCL1, and ZNF683. (B) Flow cytometry plots of PD-1 and TIM3 expression in subsets of CD103 + CD8 + and CD103 - CD8 + cells. Frequency of CD103 + and CD103 cells in eight different subsets of CD8 + TILs, defined by the expression of three inhibitory ICP molecules: PD-1, TIM3, and CD39. Data representative of six independent experiments/ patients with NSCLC (**P ≤ 0.01, two-way ANOVA test). FSC-A, forward scatter-area. (C) Representative tSNE analysis of tumorinfiltrating CD8 + cells with each cell colored based on the median fluorescence intensity of each of the markers: CD103, PD-1, TIM3, CD39, KLRG1, and GZMB. (D) Summary of the frequency of eight different CD8 + cell subsets inside CD103 + and CD103 -CD8 + cells. Data representative of six independent experiments/ patients with NSCLC. (E) Partitionbased graph abstraction of CD8 + tumor-infiltrating cells, from 11 patients with NSCLC, connects circulating (GZMK and KLF2) and resident (XCL1) CD8 + precursors with terminally differentiated CD8 + (LAYN), through GZMH cluster, without any direct connection between effector CD8 + cluster (FCGR3A). (F) Trajectory of CD8 + tumor-infiltrating cells transition state, from 11 patients with NSCLC in a two-dimensional state-space determined by Monocle3. Each dot represents a single cell. Each color represents a different CD8 + cluster. (G) Visualization of the silhouette coefficient score on the UMAP of the integrated CD8 + cells reference from 11 NSCLC tumor samples. Silhouette coefficient is calculated on the basis of the mean intracluster distance and the mean of the nearest cluster distance for each cell of each cluster. Each cell is colored on the basis of the score: 1 is the highest value, highlighted in purple and indicating robust clusters; -1 is the lower value, highlighted in gray and indicating less robust clusters; 0 is indicative of overlapping clusters. (H) Label transfer of cell type labels from GZMH cluster onto the integrated CD8 + cells reference from 11 NSCLC tumor samples. Colored cells are cells from CD8-GZMH cluster and gray cells are coming from the other clusters.

Fig. 4 .

 4 Fig. 4. Clonal expansion and TCR sharing in clusters. (A) Visualization of the clonal size on a UMAP embedding on the integrated 11 tumors reference. (B) Top: Quantification of clonal size per cluster identities. Bottom: Expansion index computed for each cluster and each patient. (C) Selected examples of cluster-specific clones from P60. (D) Normalized heatmap of clonal sharing probabilities (as defined in the subsection TCR analysis of the methods). (E) Circos plot of clonal sharing between clusters. Each line represents one shared clone with the most frequent sharing. (F) Scheme of sharing similarities between CD8 clusters as computed in (D).

Fig. 6 .

 6 Fig. 6. Tissue-specific features-integrative model. (A) UMAP representation of the integrated cells from 17 samples from tumor, juxtatumor, and blood. (B) UMAP representation of cell densities split by the three tissues. (C) Scaled proportions of cells per clusters across the three tissues. (D) Quantification of clonal expansion for each cluster across the three tissues. (E) Circos plots of clonal sharing between tumor and blood cells. Each line represents one shared clone with the most frequent sharing. (F) Top 10 shared clones of blood (top) and tumor (bottom) being shared with tumor and blood, respectively, for each CD8 cluster. (G) Summary scheme of CD8 circulation and differentiation in NSCLC. Graph created using BioRender software (BioRender.com).

DOI: 10 .OPEN 1

 101 1038/s41467-018-04985-0 Institut Curie, PSL Research University, INSERM, U932, 26 rue d'Ulm, 75005 Paris, France. 2 Sanofi, Breakthrough Laboratory, 1 Impasse des Ateliers, 94400 Vitry-sur-Seine, France. 3 Institut Curie, PSL Research University, NGS Platform, 26 rue d'Ulm, 75005 Paris, France. Correspondence and requests for materials should be addressed to E.S. (email: elodie.segura@curie.fr) NATURE COMMUNICATIONS | (2018)9:2570 | DOI: 10.1038/s41467-018-04985-0 | www.nature.com/naturecommunications 1 1234567890():,;

  and 8 displayed the highest scores. Clusters 7 and 8 had the highest score for the in vitro mo-DC signature, with clusters 9 and 11 ARTICLE NATURE COMMUNICATIONS | DOI: 10.1038/s41467-018-04985-0 2 NATURE COMMUNICATIONS | (2018)9:2570 | DOI: 10.1038/s41467-018-04985-0 | www.nature.com/naturecommunications displaying low scores. Finally, cluster 8 had the highest score for the "activated DC" signature, with some cells from cluster 9 also displaying high scores, and cluster 7 showing an intermediate score for this signature.

Fig. 1

 1 Fig. 1 Ascites DCs are distinct from classical DCs. Purified ascites DCs, ascites macrophages, tonsil DCs, and blood monocytes were analyzed by single-cell RNA-seq using a Drop-seq approach. Combined single-cell transcriptomes were analyzed. a, b t-SNE representation of cell clusters identified using unsupervised clustering. Each dot represents an individual cell. a Colors represent sample origin; don donor. b Colors represent identified clusters. Clusters are manually ordered and cell numbers for each cluster is indicated. c Heatmap of scaled expression (log values of UMI) for the top 20 differentially expressed genes of each cluster (based on log fold change)

1 Fig. 2

 12 Fig. 2 Ascites DCs are monocyte-derived cells. Purified ascites DCs, ascites macrophages, tonsil DCs, and blood monocytes were analyzed by single-cell RNA-seq using a Drop-seq approach. Combined single-cell transcriptomes were analyzed. a Signature scores in individual cells for indicated gene signatures. b Annotation of cell clusters

Fig. 3

 3 Fig. 3 Human mo-DCs and mo-Mac both cross-present efficiently. Purified DCs and macrophages from tumor ascites (a) or in vitro culture of monocytes (b, c), or DCs derived in vitro from CD34 + precursors (d) were incubated with serial concentrations of MelanA long or short peptide (a, b, d) or MelanA-coated beads (c). After washing, antigen-specific CD8 + T cells were added. After 24 h, IFN-γ secretion was assessed as a measure of T cell activation. Background level was subtracted. Mean ± SEM of three (a, d), six (b), or five (c) independent experiments

  ). Collectively, these results indicate that only m o -D C s ,b u tn o tm o -M a c ,a r ea b l et op r o v i d et h ec ostimulatory signals necessary for efficient cytotoxic CD8 + T cell differentiation.

Fig. 4 6 NATURE

 46 Fig. 4 Human mo-DCs and mo-Mac are inefficient for the transfer of exogenous proteins into their cytosol. a, b Purified DCs and macrophages from tumor ascites, derived in vitro from monocytes, or DCs derived in vitro from CD34 + precursors were loaded with a cell-permeable FRET-sensitive substrate of β-lactamase, and incubated with or without exogenous β-lactamase. After 3 h, cleavage was measured by flow cytometry. a Representative results of six (tumor ascites), ten (in vitro monocyte-derived), or eight (in vitro CD34 + cell-derived) independent experiments. b Quantification of β-lactamase transfer. Symbols represent individual donors. N = 6 for tumor ascites, N = 10 for in vitro monocyte-derived cells, and N = 8 for CD34 + cell-derived cells. **p < 0.01, Wilcoxon non-parametric test. c Purified DCs or macrophages were incubated with β-lactamase coupled to Atto dye 633 at 4 or 37 °C during 3 h. Representative results of three independent experiments

Fig. 5

 5 Fig. 5 Human mo-DCs and mo-Mac use the vacuolar pathway for cross-presentation. Purified in vitro-generated mo-DCs (a, d), mo-Mac (b, e), or CD34 + cell-derived CD1a + DCs (c, f) were incubated with serial concentrations of MelanA long or short peptide, in the absence or presence of lactacystin (a-c)or pan-cathepsin inhibitor (d-f). After washing, antigen-specific CD8 + T cells were added. After 24 h, IFN-γ secretion was assessed as a measure of T cell activation. Background level was subtracted. Mean ± SEM of three (a), five (b, d, e), or three (c, f) independent experiments

Fig. 6

 6 Fig. 6 Human mo-DCs, but not mo-Mac, are efficient inducers of effector cytotoxic CD8 + T cells. Purified DCs and macrophage from tumor ascites were cultured with allogeneic CellTrace Violet-stained naïve CD8 + T cells for 7 days, in the absence (a, b) or presence (c, d) of naïve CD4 + T cells autologous to CD8 + T cells. Expression of Granzyme A, Perforin, and IFN-γ was assessed by intracellular flow cytometry. a, c Representative results of eight independent experiments. Gated on live CD8 + T cells. b Number of proliferating CD8 + T cells is shown. Symbols represent individual donors. N = 8. Median is shown. d Number of CD8 + T cells expressing effector molecules is shown. Symbols represent individual donors. N = 8. Median is shown. *p < 0.05, **p < 0.01, Wilcoxon non-parametric test

Fig. 7

 7 Fig. 7 Human mo-DCs, but not mo-Mac, provide co-stimulatory signals for the differentiation of cytotoxic CD8 + T cells. a Heatmap of scaled expression for selected co-stimulation genes. d donor. b Purified DCs and macrophage from tumor ascites were cultured with allogeneic naïve CD4 + T cells for 7 days, in the presence of naïve CD8 + T cells. Total number of live CD4 + T cells at the end of the culture is depicted. N = 8. c Purified DCs and macrophage from tumor ascites were cultured in the absence or presence of CD40-L, IFN-γ, and R848 for 24 h. Secretion of IL-12p70 was measured in the supernatant. N = 5. b, c Symbols represent individual donors. Median is shown. ** p < 0.01, Wilcoxon non-parametric test

  COMMUNICATIONS | (2018)9:2570 | DOI: 10.1038/s41467-018-04985-0 | www.nature.com/naturecommunications Thèse de doctorat Université de Paris Ecole doctorale BioSPC ED652 Institut Curie -U932 / Équipe Amigorena / Immunité et Cancer Tumor-infiltrating CD8+ T cell differentiation paths at the single-cell level Par Paul GUEGUEN Thèse de doctorat d'Immunologie Dirigée par Sebastian Amigorena Présentée et soutenue publiquement le 29 Mars 2021 Président du jury : Eric Tartour / PU-PH / Paris Centre de Recherche Cardiovasculaire -UMR S970 Rapporteurs : Fathia Mami-Chouaib / DR1 / Gustave Roussy -Integrative Tumor Immunology and Immunotherapy -UMR 1186 Julie Déchanet-Merville / DR1 / ImmunoConcEpT lab -UMR 5164 Examinateurs : Laura Cantini / CRCN / Institut de Biologie de l'École Normale Supérieure -UMR 8197 Mickael Menager/ CRCN-HDR / Institut Imagine -SCInflaNet lab Membre invité: Vassili Soumelis / PU-PH / Hôpital Saint-Louis -Immunologie humaine, pathophysiologie et immunothérapie -UMR 976 Directeur de thèse : Sebastian Amigorena / DRCE / Institut Curie -U932

Fig. 1 .

 1 Fig.1. The updated Hallmarks of cancer.

Fig. 2 .

 2 Fig.2. Tumor cells clonal evolution and the consequences of therapeutic pressure.

Fig. 3 .

 3 Fig.3. The three classes of Tumor Micro Environment.

Fig. 4 .

 4 Fig.4. Main differences between innate and adaptive immunity.

Fig. 5 .Fig. 6 .

 56 Fig.5. The diversity of cellular subtypes of the innate and adaptive immune system Fig.6. Immune checkpoint inhibition using Anti PD-1 monoclonal antibodies.

Fig. 7 .

 7 Fig.7. Curves of progression-free survival comparing chemotherapy and a-PD1 (Pembrolizumab) treated patients with NSCLC Fig.8. Process of CAR T therapy production.

Fig. 9 .

 9 Fig.9. Recognition of cognate tumor antigen by TCR recognition of MHC-bound tumor peptide.

Fig. 10 .Fig. 11 .

 1011 Fig.10. The two models for T cell activation Fig.11. Scheme representing the key feature identified in TRM cells.

Fig. 12 .

 12 Fig.12. Epigenetics implicated in two models of CD8+ T cell differentiation Fig.13. Single-cell technologies as fruit mixtures.

Fig. 14 .

 14 Fig.14. The maximum number of sequenced single-cells profiled per study and per month.

Fig. 15 .

 15 Fig.15. Scheme presenting different scRNA-seq library preparation techniques.

Fig. 16 .

 16 Fig.16. Comparison of different scRNA-seq library preparation techniques Fig.17. The range of multimodal scRNA-seq technologies.

  Fig.17. The range of multimodal scRNA-seq technologies.

Fig. 18 .

 18 Fig.18. GEM beads from 10X Genomics Fig.19. RNA capture process by GEMs.

Fig. 20 .

 20 Fig.20. The workflow of 10x Genomics Chromium library preparation process.

Fig. 21 .

 21 Fig.21. Scheme of V(D)J recombination.

Fig. 22 .

 22 Fig.22. scRNA-seq analysis pipeline. Starting from sequencer base calls.

Fig. 23 .Fig. 24 .Fig. 25 .Fig. 26 .

 23242526 Fig.23. CellRanger preprocessing pipeline Fig.24. Number of single-cell tools according to time up to January 2021 Fig.25. Repartition of single-cell tools per category as per January 2021 Fig.26. Cell-centric and gene-centric approaches used for scRNA-seq downstream analyses.

Fig. 27 .Fig. 28 .Chapter 6 :

 27286 Fig.27. V(D)J transcript counting pipeline. Fig.28. The different subtypes of lung carcinomas.

Fig. 1 .

 1 Fig.1. The updated Hallmarks of cancer. The 6 original hallmarks of cancer are shown (grey), with two new emerging hallmarks (cyan), with two tumor-enabling factors (violet). According to this model, most cancers have acquired the same functional capabilities during their development, albeit through various mechanistic strategies. (Adapted from Hanahan and Weinberg, 2011, using Biorender).

Fig. 3 .

 3 Fig.3. The three main classes of TIME. a, Infiltrated-excluded TIME is characterized by the exclusion of CTLs from the tumor core. b, Infiltrated-inflamed TIME is defined by an abundance of PD-L1 expression on tumor and myeloid cells and highly activated CTLs characterized by expression of Grzb, IFNγ, and PD-1. c, TLS-TIMEs contain TLSs, aggregates of immune cells that resemble lymph nodes, but. These TLS include B and T cells, dendritic cells, and Treg cells. Image is taken from Binnewies et al., 2018.

Fig. 4 .

 4 Fig.4. Scheme summarizing differences between the quick and non-specific innate immunity and the long-term and adaptive immunity, which relies on specific T and B receptors. Made with BioRender.

Fig. 6 .

 6 Fig.6. Immune checkpoint inhibition using Anti PD-1 monoclonal antibodies. aPD1 ICPI allows T cell activation by inhibiting the PD1/PDL1 axis binding between T cell and tumor cells. Made with Biorender.

Fig. 7 .

 7 Fig.7. Kaplan-Meier estimates of progression-free survival, comparing chemotherapy and a-PD1 (Pembrolizumab) treated patients with NSCLC. Tick marks represent data censored at the last time the patient was known to be alive and without disease progression. Figure from[START_REF] Reck | Pembrolizumab versus Chemotherapy for PD-L1-Positive Non-Small-Cell Lung Cancer[END_REF].

  and several research groups independently showed consistent objective response rates of 40-70% (Wu et al., 2012). CAR T cells Chimeric antigen receptor (CAR) T-cell therapy represents a breakthrough in personalized cancer treatment. In this strategy, T cells (usually coming from the same patient) are genetically engineered to express a chimeric receptor that binds a tumor antigen. CAR T cells are then expanded and infused back into the patient's body to recognize and destroy tumor cells. Dramatic clinical responses and high rates of complete remission have been observed in the setting of CAR T-cell therapy in B-cell malignancies. This resulted in two recent FDA approvals of CAR T cells directed against CD19 to treat acute lymphoblastic leukemia and diffuse large B-cell lymphoma (Feins et al., 2019).

Fig. 8 .

 8 Fig.8. Process of CAR T therapy production. T cells are first acquired from the patient, transduced with the CAR construct. CAR T cells are then expanded to hundreds of millions of cells and infused back into the patient. Made with Biorender

  Infiltration by T CD8+ and transcriptomic signatures The function of immune cells is often different when inside the tumor microenvironment compared to when they circulate in the peripheral blood (Schaer et al., 2011). Therefore, the presence of various subtypes of immune cells within the tumor microenvironment could more accurately predict response to checkpoint inhibitors. Deep infiltration of the tumors by CD8+ T cells often correlates with better response to treatment and patient outcome (Geng et al., 2015). It has also been reported that CD8+ T cell infiltration of the invasive tumor margin predicts MM patient response to anti-PD-1 treatment (Tumeh et al., 2014).

Fig. 9 .

 9 Fig.9. Recognition of cognate tumor antigen by T cells. T cells recognize MHC-bound tumor peptides through the TCR. These peptides are presented by APCs, like dendritic cells or macrophages.

Fig. 10 .

 10 Fig.10. Model for T cell activation following antigen recognition. a, T cell differentiation in response to acute infection or vaccination (left) or chronic infection or cancer (right). b, Main T cell fates and phenotypes arising following T cell activation. Taken from (McLane et al., 2019).

Fig. 11 .

 11 Fig.11. Scheme representing the key features identified in human TRM cells over the years. Upregulated genes are shown in green, and downregulated genes are shown in red. Image is taken from Kumar et al., 2017.

Fig. 12 .

 12 Fig.12. In the developmental, or linear, differentiation model, the progressive acquisition of effector function during CD8+ T cell differentiation (Naive→TSCM→TCM→TEM→TEFF) depends on the strength and duration of antigenic signaling. It results in the progressive loss of memory-associated gene expression and gain of effector-associated gene expression. These transcriptional changes are accompanied by similar changes in the epigenetic landscape, which are illustrated by the gradual gain or loss of activating and repressive histone modifications. The figure is taken from (Henning et al., 2018).

Fig. 13 .

 13 Fig.13. Single-cell technologies as fruit mixtures. If single-cells are fruits, bulk genomics (preceded by flow cytometry enrichment or not) only allows us to capture a fruit smoothie. In contrast, single-cell technologies allow us to get to the single-fruit level. Image adapted from Shalek & Regev.

Fig. 14 .

 14 Fig.14. The maximum number of sequenced single-cells profiled per study and per month. The data

Fig. 18 .

 18 Fig.18. Detailed view of the surface of GEM beads. These GEMs are coated with cell barcodes (in green), Unique Molecular Identifiers (UMIs, in red), and a Poly(dT) sequence (in blue) as a way to capture the polyA tail of mRNAs. (Adapted from 10x Genomics). mRNA polyA tail is then captured thanks to the 30bp complementary poly(dT) sequence.

Fig. 19 .

 19 Fig.19. RNA capture process by gel beads. The 30bp polyA tail will bind to the poly dT sequence from the mRNA. The Reverse Transcription (RT) step will allow the transcript to be added to the barcoded sequence. (Adapted from 10x Genomics).

Fig. 20 .

 20 Fig.20. The workflow of 10x Genomics Chromium library preparation process. We have single cells in suspension and 10x Gel Beads and reagents as an input. Cells are then encapsulated in oil droplets. mRNA is then captured on barcoded gel beads. mRNA is then retro-transcribed into cDNA and amplified using PCR. Finally, barcoded cDNA are pooled together. As the output, the digital gene expression profiles from every partitioned cell are retrieved. (Adapted from 10x Genomics).

Fig. 21 .

 21 Fig.21. Scheme of V(D)J recombination. Image adapted from 10X Genomics.

Fig. 22 .

 22 Fig.22. scRNA-seq analysis pipeline. Starting from sequencer base calls. The first step is to acquire the expression matrix starting from either sequencer images or FASTQ files. In the second step, the data are cleaned and normalized. Finally, in downstream analysis, we usually look to assign names and functions to cells and groups of cells and explore gene dynamics-the figure is taken from(Zappia et al., 2018).

Fig. 23 .

 23 Fig.23. CellRanger pipeline. First barcodes (cell barcodes and UMIs) are demultiplexed. Then, STAR is used to perform alignment on the reference genome. Transcripts (UMIs) are then counted. This will lead to the generation of a gene-cell expression matrix that can be used for downstream analysis. The sample aggregation pipeline can be used when we look to merge different samples. This procedure will normalize the different samples according to their number of reads per cell. Figure Adapted from 10x Genomics.

  . It uses an unsupervised deep generative model to learn the prior distribution of gene expression profiles and estimate the background RNA profile (Fleming et al., 2019). More recently, EmptyNN (Yan et al., 2021) implements a machine learning bagging strategy, also based on the rationale that barcodes with very low total counts represent cell-free droplets. Doublet removal scRNA-seq data interpretation is confounded by technical artifacts known as doublets/multiplets, meaning droplets that have encapsulated more than one cell. They typically appear due to errors in cell sorting or capture errors, especially in droplet-based protocols (Zheng et al. 2017) involving thousands of cells. Doublets are undesirable when the aim is to characterize populations at the single-cell level. Particularly, doublets can be mistaken for intermediate populations or transitory states that do not actually exist and can confound analyses like trajectory inference. It is then desirable to identify and remove doublet libraries not to compromise the interpretation of the results. Tools like DoubletFinder (McGinnis et al., 2019), Solo (Bernstein et al., 2020), DoubletDecon (DePasquale et al., 2019) were developed in that regard. Benchmarking efforts tend to show no perfect method exist for now

  methods such as MAGIC (van Dijk et al., 2018), SAVER (Huang et al., 2018), SAVERX (Wang et al., 2019), ENHANCE (Wagner et al., 2019). However, it has been shown that imputation introduces circularity that can generate false-positive results. While they can be useful for visualization, statistical tests applied to imputed data should be treated carefully (Andrews and Hemberg, 2019b). Normalization Normalization refers to the within-sample removal of technical variability from a scRNA-seq dataset while maintaining biological variability. Technical sources of variation are stochastic loss, amplification/PCR bias, mainly due to the minimal quantity of starting material in each droplet.

  Several methodologies have been developed and implemented. Among approaches, Principal Component Analysis (PCA), t-SNE (Maaten and Hinton, 2008), Uniform Manifold Approximation and Projection (UMAP) (McInnes et al., 2020) have been used. Some methods were explicitly designed for single-cell data, like PHATE (Moon et al., 2019) or ScPhere (Ding and Regev, 2019),

  methods have recently emerged to enable automated cell type identification. Many recent tools were developed specifically for single-cell research, including scmap (Kiselev et al., 2018), SingleR (Aran et al., 2019), CHETAH (de Kanter et al., 2019), Garnett (Pliner et al., 2019) and

  al., 2020) has been developed. Cell/cell interactions Multicellular life relies on the coordination of cellular activities, which depend on cell-cell interactions (CCIs) across an organism's diverse cell types and tissues. Thus, studies on cellular functions require consideration of the relationships between each cell type. CCIs leverage diverse molecules, including ions, metabolites, integrins, receptors, junction proteins, structural proteins, ligands and secreted proteins of the extracellular matrix. Some molecules support structural CCIs (for example, cell adhesion molecules), whereas ligands such as hormones, growth factors, chemokines, cytokines and neurotransmitters mediate cell-cell communication (CCC). scRNA-seq can be used to predict in silico ligand-receptor pairs that might be used in intercellular communication, as well as autocrine communication. Approaches like CellphoneDB (Efremova et al., 2020), NicheNET (Browaeys et al., 2020), SingleCellSignalR (Cabello-Aguilar et al., 2020), ICELLNET (Noël et al., 2020) have been developed, specially designed for single-cell datasets. These tools can predict most interacting cell types and prioritize interacting ligands and receptors between senders and receivers cells. RNA Velocity RNA velocity-the time derivative of the gene expression state-can be estimated by distinguishing between unspliced and spliced mRNAs in common single-cell RNA sequencing protocols (La Manno et al., 2018). In practice, RNA velocity is a high-dimensional vector that predicts the future state of individual cells. Recently, CellRank combined velocity and trajectory modeling (Lange et al., 2020). This approach combines the robustness of trajectory inference with directional information from RNA velocity. More recently, Bergen et al. relaxed the steady-state assumption and considered the full dynamical model, thereby enabling applying the RNA velocity framework to a broader set of biological systems and states (Bergen et al., 2020). The dynamical model is implemented in their scVelo Python package. Preprocessing choices and the choice of velocity method have also been shown to affect the downstream results (Soneson et al., 2021), directly affecting biological interpretation. This tends to indicate that given the relative immaturity of the RNA velocity analyses, 2D velocity vectors shown on UMAPs and tSNEs should still be interpreted with caution. Transcript usage -isoforms Alternative splicing shapes the phenotype of cells in development and disease and have been recently probed using single-cell technologies (Arzalluz-Luque and Conesa, 2018). Different types of splicing variations can occur at different points of the mRNA molecule. Partial sequencing of the transcript, like 3' or 5' mRNA capture, will naturally overlook a fraction of the events. Still, the tool SIERRA (Patrick et al., 2020) was developed to detect differential transcript usage from data generated by commonly used polyA-captured scRNA-seq technology like 10X Genomics.Recent protocols, like Smartseq3 (Hagemann-Jensen et al., 2020), can achieve full-length transcript coverage but are still unable to assemble the complete transcript sequence and are heavily reliant on the reference annotation. Long-read sequencing can overcome this limitation and generate full-length transcript information in single cells, as illustrated in several recent studies.[START_REF] Lebrigand | High throughput error corrected Nanopore single cell transcriptome sequencing[END_REF][START_REF] Tian | Comprehensive characterization of single cell full-length isoforms in human and mouse with long-read sequencing[END_REF]. These considerations also hold true for the single-cell analysis of Transposable Elements (TEs), which are insertional mutagens that contribute to eukaryotic genomes' plasticity. Some studies demonstrated that TE expression dynamics could be identified across single preimplantation embryonic cells and cancer cell lines[START_REF] Brocks | Single cell analysis reveals dynamics of transposable element transcription following epigenetic de-repression[END_REF][START_REF] Göke | Dynamic transcription of distinct classes of endogenous retroviral elements marks specific populations of early human embryonic cells[END_REF].Machine learning applied to single-cell data Following its remarkable achievements in image recognition and natural language processing, machine learning is now effectively drawing the attention of biotechnology. In single-cell datasets, machine learning has successfully been used in almost every step of the pipeline previously described. Importantly scVI(Lopez et al., 2018) is a tool that uses deep neural networks for batch correction, visualization, clustering, and differential expression. Extensions to scVI include scANVI for cell annotation of scRNA-seq data using semi-labeled examples (Xu et al., 2020), totalVI for analysis of CITE-seq (Cellular Indexing of Transcriptomes and Epitopes by sequencing) data (Gayoso et al., 2020), among many other extensions. As single-cell datasets continue to grow in size and scope, it is expected that machine learning tools will be of growing importance and could further bridge the gap between biologists and artificial intelligence researchers. Integrating modalities Multimodal single-cell technologies, which simultaneously profile multiple data types in the same cell, represent a new frontier for cell states' discovery and characterization. It is now possible to Seurat v4 uses a 'weighted-nearest neighbor' (WNN) framework to integrate multiple data types measured within a cell and to obtain a joint definition of cellular state. (Hao et al., 2020). They demonstrate better resolution of cell types integrating RNA / protein and RNA/ATAC. The tool TESSA and CoNGA can combine scRNA and TCR modalities (Schattgen et al., 2020; Zhang et al., 2021). It identifies correlations between GEX profile and TCR sequence through statistical analysis of a pair of T cell similarity graphs. Applying CoNGA and TESSA across diverse human and mouse T cell datasets uncovered known and novel associations between TCR sequences and cellular phenotype.

Fig. 27 .

 27 Fig.27. V(D)J transcript counting pipeline. Image taken from 10X Genomics.

  stem cell/memory gene expression program is under the control of Suv39h1 by imposing the H3K9me3 modification on chromatin at the corresponding loci. In doing so, Suv39h1/H3K9me3 would establish an epigenetic barrier on the stem/memory gene expression program, preventing effector reprogramming into memory cells. It is most likely that the possibly reversible silencing of effector gene expression in memory cells occurs through other mechanisms, as memory cells do effectively reprogram into effectors upon rechallenge. These results open new perspectives for the manipulation of epigenetic programming of T lymphocyte identity in the context of vaccination, checkpoint-based immunotherapies, and adoptive T cell therapies. A greater understanding of the CD8+ T cell epigenome may therefore provide essential clues for how to unlock the potential of highly differentiated, tumor-antigen-specific T cells infiltrating tumors. In the context of chronic TCR activation, cells have been shown to express an exhaustion phenotype (Wherry and Kurachi, 2015), and being hyporesponsive. CD8+ T cells in Suv39h1-KO mice are enriched for genes associated with naïve and memory signatures and showed enhanced memory potential and increased survival capacity. It is possible that epigenetic modifying drugs may reverse the repression of stem and memory genes in differentiated T cells and improve T cell-based immunotherapies. Thus, it still remains to know what would happen in the context of tumor immunology and immunotherapy. Introduction In the absence of Suv39h1 expression, CD8+ T cell effectors co-express memory genes, and display higher reconstitution potential after adoptive transfer (Pace et al., 2018). These results show that Suv39h1-induced heterochromatin dynamics plays a critical role during lymphocyte differentiation, and that chromatin modifiers control specific gene expression programs at precise stages of immune responses. Other epigenetic regulators, including EZH2 (Gray et al., 2017), DNA methyl-transferases (Youngblood et al., 2017) are also required to silence stem cell/memory programing during effector differentiation, indicating a general requirement for coordinated changes in chromatin dynamics. Several recent studies have investigated chromatin accessibility changes during CD8+ T cell differentiation into memory, effector or exhausted phenotypes. For the latter, transcriptomic and epigenomic changes that occur while T cells become exhausted in the context of chronic viral infection and cancer have been analyzed. A critical question in the field is to understand how blockade of checkpoints, such as PD-1, reverses the dysfunction of exhausted CD8+ T cells. Recent studies have proposed that reversion of the exhausted phenotype can occur early during the onset of the exhaustion program, and that once this program is fully established it becomes irreversible, and that fully exhausted T cells can no longer respond to anti-PD-1 (Philip et al., 2017). Understanding the molecular basis and underlying chromatin dynamics for the establishment of the exhaustion program is a cornerstone with direct consequences for the development of therapeutic strategies in cancer and chronic infections. Paper Conclusion Our results identify a critical epigenetic regulator of exhaustion. They also establish a clear link between IFN-I, early activation, cytolytic effectors and re-programing of progenitor exhausted cells. This link may be particularly relevant to the multiple immunotherapy approaches based on inducing IFN-I responses in the tumor environment, including cytolytic viruses, STING agonists, and DNA-demethylating agents. Since Suv39h1 genetic deficiency or pharmacological inhibition both re-activate anti-tumor immune responses, and promote re-programing by anti-PD-1, the role of Suv39h1 in this process must be at least partially non-redundant. By analogy to other non-redundant immunosuppressive proteins whose inhibition unleashes anti-tumor immune responses (often referred to as "immune checkpoints"), Suv39h1 can be considered as an "epigenetic immune checkpoint". Should the role of Suv39h1 be conserved in human T cells, its blockade would open new perspectives for the epigenetic manipulation of anti-tumor T cell responses in the clinic. Mice models represent attractive tumor models because they are low cost, have a short reproductive cycle, exhibit high tumor growth rates, and can be easily genetically modified. However, the obvious problem of these models is the high failure rate observed in human clinical trials after promising results obtained in mouse models. Efforts are ongoing to confirm these results in human tumors. Introduction Human primary samples are the closest we have from what happens in real patients. Lung cancer is the leading cause of cancer-related death worldwide. NSCLC represents 85% of lung cancer diagnoses, and tumor infiltration by lymphocytes is associated with a favorable survival prognosis and a better clinical response to immune checkpoint blockade (ICB).

Fig. 28 .

 28 Fig.28. The different subtypes of lung carcinomas. Lung cancers are broadly segregated into two subtypes : Small cell lung cancer (SCLC) and NSCLC. Our study was composed of 11 NSCLC patients, with nine adenocarcinomas, one squamous cell carcinoma and one undetermined tumor).

( 1 :

 1 50, clone 2ST8SH7, #9259590, BD Biosciences), PD-1-BV421 (1:40, clone EH12.2H7, #B268454, Biolegend), KLRG1-PE-Cy7 (1:40, clone SA231A2,#B256810, Biolegend), TIM3-BV786 (1:40, clone 7D3, #9259393, BD Biosciences), CD39-APC (1:50, clone eBioA1, #2071264, Thermofischer Scientific), PD-L1-APC-R700 (1:40, clone M1H1, #9010854, BD Biosciences), VCAM1-BB515 (1:20, clone S1-10C9, #9093693, BD Biosciences), CD103-PerCP-Cy5.5 (1:50, clone Ber-ACT8, #B246047, Biolegend). Intracellular markers staining was performed by using a combination of monoclonal antibodies from the following: GZMB-PE (1:20, clone QA16A02, #372207, BD Biosciences), Ki67-BV711(1:20, clone Ki67, #B291484, Biolegend).

FigFigure S3 .

 S3 Fig.S1.

Fig

  Fig.S4.

  In the absence of Suv39h1 expression, CD8+ T cell effectors co-express memory genes, and display higher reconstitution potential after adoptive transfer(Pace et al., 2018). These results show that Suv39h1-induced heterochromatin dynamics plays a critical role during lymphocyte differentiation, and that chromatin modifiers control specific gene expression programs at precise stages of immune responses. Other epigenetic regulators, including EZH2(Gray et al., 2017), DNA methyl-transferases(Youngblood et al., 2017) and histone chaperone CAF-1[START_REF] Ng | The histone chaperone CAF-1 cooperates with the DNA methyltransferases to maintain Cd4 silencing in cytotoxic T cells[END_REF] are also required to silence stem cell/memory programing during effector differentiation, indicating a general requirement for coordinated changes in chromatin dynamics.Several recent studies have investigated chromatin accessibility changes during CD8+ T cell differentiation into memory, effector or exhausted phenotypes[START_REF] Pauken | Epigenetic stability of exhausted T cells limits durability of reinvigoration by PD-1 blockade[END_REF][START_REF] Sen | The epigenetic landscape of T cell exhaustion[END_REF] Philip et al., 2017). For the latter, transcriptomic and epigenomic changes that occur while T cells become exhausted in the context of chronic viral infection and cancer have been analyzed. A critical question in the field is to understand how blockade of checkpoints, such

expression, Th2 cells differentiate, but fail to fully commit to the Th2 fate, and start producing IFNγ i n vi vo. In CD8+ T cells, Suv39h1 silences the memory/stem cell gene expression program during differentiation of cytotoxic effectors. as anti-PD-1 treatment, reverses the dysfunction of exhausted CD8+ T cells. Recent studies have proposed that reversion of the exhausted phenotype can occur early during the onset of the exhaustion program, and that once this program is fully established it becomes irreversible, and that fully exhausted T cells can no longer respond to PD-1 blockade

(Philip et al., 2017)

. Understanding the molecular basis and underlying chromatin dynamics for the establishment of the exhaustion program is a cornerstone with direct consequences for the development of therapeutic strategies in cancer and chronic infections.

RESEARCH | RESEARCH ARTICLE on January 21, 2021 http://science.sciencemag.org/ Downloaded from Pace et al., Science 359, 177-186 (2018)

of 10 RESEARCH | RESEARCH ARTICLE on January 21, 2021 http://science.sciencemag.org/ Downloaded from

at INSERM on January 29, 2021 http://immunology.sciencemag.org/ Downloaded from Gueguen et al., Sci. Immunol. 6, eabd5778 (2021) 29 January 2021

at INSERM on January 29, 2021 http://immunology.sciencemag.org/ Downloaded from

NATURE COMMUNICATIONS | (2018)9:2570 | DOI: 10.1038/s41467-018-04985-0 | www.nature.com/naturecommunications

Paul Gueguen -Doctoral thesis -2021

ACKNOWLEDGMENTS

We thank all the U932 and UMR3664 members; C. Maison, N. Lacoste, and F. X. Gobert for helpful discussions; A. Cros and V. Mondin for technical help; B. Amati and J. Sedat for critical reading of

Acknowledgments

We thank Virginie Dangles-Marie, Celine Daviaud, Isabelle Grandjean, Mikael Garcial and Cedrik Paurchad, the mouse facility technicians and flow cytometry core at Institut Curie. We thank Nathalie Amzallag for technical assistance. This work was supported by the ANR-14-CE16-009 (Programme ANR Innovation biomedical DS0404, 2014), Institut National de la Santé et de la Recherche Médicale, Labex DCBIOL (ANR-10-IDEX-0001-02 PSL and ANR-11-LABX0043), SIRIC INCa-DGOS-Inserm_12554. P.G. is supported by a PhD fellowship from Ligue Nationale contre le Cancer at Université Paris Descartes (Paris V). P.C. was Acknowledgments: We wish to thank the flow cytometry platform of Institut Curie and L. Guyonnet for helpful discussions. We thank W. Richer for fruitful discussions on label transfer and TCR analyses. We thank all the thoracic surgery, respiratory medicine, pathology, and tissue procurement staff of Institut Mutualiste Montsouris, as well as all the patients and their families. Funding: C.M. received a PhD fellowship from Servier Laboratories and funds by the Acknowledgments This work was supported by INSERM, Agence Nationale de la Recherche (ANR-10-LABX-0043, ANR-CHIN-0002, and ANR-10-IDEX-0001-02 PSL), Institut Curie (CIC IGR-Curie 1428 and ICGex-2017), the European Research Council (2013-AdG No. 340046 DCBIOX), La Ligue Nationale contre le Cancer (EL2014.LNCC/SA). T.-L.T.-H. is supported by an industry-funded PhD scholarship (CIFRE) at Université Paris Descartes (Paris V). P.G. is supported by a PhD fellowship from Ligue Nationale contre le Cancer at Université Paris Descartes (Paris V). High-throughput sequencing performed by the ICGex NGS platform of the Institut Curie was supported by grants ANR-10-EQPX-03 (Equipex) and ANR-10-INBS-09-08 (France Génomique Consortium), by Cancéropôle Ile-de-France, by SESAME program from Région Ile-de-France, and by SiRIC-Curie program (SiRIC Grant "INCa-DGOS-4654"). The authors wish to thank the Flow Cytometry Platform of Institut Curie for cell sorting.

Acknowledgments

PhD program "FIRE -Programme Bettencourt". P.G. is supported by a PhD fellowship from Ligue Nationale contre le Cancer at Université Paris Descartes (Paris V). This study was supported by grants from "Janssen Horizon," "Fondation ARC pour la recherche sur le cancer," "Fondation chercher et trouver," "European Research Council COG 724208" (T.D., (ANR-10-IDEX-0001-02 PSL and ANR-11-LABX-0043), Center of Clinical Investigation (CIC IGR-Curie 1428), Equipex (ANR-10-EQPX-03), France Génomique Consortium (ANR-10-INBS-09-08) from the Agence Nationale de la Recherche ("Investissements d'Avenir" program), and the SiRIC-Curie program (SiRIC grant INCa-DGOS-4654) (J.J.W.).

The data for this study are being deposited in the Gene Expression Omnibus database (GSE162500). Code used to analyze the data is available at https://github.com/p-gueguen/ single_cell_NSCLC.

Data availability. Single-cell RNA-seq data that support the findings of this study have been deposited in GEO with accession codes: GSE115007 for ascites DCs and ascites macrophages and GSE115006 for tonsil DCs.

Author contributions

L.L.N., C.S., L.P., S.A., and E.P. conceptualization. L.L.N., J.C., C.S., S.A., and E.P. designed experiments. L.L.N., E. Z., A.T., P.C., R.R., J.D., Q.Z., Z.L., B.K., and X.H. performed the experiments. M.Y., E. Z., L.L.N. analyzed ATACseq data. L.L.N., P.G., C.G., and J.W. analyzed scRNAseq data. J.-M.L., A.S., B.B., P.B., performed ETP-69 compound synthesis and pharmacokinetics studies. L.L.N., E. Z., S.A., and E.P. wrote the manuscript, with all authors contributing to writing and providing feedback. S.A., and E.P. supervision. J.-P. Q, G.A., S.A., and E.P. funding acquisition.

Declaration of conflict of interest

S. Amigorena, J. Waterfall and E. Piaggio are shareholders in Mnemo Therapeutics. P. Broqua and Anne Soude are employees and shareholders in Inventiva Pharma. The other authors declare no competing interests. matrix to preform preprocessing and UMAP reduction. DDRTree algorithm was then used to reconstruct the tree embedding. PAGA output was generated using Scanpy v. 1.4.3 (40) with default values and a threshold of 0.15.

TCR analysis

TCR-seq data for each sample were processed using Cell Ranger software (versions as above), with the command "cellranger vdj" using the human reference genome GRCh38.

Statistical analysis

All flow cytometry experiments were analyzed using Prism 8 (GraphPad software, version 8.4.2). All data were presented as means ± SEM. Statistical differences were assessed using two-way analysis of variance (ANOVA) test or two-tailed paired Student's t test.

SUPPLEMENTARY MATERIALS

immunology.sciencemag.org/cgi/content/full/6/55/eabd5778/DC1 Materials and Methods Fig. S1. Comparison of the 3′ versus 5′ single-cell chemistries and general characterization of total CD3 + T cell subsets. Fig. S2. Gene signatures and biological pathways differentiating precursors and terminally differentiated CD8 + T cells. Fig. S3. Flow cytometry gating strategy and quantification of silhouette score and maximum prediction score of the label transfer method in CD8 + T cell subsets. Fig. S4. Quantification of shared clones among the CD8 + clusters and transition index. Fig. S5. Characterization of the cycling CD8 + T cells. Fig. S6. Analysis of tissue specificity of tumor, juxtatumor, and blood tissue in patients with NSCLC. Table S1. Clinical data. Table S2. Differentially expressed genes per cluster. Table S3. Gene signatures. Table S4. TCR information. Table S5. Raw data table.

View/request a protocol for this paper from Bio-protocol.

Author contributions

Additional information

Supplementary Information accompanies this paper at https://doi.org/10.1038/s41467-018-04985-0.

Competing interests: The authors declare no competing interests.

Reprints and permission information is available online at http://npg.nature.com/ reprintsandpermissions/ Publisher's note: Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Preface

This preface provides a summary of the chapters in this thesis and describes my contribution to them. This thesis is with publications, where each publication is corresponding to a chapter. The following publications are included as part of this thesis:

1. Pace, L., Goudot, C., Zueva, E., Gueguen, P., Burgdorf, N., Waterfall, J.J., Quivy, J.-P., Almouzni, G., and Amigorena, S. (2018). The epigenetic control of stemness in CD8+ T cell fate commitment. Science 359, 177-186.

2. Niborski, L.L., Gueguen, P., Ye, M., Thiolat, A., Ramos, R., Caudana, P., Denizeau, J., These publications are as they were published or submitted and should be read as independent papers. References are available at the end of each publication rather than in the reference list for this thesis. Supplementary figures and methods are presented in annexes. The contributions of authors to these papers are explained below.

Chapter 1: Introduction provides background and motivation for this work. This chapter briefly summarizes the current immunotherapy landscape and our current need to understand basic immunological systems at the single-cell level to inform therapies. Contributions to the work in this chapter:

• I carried out bioinformatics work regarding the scRNA-seq data: data processing and analysis, with input from Christel Goudot, Josh Waterfall, Sebastian Amigorena, and Luigia Pace.

• Luigia Pace conceived and designed the project, carried out experimental work, and wrote the manuscript.

Chapter 4: CD8+ differentiation in mice in response to immunotherapy describes the impact of Suv39h1 in the context of immunotherapy. This work is currently under revision.

Contributions to the work in this chapter:

• I carried out bioinformatics work regarding the scRNA-seq data: data processing and analysis, with input from Josh Waterfall, Eliane Piaggio, Sebastian Amigorena, and Leticia Niborski.

• I created the figures shown in this paper with input from Josh Waterfall, Eliane Piaggio, Sebastian Amigorena, and Leticia Niborski.

• Leticia Niborski carried out experimental work.

Chapter 5: CD8+ differentiation in human NSCLC patients describes a paired scRNA-seq / scTCR-seq analysis of tumor, juxta tumor, and blood, where we study circulation and filiation of CD8+ T cells in patients with Non-Small Cell Lung Cancer (NSCLC). This work has been published in Science Immunology.

Contributions to the work in this chapter:

• I carried out bioinformatics work regarding the scRNA-seq data: data processing and analysis, with input from Josh Waterfall and Sebastian Amigorena.

• I performed the scRNA-seq and part of the scTCR-seq analysis.

• I designed and created the figures shown in this paper with input from Sebastian Amigorena and Christina Metoikidou.

• Thomas Dupic performed scTCR analyses, namely regarding cluster TCR sharing.

• Myriam Lawand and Christina Metoikidou carried out experimental work.

subjected to various pressures, including hypoxia and metabolic and nutritional limitations; consequently, some clones containing certain mutations may be more adept at coping with this pressure and grow in clonal size. Similarly, when challenged with therapy, some clones are eradicated, whereas others have mutations that confer a clonal advantage in the context of therapeutic pressure.

Such ecology of tumor clones is usually represented on so-called fish plots, where tumor clonal size can be appreciated according to time. We can appreciate how tumor clonal size and heterogeneity changes as different lines of treatments are applied, and eventually, as resistance occurs (Fig. 2).

Fig. 2. Illustrated here is an often-used schematic, known as a fish plot, which is used to describe tumor evolution and the consequences of therapeutic pressure. This evolutionary trajectory in a cancer cell population can be interpreted as genetic evolution due to selection pressure. Adapted from Marine et al., 2020.

Cancer immune infiltration

Technological and conceptual advances recently highlighted how tumors are complex mixtures of malignant, immune, and stromal cells, among many other cell types, which often have substantial levels of intratumor and interpatient heterogeneity. The tumor microenvironment (TME) comprises an amalgam of tumor-promoting and anti-tumor cells, producing or reacting to signals that can modulate tumor growth and influence tumor evolution.

ranges from the first study, which profiled 5 cells, to the latest and biggest study, published in November 2020 with over 4 million cells. Data was collected from www.nxn.se/single-cell-studies.

A bit of history: Early single-cell capture technologies

The first scRNA-seq protocol was published in 2009 [START_REF] Tang | mRNA-Seq whole-transcriptome analysis of a single cell[END_REF]. While this original approach allowed the transcriptome measurements in five individual cells, it required manual manipulation and was limited to probe few cells. Further studies showed that cell types could be identified with this technique, and approaches were developed to allow unbiased capture of the transcriptome. Since then, many plate-based scRNA-seq protocols have been developed When comparing different library generation techniques (Fig. 16), the authors found that Quartz-seq2 is the most sensitive method as it detects the most genes per cell and allowed better detection of cell-type-specific markers [START_REF] Mereu | Benchmarking single-cell RNA-sequencing protocols for cell atlas projects[END_REF]. The choice of the method should consider cost-efficiency, if polyA mRNA capture is enough or if we want full-length RNA capture, the number of cells desired, and sequencing depth needed according to the biological system analyzed.

CRISPR perturbations).

A more recent approach allowed simultaneous single-cell measurements of intranuclear proteins and gene expression, measuring multiplexed intranuclear protein levels and the transcriptome in parallel in thousands of cells, enabling joint analysis of TF levels and gene expression in vivo. [START_REF] Chung | Simultaneous single cell measurements of intranuclear proteins and gene expression[END_REF].

Protocols like

Raw data processing: Cellranger preprocessing workflow

After the sequencing of a single-cell experiment is done, the result of a sequencing experiment is usually a set of image files that came from the sequencer. For most analyses, we are using the expression matrix where each column is a cell, each row is a gene, and the values indicate expression levels. Before working on this matrix, we need to perform a series of pre-processing steps, typically beginning with some quality control of the raw reads. Reads are then aligned to a reference genome, and the number of reads overlapping annotated features (genes or transcripts) is counted.

Calling cell Barcodes

Cell Ranger 3.0 introduces an improved cell-calling algorithm that is better to identify populations of low RNA content cells, especially when low RNA content cells are mixed into a population of high RNA content cells. For example, tumor samples often contain large tumor cells mixed with smaller tumor-infiltrating lymphocytes (TIL), and researchers may be particularly interested in the TIL population. The new algorithm is based on the EmptyDrops method [START_REF] Lun | EmptyDrops: distinguishing cells from empty droplets in droplet-based single-cell RNA sequencing data[END_REF].

The algorithm has two key steps:

-It uses a cutoff based on the total UMI counts of each barcode to identify cells. This step identifies the primary mode of high RNA content cells.

-Then the algorithm uses the RNA profile of each remaining barcode to determine if it is an "empty" or a cell containing partition. This second step captures low RNA content cells whose total UMI counts may be similar to empty GEMs.

Cleaning the data: Pre-processing and quality control

Cell filtering

Quality control of individual cells is important as experiments that contain low-quality cells will eventually lead to misleading results. It may be cells that were stressed or lysed during the experimental handling and processing, or it could be that droplets that didn't encapsulate any cells, or even droplets that encapsulated more than one cell (referred to as doublets or multiplets).

Empty / bad quality doublet removal

The first step is to select which droplet or well did encapsulate a cell or not. Since ambient RNA in the cell-free droplets may be reverse transcribed into cDNA during the library preparation, which will produce unique molecular identifier (UMI) counts in the resulting gene expression matrices. Therefore, cell-free droplets are difficult to distinguish from cell-containing droplets.

Failure to remove cell-free droplets may introduce spurious biological signals in the downstream analysis (Lun et In practice, PCA and t-SNE/UMAP are usually combined. To reduce the importance of less interesting genes and be more efficient, it is standard to perform PCA (using highly variable genes), retain only principal components that explain significantly, and then compute tSNE or UMAP into a two-dimensional representation, taking the first (usually from 10 to 50), principal components.

Getting to insights: Downstream analyses

As scRNA-seq datasets have become more widely available, athe single-cell community has developed a standard analysis workflow. This pipeline is broad enough in scope to be applied to a wide range of experiments. The number of tools has rapidly increased over the last few years (Fig. 24) and are diverse enough to meet all the field's current needs (Fig. 25). As the number of cells and modalities surveyed keeps on increasing, an active survey of tools being developed seems like good practice. PROGENy estimates the activity of signaling pathways by combining corresponding gene sets with a linear model. On the other hand, DoRothEA is a collection of resources of TF's targets (called regulons) that can serve as gene sets for TF activity inference. Some tools were specifically designed for application on scRNA-seq data SCENIC [START_REF] Aibar | SCENIC: single-cell regulatory network inference and clustering[END_REF]. SCENIC is a computational workflow that comprises the construction of gene regulatory networks (GRNs) from scRNA-seq data that are subsequently interrogated to infer TF activity with the statistical method AUCell. A recent benchmark noted that the overall performance of these approaches was less than ideal. [START_REF] Pratapa | Benchmarking algorithms for gene regulatory network inference from single-cell transcriptomic data[END_REF]. One possible reason is that single-cell RNA-seq techniques may not still provide sufficient resolution and variation in expression for the reliable inference of GRNs despite rapid advances both in the number of cells and coverage. objective was to construct high-quality datasets to work with. This was done mainly using lenient cell calling and stringent quality control, ensuring that we worked with the highest number of high-quality cells. This step is among the most important since single-cell analyses are also subject to the "garbage in, garbage out" concept. When multiple samples or conditions were present, it was also critical to assess sample integration. Integration was examined by comparing single samples annotation to the merged samples, ensuring no bias was introduced.

When possible, integration was also closely monitored with scTCR-seq information used as a ground truth. scTCR-seq was able to detect when Seurat v3 integration was badly merging different datasets. This happened when merging contrasting datasets, like blood and tumor samples, especially when there was an unbalance of the number of cells between the conditions. showing an over smoothed picture after batch-correction. In this case, we chose to rely on supervised reference-based annotation and not on unsupervised clustering, which led to a much clearer picture. System using a High Sensitivity DNA chip. At the end, indexed libraries were tested for quality, equimolarly pooled and sequenced on an Illumina HiSeq2500 using paired-end 26x98bp as sequencing mode, using a full Rapid flow cell, with a coverage around 75M reads per sample.

Approximately 500-1,000 cells were obtained corresponding to 86,000 reads/cell. Each pool of cells was tested for library quality and concentration.

scRNA-seq analysis.

Single-cell sequencing files (basecalls) were processed using the Cell Ranger Single Cell Software Suite (v1.3.1) to perform quality control, sample de-multiplexing, barcode processing, and single-cell 3′gene counting (http://software.10xgenomics.com/singlecell/overview/welcome) (35). Samples were first demultiplexed and then aligned to the UCSC mouse (mm10) transcriptome and genome using "cellranger mkfastq" with default parameters for both WT and Suv39h1-KO samples. UMI were counted using "cellranger count". Each sample, was run twice, WT and Suv39h1-KO naive and Kb -OVA+ CD8+ T were merged using "cellranger aggregate" procedure. After UMI counting, 4,959 WT (naive, and two technical replicates) and 4,074 Suv39h1-KO (naive, and two technical replicates) were analyzed. For the second biological replicates: 1,351 wt and 1,072 Suv39h1-KO. Approximately 36,000 and 43,000 raw reads per cell were obtained in average, respectively. Further analysis was performed in R (v3.3) using the Seurat package (v1.4.0.14) (36). The gene-cell-barcode matrix of the samples was log-transformed and filtered based on the number of genes detected per cell (any cell with less than 400 genes per cell was filtered out). Regression in gene expression was performed based on the number of UMI and the percentage of mitochondrial genes. Only those genes detected in at least ten cells were included in the analysis. Cells were then scaled to a total of 1 × e4 molecules (37). Any cell with more than 30% of mitochondrial UMI counts was filtered out. For the first biological replicate, 4,619 (WT) and 3,946 (Suv39h1-KO) cells were kept for statistical analysis. To reduce data dimensionality, 2,663 (WT) and 2,385 (Suv39h1-KO) variable genes were selected based on their expression and dispersion (expression cutoff = 0, and dispersion cutoff = 0.5), PCA was run on the normalized gene-barcode matrix. Using the Elbow Plot 7 approach, the first 15 principal components were selected. Barnes-hut approximation to tSNE was then performed on the selected principal components to visualize cells in a twodimensional space used for the t-SNE projection and clustering analysis. Clusters were identified using the "Find_Clusters" function in Seurat with default parameters and a resolution of 0.8 for the first biological replicate and 1,0 for the second one. This graph-based clustering method relies on a clustering algorithm based on shared nearest neighbor (SNN) modularity optimization. Clusters were grouped into semi-unsupervised immunological categories, according to the differential expression of hallmark genes. Genes that were differentially expressed between clusters were identified by running the Seurat "Find_All_Markers" function using default parameters. Sensitivity DNA chip. Indexed libraries were tested for quality, equimolarly pooled and sequenced on an Illumina HiSeq2500 using paired-end 26x98bp as sequencing mode. By using a Unique cluster-specific genes were identified by running the Seurat "Find_All_Markers" function using the MAST framework (Version 3alpha) (Finak et al., 2015). Contaminating cells were removed using the expression of Tyrobp, H2-Eb1, Apoe and Hbb-b.

Signature scores were computed using the Seurat function "AddModuleScore" using the gene signature of interest. This function calculates for each individual cell the average expression of each gene signature, subtracted by the aggregated expression of control gene sets. All analyzed genes are binned into 25 bins based on averaged expression, and for each gene of the gene signature, 100 control genes are randomly selected from the same bin as the gene. Featureplots were plotted using minimum and maximum cutoff values for each feature were respectively quantile 3 and quantile 97. Cells were then selected using exhaustion signature score. Top 30% cells on the whole dataset were selected for Fig. 6. Trajectory analysis were performed using Monocle3 using proposed parameters [START_REF] Qiu | Single-cell mRNA quantification and differential analysis with Census[END_REF]. UMAP was used for dimension reduction and SimplePPT was used to learn the graph. Data is available at GEO (accession numbers XXX). Scripts used to perform this analysis are available on GitHub (https://github.com/p-gueguen/Suv39h1_aPD1).

Appendix 3: Gueguen et al. 2021

Materials and Methods

Single-Cell RNA sequencing and TCR (VDJ) profiling

All samples for a given donor (blood, tumor and normal adjacent lung tissue) were processed simultaneously with the Chromium Controller (10X Genomics) and the resulting libraries were prepared in parallel in a single batch. After Gel Bead-in-Emulsion reverse transcription (GEM-RT) reaction and clean-up, a total of 14 cycles of PCR amplification were performed to obtain sufficient cDNAs used for both RNA-seq library generation and TCR V(D)J targeted enrichment followed by V(D)J library generation. TCR V(D)J enrichment was done per manufacturer's user guide using Chromium Single Cell V(D)J Enrichment Kit, Human T cell (10X Genomics). cDNA before and after TCR enrichment was profiled using both Qubit (Thermofischer scientific) and Bioanalyzer High Sensitivity DNA kit (Agilent Technologies). Libraries for RNA-seq and V(D)J were prepared following the manufacturer's user guide (10X Genomics), then profiled using Kapa Library Quantification kit (Kapa Biosystems) and quantified with Qubit (Thermo Fisher Scientific). Single-cell 3' RNA-seq libraries were sequenced on an HiSeq2500 (Illumina).

Single-cell 5' RNA-seq libraries were sequenced on an Novaseq600 (Illumina). Single-cell TCR V(D)J libraries were sequenced either on an MiSeq (Illumina) or Hiseq2500 (Illumina), according to the number of samples processed each time. All sequencing was done according to the manufacturer's specification (10X Genomics).

Flow cytometry and antibodies

Surface markers staining was performed by using a combination of monoclonal antibodies from the following: CD3-AF532 (1:50, clone UCHT1, #2067979, Thermofischer Scientific), CD8-BV510

Label transfer was performed using Seurat v3 FindTransferAnchors and TransferData functions.

Briefly FindTransferAnchors finds a set of anchors between a reference and query object while

TransferData uses the set of anchors computed to transfer categorical or continuous data from the reference to the query object. For each label transfer, default parameters were used.

Namely 30 PCs were selected for FindTransferAnchors and TransferData functions. For the label transfer between 11 tumors and 4 blood and 2 juxta samples, the uncorrected original matrix was used in the query object since we did label transfer directly to classify our juxta-tumor and blood cells. We chose to term "unassigned" cells that had a prediction score below the bottom 5% of the positive control represented by the merging of 4 tumors in the matching 11 tumors reference. For the label transfer between cycling / non-cycling cells, and for label transfer between CD8-GZMH and the rest of the CD8 subsets, the integrated matrix was used for both the reference and the query objects.

Denoising using ENHANCE

We used ENHANCE v0.1 (42), a denoising strategy based on PCA for visualization when genes were dropped-out (TOX and TCF7). This algorithm is base PCA to separate biological heterogeneity from technical noise and then mitigating the bias towards highly expressed genes by aggregating closely related cells. The algorithm was run with default parameters as described in the preprint (42). Since such denoising algorithms might lead to overcorrection, we still used the original matrix for the rest of the analysis.

Silhouette score

Silhouette widths were computed using the first 30 PCs using the silhouette function from the R package cluster v2.1.0. The silhouette value is a measure of how similar a cell is to its own clusters compared to other clusters. The silhouette ranges from -1 to +1. A high positive value

shows that the cell is well matched to its own cluster and poorly matched to neighboring clusters. A value around 0 indicates that the sample is close to the decision boundary between two neighboring clusters and negative values indicate that those samples might have been wrongly assigned to a cluster. Emptydrops we notice a significant gain of cells retrieved, even after subsequent drastic quality control. This means that cells retrieved were most likely of good quality.

Supplementary figures

TCR clone merging

Assessment of TCR clones as a connected graph. Each node represents one clonotype, defined by a TRA and/or a TRB sequence. A link is drawn when either TRA or TRB sequences are shared between two clonotypes. The first plot shows that clusters of TCR appear if we were to apply a greedy TCR merging approach. After showing the link with only the TRB sequences, the picture is much. Merging with TRB sequence can be reasoned because there are ~4x more TRB reads than TRA, and TRB recombination happens before TRA recombination.

Assessing merging (Seurat v2 vs v3) with TCR clones as positive control

Assessing Seurat v3 merging approach, considering TCR information as a ground truth. On top, we merge only 4 samples of T cells (two patients with matched tumor and blood samples). On the bottom, we merge 11 samples from 9 patients (9 tumors and 2 blood samples). We can see that the unbalance between conditions and/or patients led to an over smoothing of the data on the bottom. TCR clones were nicely segregated in their respective cluster on the top configuration, while TCR clones were scattered all over, and not cluster-specific anymore on the bottom representation.